Lucas Bodelle 1, Marc Lavielle 3,4, Donato Teutonico 2, Emmanuelle Comets 1,5
1 Université Paris Cité and Université Sorbonne Paris Nord, Inserm, IAME (75018 Paris, France), 2 Sanofi, Translational Medicine Unit (94400 Vitry sur Seine , France), 3 Inria (91400 Saclay , France), 4 Ecole Polytechnique (91120 Palaiseau , France), 5 Univ Rennes, Inserm, EHESP, Irset (Institut de recherche en santé, environnement et travail) - UMR_S 1085 (F-35000 Rennes, France )
Introduction:
Quantifying the sources of individual variability in drug concentrations through non-linear mixed effect models (NLMEM) has proven instrumental in drug development [1]. Modelling software includes the commercial Nonmem and Monolix software, as well as open-source implementations such as the saemix package [2] for R. saemix implements the SAEM (Stochastic Approximation Expectation Maximization) algorithm and contains tools for estimation, evaluation and diagnostics. Among all available options, saemix has the advantage of being highly flexible since it is primarily coded in R, making it easier to tweak, as shown by a recent extension to joint models [3]. It can also be easily coupled with other software, such as the PKSim software in the OSP suite [4] through the saemixPBPK package [5], which leverages NLMEM estimation.
Objectives:
The current version of the saemix package was designed for single response models and closed-form solutions benefitting from R’s vectorisation approach. We first extended and evaluated saemix to handle multiple responses. The second objective was computational optimisation of the grid approximation introduced in saemixPBPK and used to limit the number of evaluations of the ODE system. The third objective was to evaluate these extensions using simulations.
Methods:
Starting from saemix version 3.4, we extended every function in the source code to handle multi-responses, from estimation tools to diagnostic tools. We performed a simulation study to evaluate the estimation performance, using a closed-form PK/PD model. Simulations were conducted under four different scenarios: rich (N=100, n=10 samples for PK and for PD) or sparse design (N=10, n=6), adding correlation and covariate effects, and increasing inter-individual variability (IIV) and residual error. In total, 1000 datasets were simulated for each scenario. The parameters and their associated standard errors (SE) were estimated using the extended saemix, and we computed the Relative Estimation Error (REE) to assess bias, and the Relative Standard Error (RSE) to verify whether the predicted SE distribution aligns with the empirical SE.
The second objective of this work was to optimise and evaluate the grid interpolation approach in saemixPBPK. To accelerate computation, the source code was improved: critical sections were recoded in C++, and the data structure was modified to drastically reduce memory costs. The grid interpolation method uses a grid error threshold that defines the deviation from linearity, a lower error threshold leads to more accurate estimation but slower computation. We evaluated the impact of error thresholds on estimation performance and computational time, applying thresholds of 1%, 5%, and 10%. This evaluation was performed on 200 datasets from the rich sampling scenario.
Finally, a comparison was conducted on the warfarin dataset [6]to evaluate the number of model evaluations required when estimating 1 to 9 fixed parameters, both with and without grid interpolation. We used deSolve to integrate the PKPD model involving a turnover response, and the analysis was performed using 500 saemix iterations with two Markov chains.
Results:
In the multi-response evaluation approach, the REE showed no bias across all scenarios, except for the variance parameters in the sparse scenario, which was expected to be less informative due to fewer individuals and observations. The RSE were well predicted but slightly underestimated for EC50.
Regarding the evaluation of the grid interpolation approach, the bias in REE increases as the grid error threshold increases. The lower the error threshold, the more accurate the estimation. However, a more precise grid leads to exponentially longer runtime, therefore, the trade-off between approximation accuracy and computational efficiency should be carefully considered. In the settings used for the evaluation, the grid interpolation method proved efficient when estimating 5 or fewer fixed parameters, with respectively 3, 8, 28, 95 and 826-fold reductions in computational time for 5 estimated parameters to 1.
Conclusion:
The saemix package was successfully extended to handle multiple responses and will be available in saemix version 3.5 soon to be released on CRAN. A tutorial notebook will be added to the GitHub repository. Furthermore, the new grid interpolation method demonstrates satisfactory estimation results and can be used to accelerate computation with complex models such as ODE systems or PBPK models when interfaced with PKSim through the saemixPBPK package [5].
Funding Statement:
During this work, Lucas Bodelle was employed by Sanofi.
References:
References:
[1] Lavielle, M. (2014). Mixed effects models for the population approach: models, tasks, methods and tools. CRC press.
[2] Comets, E. et al (2017) Parameter Estimation in Nonlinear Mixed Effect Models Using saemix, an R Implementation of the SAEM Algorithm. Journal of statistical software
[3] Lavalley-Morelle (2024) Extending the code in the open-source saemix package to fit joint models of longitudinal and time-to-event data, Computer Methods and Programs in Biomedicine
[4] Lippert, J., et al. (2019). Open systems pharmacology community—an open access, open source, open science approach to modeling and simulation in pharmaceutical sciences. CPT: pharmacometrics & systems pharmacology
[5] Teutonico, D. (2026) Integrating Population Approaches With Physiologically Based Pharmacokinetic Models: A Novel Framework for Parameter Estimation, CPT: pharmacometrics & systems pharmacology
[6] ROBERT A O’REILLY and Paul M Aggeler (1968) . Studies on coumarin anticoagulant drugs: initiation of warfarin therapy without a loading dose, Circulation
Reference: PAGE 34 (2026) Abstr 12274 [www.page-meeting.org/?abstract=12274]
Poster: Software Demonstration