(M-097) Automatic differentiation, a possible solution to numerical instability in Population PK parameter estimation

Rong Chen, PhD – Research Scientist, Certara; James Craig, MS – Software engineer, Certara; Shuhua Hu, PhD – Principal Research Scientist, Certara; Alex Mazur, BS – engineer, Certara; Keith Nieforth, PharmD – Sr. Director, Product Management, Certara; Michael Tomashevskiy, MD – engineer, Certara

Objective: Evaluate Automatic Differentiation (AD) vs Finite Difference (FD) for fitting population PK (PopPK)

Methods: Imprecision in estimation of PopPK model may result in numerical problems such as failure to converge to a sufficiently precise minimum and non-positive definite covariance matrices as well as slower convergence. The standard method for parameter estimation in PopPK FOCE models is quasi Newton, typically BFGS [1], with gradients estimated by FD [2]. We propose AD [3] as an alternative to FD for the “inner loop”, the estimation of ETA effects. AD may provide a more precise estimate of the gradient than FD, potentially resulting in both improved speed and precision. To assess the performance of FD vs AD methods, simulated data sets using ordinary differential equation (ODE) models with saturable elimination were developed, consisting of 60 subjects, 7 samples each. Models (number of models) included:

• First order absorption (48)


• 1 (24) or 2 (24) compartment


• Volume an allometric function of weight (exponent 1.1) (24) or independent of weight (24)


• Without the gamma parameter for the saturable term (24) or fixed to a value of 1.6 (24)

OMEGA structures (block [32] and diagonal [16]) with 1 (8), 2 (8), 3 (16), 4 (8) or 6 (8) estimated OMEGA values. The “true” models (those used for simulation) were fitted using NLME with both FD (NLME FD) and AD (NLME AD) methods and the DVERK ODE solver. For reference, the same “true” models were also fitted in NONMEM 7.5, with ADVAN6 (DVERK ODE solver) and the same precision (TOL=6). Identical initial estimates (changed ~20% from simulation values) were used for all models. Evaluation criteria included:

• Estimation time (ESTtime, seconds)


• Covariance time (COVtime, seconds)


• Fraction with successful convergence (fSuccess, %)


• Fraction with successful covariance (fCOV, %)

All benchmarks were run on the same 4 core computer using a single core (no parallel execution) and was rebooted between runs. No other computationally intensive software was running. during the benchmark runs.

Results:	NLME FD	NLME AD	NONMEM
GeoMean ESTtime	133	265	405*
GeoMean COVtime	15.1	34.0	274*
fSuccess	98	100	75*
fCOV	100	100	88*

*p < 0.05
NLME AD converged in all cases, compared to 47 of 48 for NLME FD (NS) and 75% (36 of 48) for NONMEM (p < 0.05). NLME FD and NLME AD had successful covariance steps in all cases, compared to 88% (42 of 48, p< 0.05) for NONMEM. A larger sample size with more models that fail convergence or covariance step may be necessary for an adequate statistical test for FD vs AD.

Conclusions: CPU time for fitting was generally greater for AD than for FD for NLME. Compared to traditional FD methods of calculating the gradient AD may improve numerical stability. NLME FD and NLMD AD were faster than NONMEM for both estimation and covariance. R code for execution of the benchmark runs for NLME FD and NONMEM is posted to https://github.com/certara/nlme-performance-demo. Work is ongoing to improve execution speed.

Citations: [1] Shanno D, “Conditioning of quasi-Newton methods for function minimization”, Math. Comp. 24, 647–656 (1970)

[2] Broyden, C. G. (1972). "Quasi-Newton Methods". In Murray, W. (ed.). Numerical Methods for Unconstrained Optimization. London: Academic Press. pp. 87–106. ISBN 0-12-512250-0.

[3] Baydin AG, Pearlmutter BA, Radul AA, Siskind JM. “Automatic differentiation in machine learning: a survey” The Journal of Machine Learning Research, 18(153):1--43, (2018)