Anna Mikhailova 1,2, Huijuan Hu 3, Yuan Xiong 4, Kirill Peskov 1,5, Wenping Wang 3, Victor Sokolov 1,2
1 M&S Decisions (Dubai, United Arab Emirates), 2 Marchuk Institute of Numerical Mathematics of the Russian Academy of Sciences (INM RAS), (Moscow, Russia), 3 Takeda Development Center Americas, Inc (Cambridge, USA), 4 Johnson & Johnson (Raritan, USA), 5 Research Center of Model-Informed Drug Development, Sechenov First Moscow State Medical University (Moscow, Russia)
Introduction
Most model-based analyses in pharmacometrics (PMx) and the associated statistical inference follow the frequentist paradigm, established in the early XX century by Ronald Fisher in opposition to Bayesian statistics, which was considered too subjective due to its reliance on priors [1]. A limitation in one context can be an advantage in another: the Bayesian approach naturally quantifies uncertainty by incorporating prior knowledge – which is the essence of quantitative systems pharmacology (QSP) modeling [2]. Application of the Bayesian approach to QSP or MBMA, however, remains limited because few tools provide out-of-the-box support for Bayesian inference in custom systems of ordinary differential equations (ODEs) with dosing events. Several R packages have been developed to address this gap, including MCSim, NIMBLE, Torsten, and stanette. In this work, we present a Bayesian workflow for aggregated data analysis and demonstrate its application through the development of a generalized dapagliflozin pharmacokinetic model, implemented and benchmarked across these four packages.
Methods
A generalized ODE-based dapagliflozin PK model was implemented in Torsten (v0.89.0), NIMBLE (v1.1.0), MCSim (v6.2.0), and stanette (v2.21.4), following case example described in [6]. Stanette package was updated to ensure compatibility with recent versions of StanHeaders (>2.36.0), Rcpp (>= 1.0.7), RcppParallel (>= 5.1.4) and RcppEigen (>= 0.3.3.9.3) libraries. Functionality was extended to include transit-compartment absorption models, optional exclusion of random effects for selected parameters, and support for log-additive relative error models.
Performance was assessed using both statistical and computational criteria. Statistical criteria included posterior parameter estimates, 95% credible intervals (CrIs), convergence diagnostics (R̂), and effective sample size (ESS), evaluated in both central and tail posterior regions. Predictive adequacy was assessed using the Watanabe–Akaike information criterion (WAIC) [3] and standard graphical diagnostics (observed-versus-predicted plots, time-course profiles, and residual distributions). Computational performance was evaluated based on runtime per iteration per chain. Sampling efficiency was interpreted in relation to ESS and runtime.
Results
The 95% credible intervals (CrI) obtained with stanette overlapped with those produced by MCSim, NIMBLE and Torsten for most of the parameters. Estimated posterior means (95% CrI) were: ka = 5.49 1/h (4.81, 6.03), ktr = 5.30 1/h (3.39, 9.35), Cl = 14.86 L/h (14.37, 15.31), Vd = 65.59 L (58.76, 69.3), Q = 6.82 L/h (5.09, 8.99) and Vp = 62.55 L (47.82, 77.56).
Stanette demonstrated sampling efficiency comparable to Torsten. The mean relative effective sample size (ESS), averaged across all parameters, was 0.58 for the tail and 0.76 in central posterior regions with stanette, compared with 0.75 in the tail region and 0.95 in the central posterior region for Torsten. In contrast, stanette showed superior sampling efficiency relative to MCSim (mean relative ESS: 0.13 for the tail and 0.10 for the central posterior region) and NIMBLE (mean relative ESS: 0.11 for the tail and 0.17 for the central posterior region).
Computational time was slightly higher for stanette (14.9 seconds for one sample in a chain) compared with Torsten (12 seconds for one sample in a chain). The observed difference in computational runtime between Torsten and stanette is primarily due to the integration method used: Torsten employed a specialized solver for linear ODE systems (based on a matrix exponential approach), whereas stanette relied on Stan’s built-in general-purpose ODE solvers. In future iterations of the package, we plan to implement dedicated, model-specific solvers to improve computational efficiency and numerical performance. Other programs have shown better runtime efficacy (MCSim – 0.2 seconds for one sample in a chain, NIMBLE – 0.4 seconds), due to the use of Gibbs sampler.
Conclusion
Stanette produced posterior estimates consistent with established Bayesian ODE-based tools and demonstrated robust convergence properties. Its sampling performance was comparable to Torsten, with moderately lower ESS and slightly longer runtime per iteration. Although MCSim and NIMBLE achieved faster per-iteration runtimes, their sampling efficiency was lower in this application. Overall, Stanette provides a flexible and competitive framework for Bayesian pharmacometric modeling within a PMx-oriented workflow, with performance characteristics similar to state-of-the-art Hamiltonian Monte Carlo–based tools.
References:
References
[1] A. Gelman, J. B. Carlin, H. S. Stern, D. B. Dunson, A. Vehtari, and D. B. Rubin, Bayesian data analysis, Third edition. in Texts in statistical science series. Boca Raton London New York: CRC Press, Taylor and Francis Group, 2014.
[2] R. Van De Schoot, D. Kaplan, J. Denissen, J. B. Asendorpf, F. J. Neyer, and M. A. G. Van Aken, ‘A Gentle Introduction to Bayesian Analysis: Applications to Developmental Research’, Child Development, vol. 85, no. 3, pp. 842–860, May 2014, doi: 10.1111/cdev.12169.
[3] S. Watanabe, ‘Asymptotic Equivalence of Bayes Cross Validation and Widely Applicable Information Criterion in Singular Learning Theory’, 2010, arXiv. doi: 10.48550/ARXIV.1004.2316.
[4] A. Vehtari, A. Gelman, D. Simpson, B. Carpenter, and P.-C. Bürkner, ‘Rank-Normalization, Folding, and Localization: An Improved Rˆ for Assessing Convergence of MCMC (with Discussion)’, Bayesian Anal., vol. 16, no. 2, Jun. 2021, doi: 10.1214/20-BA1221.
[5] A. Gelman et al., ‘Bayesian Workflow’, 2020, arXiv. doi: 10.48550/ARXIV.2011.01808.
[6] A. Mikhailova et al., ‘A Bayesian Workflow for Minimal PBPK Modeling: A Case Study with Dapagliflozin‘, PAGE 2025, Abstract #11574
Reference: PAGE 34 (2026) Abstr 11952 [www.page-meeting.org/?abstract=11952]
Poster: Methodology - Other topics