Abstract
In this paper, we derive a Markov chain Monte Carlo (MCMC) algorithm supported by a neural network. In particular, we use the neural network to substitute derivative calculations made during a Metropolis adjusted Langevin algorithm (MALA) step with inexpensive neural network evaluations. Using a complex, high-dimensional blood coagulation model and a set of measurements, we define a likelihood function on which we evaluate the new MCMC algorithm. The blood coagulation model is a dynamic model, where derivative calculations are expensive and hence limit the efficiency of derivative-based MCMC algorithms. The MALA adaptation greatly reduces the time per iteration, while only slightly affecting the sample quality. We also test the new algorithm on a 2-dimensional example with a non-convex shape, a case where the MALA algorithm has a clear advantage over other state of the art MCMC algorithms. To assess the impact of the new algorithm, we compare the results to previously generated results of the MALA and the random walk Metropolis Hastings (RWMH).
A Appendix: A blood coagulation model
What follows is the mathematical formulation of the ODE system of the blood coagulation model we use in Section 4.1. For a detailed, biological description of the parameters, we refer to [6] and the added supplementary materials.
The initial values of the system are given as
List of fixed parameters (1/2).
Dilutionfactor | 0.6667 | S_TF | 7.658989e | S_VII | 1.561404e | |||
S_VIIa | 5.326683e | S_X | 2.039697e | S_IX | 9.096479e | |||
S_II | 2.017022e | S_VIII | 7.518441e | S_V | 1.886433e | |||
S_TFPI | 2.307784e | S_ATIII | 2.661478e | S_alpha1AT | 2.37676e | |||
S_Tm | 9.79467e | Kd2 | Kd8 | 3.489801e | ||||
Kd9 | 2.426951e | Kd10 | 2.535596e | Kd11 | ||||
Kd11a | 1.602839e | k4 | 26364809 | k6 | 16759672 | |||
k8 | 1.578925 | k9 | 32443880 | k10 | 1.828158 | |||
k11 | 2.486723 | k17 | 15.61378 | k19 | 17.88283 | |||
k21 | 811283301 | k26b | 15.5222 | k27 | 0.3163996 | |||
k28 | 21.5559 | k29 | 136.4629 | k30 | 30080707 | |||
k31 | 104.459 | k32 | 53068894 | k36 | 240193650 | |||
k39 | 88172.98 | Kog10 | Kog8 | 0.2416168 | ||||
Kog8m | 4.0979e | kbu50 | 14.02588 | Va_pre_ActofV | ||||
WarfarinFactor | 0.5978258 | Veer1 | 1537.89 | Down2 | 89.12773 | |||
Bourin1 | 10.3116 | S_HCII | 1.646771e | S_PCI | 8.728542e | |||
S_vWF | 3.707884e | kvWF1 | 583672.2 | S_PL | 2.282184e |
List of fixed parameters (2/2).
0 | 0 | 7.091038e | 1.201881e | 3.287149e | |||||
0 | 9.586261e | 4.885827e | 0 | 5.867748e | |||||
5.268007e | 0 | 1.381180e | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 1.233638e+07 | 1.932636e+06 | |||||
6.253566e | 5.787158e+07 | 5.551706e+07 | 1.350195e+06 | 6.886843e+06 | |||||
4.219812e+06 | 3.486538e | 9.838915e+04 | 3.343605e | 4.907052e+05 | |||||
1.659267e+04 | 2.034925e+07 | 2.463716 | 1.453219e+07 | 3.667853e | |||||
6.113555e+04 | 7.604036e+05 | 5.642855e | 5.850895e | 9.693666e+01 | |||||
2.333402e+04 | 6.103863e | 4.875274e | 2.475065e+06 | 5.122565e+04 | |||||
3.969542e | 8.620074e+05 | 7.423331e | 4.225866e+07 | 2.037410e+03 | |||||
4.950986e+02 | 3.683459e+02 | 2.290561e+02 | 6.109323e+07 | 1.717579e+06 | |||||
0 | 0 | 3.370790e | 5.122001e | 1.745986e | |||||
3.209826e | 3.498565e | 5.679843e+02 | 2.771999e | 1.757158e | |||||
1.783598e | 5.259668e | 3.746471e | 2.113415e | 1.007405e+01 | |||||
5.344848e+01 | 1.395878e+02 | 5.080033e+02 | 1.981534e+05 | 8.094526e+05 | |||||
0 | 5.222449e | 1.380312e+07 | 1.210472e+08 | 7.532941e | |||||
8.460765e | 3.101805e | 7.688052e | 3.924507 | 1.159746e | |||||
2.106957e | 1.242754e+07 | 3.635741e | 1.358760e+07 | 4.299456e | |||||
2.809204e+07 | 4.586732e | 4.930859e | 8.859077 | 4.495039e+06 | |||||
2.208813e | 4.651777e+06 | 1.259095e | 8.571367e+02 | 5.330477e+06 | |||||
1.940336e | 5.136096e | 5.328625e | 0 | 0 | |||||
0 | 0 | 6.951241e+02 | 4.685068e+02 | 8.218390e+01 | |||||
1.208968e+02 | 9.538816e+03 | 3.820788e | 3.268536e+04 | 5.697953e | |||||
8.439507e+02 | 4.546240e | 1.418792e+06 | 1.333866e+06 | 3.369624e | |||||
4.492474e+03 | 7.736602e | 3.181825e+02 | 1.850638e+01 | 1.866119e | |||||
3.006009e | 1.994891e | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 0 | 0 | 0 | 0 | |||||
0 | 1.884206 | 0 | 0 | 2.268764e | |||||
6.094237e+04 | 3.000000e | 2.686606e+04 | 5.890622 | 2.927487e | |||||
0 | 1.616306 | 6.691869e+07 | 9.311524e | 1.575627e+06 | |||||
5.279387e+06 | 1.404290e+02 | 1.741800e+04 | 4.828415e+06 | 5.983991e | |||||
1.879658e+01 | 0 | 0 | 0 | 0 |
And lastly, in Tables 3 and 4, we give a list of values for fixed parameters u and start values for parameters p to be varied, based on [6].
The measurements we use throughout the paper were taken from [15, Figure 7a]. We use the curve corresponding to the case of platelet poor plasma with 10 nM Thrombomodelin initial dosage. The measured output is
The total Thrombin is given as the weighted sum
References
[1] A. Arenas, E. Omodei and M. Brashears, A mechanistic model of human recall of social network structure and relationship affect, Sci. Rep. 7 (2017), Article No. 17133. Search in Google Scholar
[2] F. Baermann, NN-Tool, http://www.nntool.de/Englisch/index_engl.html. Search in Google Scholar
[3] A. R. Barron, Approximation and estimation bounds for artificial neural networks, Mach. Learn. 14.1 (1994), 115–133. 10.1016/B978-1-55860-213-7.50025-0Search in Google Scholar
[4] G. Bellu, M. P. Saccomani, S. Audoly and L. D’Angiò, DAISY: A new software tool to test global identifiability of biological and physiological systems, Comp. Methods Programs Biomed 88 (2007), 52–61. 10.1016/j.cmpb.2007.07.002Search in Google Scholar PubMed PubMed Central
[5] S. Brooks, A. Gelman, G. L. Jones and X.-L. Meng, Handbook of Markov chain Monte Carlo, Chapman & Hall/CRC Handb. Mod. Stat. Methods, CRC Press, Boca Raton, 2011. 10.1201/b10905Search in Google Scholar
[6] R. Burghaus, K. Coboeken, T. Gaub, C. Niederalt, A. Sensse, H.-U. Siegmund, W. Weiss, W. Mueck, T. Tanigawa and J. Lippert, Computational investigation of potential dosing schedules for a switch of medication from warfarin to rivaroxaban–an oral, direct Factor Xa inhibitor, Frontiers Physiol. 5 (2014), Article No. 417. Search in Google Scholar
[7] S. D. Cohen and A. C. Hindmarsh, CVODE, A stiff/nonstiff ODE Solver in C., Comput. Phys. 10 (1996), 138–143. 10.1063/1.4822377Search in Google Scholar
[8] M. K. Cowles and B. P. Carlin, Markov chain Monte Carlo convergence diagnostics: A comparative review, J. Amer. Statist. Assoc. 91 (1996), no. 434, 883–904. 10.1080/01621459.1996.10476956Search in Google Scholar
[9] G. Cybenko, Approximation by superpositions of a sigmoidal function, Math. Control Signals Syst. 2 (1989), no. 4, 303–314. 10.1007/BF02551274Search in Google Scholar
[10] A. Gelman and D. Rubin, Inference from iterative simulation using multiple sequences, Statist. Sci. 7 (1992), 457–472. 10.1214/ss/1177011136Search in Google Scholar
[11] M. Girolami and B. Calderhead, Riemann manifold Langevin and Hamiltonian Monte Carlo methods, J. R. Stat. Soc. Ser. B Stat. Methodol. 73 (2011), no. 2, 123–214. 10.1111/j.1467-9868.2010.00765.xSearch in Google Scholar
[12] A. Griewank and A. Walther, Evaluating Derivatives. Principles and Techniques of Algorithmic Differentiation, 2nd ed., Society for Industrial and Applied Mathematics, Philadelphia, 2008. 10.1137/1.9780898717761Search in Google Scholar
[13] H. Haario, E. Saksman and J. Tamminen, An adaptive Metropolis algorithm, Bernoulli 7 (2001), no. 2, 223–242. 10.2307/3318737Search in Google Scholar
[14] W. K. Hastings, Monte Carlo sampling methods using Markov chains and their applications, Biometrika 57 (1970), no. 1, 97–109. 10.1093/biomet/57.1.97Search in Google Scholar
[15] H. C. Hemker, P. Giesen, R. Al Dieri, V. Regnault, E. de Smedt, R. Wagenvoord, T. Lecompte and S. Béguin, Calibrated automated thrombin generation measurement in clotting plasma, Pathophysiol. Haemost Thromb 33 (2003), 4–15. 10.1159/000071636Search in Google Scholar PubMed
[16] M. Krauss, Using Bayesian-PBPK modeling for assessment of inter-individual variability and subgroup stratification, Silico Pharmacology 1 (2013), 10.1186/2193-9616-1-6. 10.1186/2193-9616-1-6Search in Google Scholar PubMed PubMed Central
[17] M. Krauss, Bayesian population physiologically-based pharmacokinetic (PBPK) approach for a physiologically realistic characterization of interindividual variability in clinically relevant populations, PLoS One 10 (2015), Article No. e0139423. 10.1371/journal.pone.0139423Search in Google Scholar PubMed PubMed Central
[18] C. Müller, F. Weysser, T. Mrziglod and A. Schuppert, Markov-Chain Monte-Carlo methods and non-identifiabilities, Monte Carlo Methods Appl. 24 (2018), no. 3, 203–214. 10.1515/mcma-2018-0018Search in Google Scholar
[19] A. Raue, C. Kreutz, T. Maiwald, J. Bachmann, M. Schilling, U. Klingmüller and J. Timmer, Structural and practical identifiability analysis of partially observed dynamical models by exploiting the profile likelihood, Bioinformatics 25 (2009), no. 15, 1923–1929. 10.1093/bioinformatics/btp358Search in Google Scholar PubMed
[20] G. O. Roberts and J. S. Rosenthal, Optimal scaling for various Metropolis-Hastings algorithms, Statist. Sci. 16 (2001), no. 4, 351–367. 10.1214/ss/1015346320Search in Google Scholar
[21] S. Shammas, A mechanistic model of tau amyloid aggregation based on direct observation of oligomers, Nature Commun. 6 (2015), Article No. 7025. 10.1038/ncomms8025Search in Google Scholar PubMed PubMed Central
© 2020 Walter de Gruyter GmbH, Berlin/Boston