Perspective: Stochastic algorithms for chemical kinetics

Journal of Chemical Physics

Volumn 138, Issue 17, 2013, Pages

  Gillespie, Daniel T ^a   Hellander, Andreas ^b   Petzold, Linda R ^b  

^a Dan T Gillespie Consulting   (United States)

^b University of California   (United States)

Author keywords

[No Author keywords available]

Indexed keywords

NUMERICAL SIMULATION ALGORITHMS; RARE EVENT; REACTANT MOLECULES; STIFF SYSTEMS; STOCHASTIC ALGORITHMS; STOCHASTIC CHEMICAL KINETICS; SUBVOLUMES;

ALGORITHMS; KINETICS; ORDINARY DIFFERENTIAL EQUATIONS;

STOCHASTIC SYSTEMS;

ALGORITHM; ARTICLE; CHEMICAL MODEL; COMPUTER SIMULATION; DIFFUSION; KINETICS; STATISTICS;

ALGORITHMS; COMPUTER SIMULATION; DIFFUSION; KINETICS; MODELS, CHEMICAL; STOCHASTIC PROCESSES;

EID: 84877749614     PISSN: 00219606     EISSN: None     Source Type: Journal    
DOI: 10.1063/1.4801941     Document Type: Article

References (110)

1. L. Wilhelmy, Ann. Phys. Chem. 81, 413 (1850);
2. L. Wilhelmy, Ann. Phys. Chem. 81, 499 (1850).
3. M. Delbrück, J. Chem. Phys. 8, 120 (1940). 10.1063/1.1750549
4. D. McQuarrie, J. Appl. Probab. 4, 413 (1967). 10.2307/3212214
5. H. McAdams and A. Arkin, Proc. Natl. Acad. Sci. U.S.A. 94, 814 (1997). 10.1073/pnas.94.3.814
6. A. Arkin, J. Ross, and H. McAdams, Genetics 149, 1633 (1998).
7. M. Elowitz, A. Levine, E. Siggia, and P. Swain, Science 297, 1183 (2002). 10.1126/science.1070919
8. L. Weinberger, J. Burnett, J. Toettcher, A. Arkin, and D. Schaffer, Cell 122, 169 (2005). 10.1016/j.cell.2005.06.006
9. 0.
10. D. Gillespie, J. Comput. Phys. 22, 403 (1976). 10.1016/0021-9991(76) 90041-3
11. D. Gillespie, J. Phys. Chem. 81, 2340 (1977). 10.1021/j100540a008
12. BT), which is the famous Arrhenius factor. 10.1016/0378-4371(92)90283-V
13. D. Gillespie, J. Chem. Phys. 131, 164109 (2009). 10.1063/1.3253798
14. As discussed in Sec. VI of this reference, the analysis leading to the result can be regarded as a refined, corrected, and stochastically extended version of the analysis of F. Collins and G. Kimball, J. Colloid Sci. 4, 425 (1949). 10.1016/0095-8522(49)90023-9
15. A derivation of Eq. that is further improved over the one given in the 2009 paper can be found in Secs. 3.7 and 4.8 of D. Gillespie and E. Seitaridou, Simple Brownian Diffusion (Oxford University Press, 2012).
16. M. Smoluchowski, Z. Phys. Chem. 92, 129 (1917).
17. T. Nakanishi, J. Phys. Soc. Jpn. 32, 1313 (1972); 10.1143/JPSJ.32.1313
18. T. Nakanishi, J. Phys. Soc. Jpn. 40, 1232 (1976). 10.1143/JPSJ.40.1232
19. M. Šolc and I. Horsák, Collect. Czech. Chem. Commun. 37, 2994 (1972); 10.1135/cccc19722994
20. M. Šolc and I. Horsák, Collect. Czech. Chem. Commun. 38, 2200 (1973); 10.1135/cccc19732200
21. M. Šolc and I. Horsák, Collect. Czech. Chem. Commun. 40, 321 (1975). 10.1135/cccc19750321
22. D. Bunker, B. Garrett, T. Kleindienst, and G. Long III, Combust. Flame 23, 373 (1974). 10.1016/0010-2180(74)90120-5
23. 0(x) reveals that τ and j are statistically independent random variables with those respective PDFs. Integrating (summing) those PDFs over their respective arguments gives their CDFs, and Eqs. then follow from the inversion method. The generic forms of Eqs., being such straightforward consequences of applying the inversion Monte Carlo method to the aforementioned PDFs of τ and j, were known long before Gillespies 1976 paper; e.g., see page 36 of J. Hammersley and D. Handscomb, Monte Carlo Methods (Methuen, 1964). The main contributions of Gillespies 1976 paper were: (i) proving from simple kinetic theory that bimolecular reactions in a dilute gas, like unimolecular reactions, are describable by propensity functions as defined in; and (ii) proving that propensity functions as defined in imply that the time to next reaction and the index of next reaction are random variables distributed according to the joint PDF.
24. M. Gibson and J. Bruck, J. Phys. Chem. 104, 1876 (2000). 10.1021/jp993732q
25. D. Gillespie, Annu. Rev. Phys. Chem. 58, 35 (2007). 10.1146/annurev. physchem.58.032806.104637
26. Y. Cao, H. Li, and L. Petzold, J. Chem. Phys. 121, 4059 (2004). 10.1063/1.1778376
27. J. McCollum, G. Peterson, C. Cox, M. Simpson, and N. Samatova, Comput. Biol. Chem. 30, 39 (2006). 10.1016/j.compbiolchem.2005.10.007
28. D. Anderson, J. Chem. Phys. 127, 214107 (2007). 10.1063/1.2799998
29. T. Kurtz, Ann. Probab. 8, 682 (1980); 10.1214/aop/1176994660
30. S. Ethier and T. Kurtz, Markov Processes: Characterization and Convergence (Wiley, 1986).
31. D. Anderson and T. Kurtz, Continuous time Markov chain models for chemical reaction networks., in Design and Analysis of Biomolecular Circuits: Engineering Approaches to Systems and Synthetic Biology, edited by, H. Koeppl, (Springer, 2011).
32. A. Slepoy, A. Thompson, and S. Plimpton, J. Chem. Phys. 128, 205101 (2008). 10.1063/1.2919546
33. S. Mauch and M. Stalzer, IEEE/ACM Trans. Comput. Biol. Bioinf. 8, 27 (2011). 10.1109/TCBB.2009.47
34. D. Gillespie, J. Chem. Phys. 115, 1716 (2001). 10.1063/1.1378322
35. An early conjecture was that the negativity problem in tau-leaping was caused by the unbounded Poisson random variables in the tau-leaping formula occasionally allowing too many reaction firings in a single leap. That conjecture led to proposals to replace the Poisson random variables in Eq. with binomial random variables, which are bounded. But it was subsequently determined that the principal causes of the negativity problem lay elsewhere. First was the flawed implementation of the leap condition that was used in Ref.: the small-valued propensity functions that it failed to protect often have a reactant with a small molecular population which can easily be driven negative. Second, since the firing numbers of the individual reactions in the tau-leaping formula are generated independently of each other, two or more reaction channels that decrease the population of a common species could inadvertently collude to overdraw that species. Neither of these two problems is fixed by the ad hoc substitution of binomial random variables for the Poisson random variables in Eq.. But both problems are effectively dealt with in the heavily revised tau-leaping procedure described in Ref., which uses the theoretically appropriate Poisson random variables.
36. Y. Cao, D. Gillespie, and L. Petzold, J. Chem. Phys. 124, 044109 (2006). 10.1063/1.2159468
37. D. Gillespie, Simulation methods in systems biology., in Formal Methods for Computational Systems Biology, edited by, M. Bernardo, P. Degano, and, G. Zavattaro, (Springer, 2008), pp. 125-167, Sec. 3 gives a tutorial on tau-leaping, Secs. 4 and 5 give a tutorial on the ssSSA.
38. M. Rathinam, L. Petzold, Y. Cao, and D. Gillespie, J. Chem. Phys. 119, 12784 (2003); 10.1063/1.1627296
39. M. Rathinam, L. Petzold, Y. Cao, and D. Gillespie, J. Chem. Phys. 121, 12169 (2004); 10.1063/1.1823412
40. M. Rathinam, L. Petzold, Y. Cao, and D. Gillespie, Multiscale Model. Simul. 4, 867 (2005). 10.1137/040603206
41. A. Auger, P. Chatelain, and P. Koumoutsakos, J. Chem. Phys. 125, 084103 (2006). 10.1063/1.2218339
42. D. Anderson, J. Chem. Phys. 128, 054103 (2008). 10.1063/1.2819665
43. t in Eq. and Gaussian white noise, see D. Gillespie, Am. J. Phys. 64, 225 (1996); 10.1119/1.18210
44. D. Gillespie, Am. J. Phys. 64, 1246 (1996). 10.1119/1.18387
45. D. Gillespie, J. Chem. Phys. 113, 297 (2000). The CLE can be shown to be mathematically equivalent to the Fokker-Planck equation that is obtained by first making a formal Taylor series expansion of the right side of the CME, which yields the so-called Kramers-Moyal equation, and then truncating that expansion after the second derivative term. However, that way of obtaining the CLE, which was known well before this 2000 paper, does not qualify as a derivation, because it does not make clear under what conditions the truncation will be accurate. In contrast, the derivation of the CLE given here provides a clear and testable criterion for the accuracy of its approximations, namely, the extent to which both leap conditions are satisfied. But see also the proviso in Ref.. 10.1063/1.481811
46. The approximation P(m)≈N(m,m) that was made in deriving the CLE from the tau-leaping formula, while accurate for likely values of those two random variables, is very inaccurate in the tails of their probability densities, even when m 1. Although both tails are very near 0, they differ by many orders of magnitude. Since rare events arise from the unlikely firing numbers under those tails, it follows that the CLE will not accurately describe the atypical behavior of a chemical system, even if both leap conditions are well satisfied.
47. D. Gillespie, J. Phys. Chem. B 113, 1640 (2009). 10.1021/jp806431b
48. The foregoing chain of reasoning can be summarized from the perspective of computational mathematics as follows: The definition of the propensity function implies, for any time step τ that is small enough to satisfy the first leap condition, the tau-leaping formula. If the second leap condition is also satisfied, the tau-leaping formula becomes Eq., which is the forward Euler formula for a stochastic differential equation. And that formula becomes in the thermodynamic limit, where its diffusion term will be negligibly small compared to its drift term, the forward Euler formula for an ordinary differential equation.
49. A common misconception is that, while the molecular population X is obviously a discrete variable, the molecular concentration Z XΩ is a continuous variable. The error in that view becomes apparent when one realizes that simply by adopting a unit of length that gives Ω the value 1, Z becomes numerically equal to X. But even if that is not done, a sudden change in X, say from 10 to 9, will always result a discontinuous 10 decrease in Z. The molecular concentration Z is no less discrete, and no more continuous, than the molecular population X.
50. T. Kurtz, Stochastic Proc. Appl. 6, 223 (1978). 10.1016/0304-4149(78) 90020-0
51. N. van Kampen, Adv. Chem. Phys. 34, 245 (1976); 10.1002/9780470142530.ch5
52. Stochastic Processes in Physics and Chemistry (North-Holland, 1992).
53. E. Wallace, D. Gillespie, K. Sanft, and L. Petzold, IET Syst. Biol. 6, 102 (2012). 10.1049/iet-syb.2011.0038
54. R. Grima, P. Thomas, and A. Straube, J. Chem. Phys. 135, 084103 (2011). 10.1063/1.3625958
55. jX(t′)dt′, see Refs. and. Note that all the integrals here can be written as finite algebraic sums, since X(t′) stays constant between successive reactions.
56. Y. Cao, D. Gillespie, and L. Petzold, J. Chem. Phys. 122, 014116 (2005); 10.1063/1.1824902
57. Y. Cao, D. Gillespie, and L. Petzold, J. Chem. Phys. 123, 144917 (2005). 10.1063/1.2052596
58. A more concise tutorial presentation of the ssSSA, featuring illustrative applications to reactions and the classical enzyme-substrate reactions E S ES → E P, is given in Ref.. A more thorough examination of the relation between the ssSSAs treatment of the enzyme-substrate system and the classical Michaelis-Menten approximation is given in K. Sanft, D. Gillespie, and L. Petzold, IET Syst. Biol. 5, 58 (2011). 10.1049/iet-syb.2009.0057
59. 1 must be also. This illustrates the important fact that fast and slow reaction channels cannot always be identified solely on the basis of the magnitudes of their rate constants.
60. Y. Cao, D. Gillespie, and L. Petzold, J. Comput. Phys. 206, 395 (2005). 10.1016/j.jcp.2004.12.014
61. L. Ferm, P. Lotstedt, and A. Hellander, J. Sci. Comput. 34, 127 (2008). 10.1007/s10915-007-9179-z
62. A. Singh and J. Hespanha, IEEE Trans. Autom. Control. 56, 414 (2011). 10.1109/TAC.2010.2088631
63. W. E, D. Liu, and E. Vanden-Eijnden, J. Chem. Phys. 123, 194107 (2005); 10.1063/1.2109987
64. D. Gillespie, L. Petzold, and Y. Cao, J. Chem. Phys. 126, 137101 (2007). 10.1063/1.2567036
65. K. Sanft, Ph.D. thesis, University of California, Santa Barbara, 2012.
66. H. Kuwahara and I. Mura, J. Chem. Phys. 129, 165101 (2008). 10.1063/1.2987701
67. D. Gillespie, M. Roh, and L. Petzold, J. Chem. Phys. 130, 174103 (2009). 10.1063/1.3116791
68. M. Roh, D. Gillespie, and L. Petzold, J. Chem. Phys. 133, 174106 (2010). 10.1063/1.3493460
69. B. Daigle, Jr., M. Roh, D. Gillespie, and L. Petzold, J. Chem. Phys. 134, 044110 (2011). 10.1063/1.3522769
70. M. Roh, B. Daigle, Jr., D. Gillespie, and L. Petzold, J. Chem. Phys. 135, 234108 (2011). 10.1063/1.3668100
71. R. Rubinstein and D. Kroese, The Cross-Entropy Method (Springer, 2004).
72. M. Rathinam, P. Sheppard, and M. Khammash, J. Chem. Phys. 132, 034103 (2010). 10.1063/1.3280166
73. D. Anderson, SIAM J. Numer. Anal. 50, 2237 (2012). 10.1137/110849079
74. K. Takahashi, S. Tǎnase-Nicola, and P. ten Wolde, Proc. Natl. Acad. Sci. U.S.A. 107, 2473 (2010). 10.1073/pnas.0906885107
75. K. Doubrovinski and M. Howard, Proc. Natl. Acad. Sci. U.S.A. 102, 9808 (2005). 10.1073/pnas.0500529102
76. J. Elf and M. Ehrenberg, Syst. Biol. 1, 230 (2004). 10.1049/sb:20045021
77. R. Metzler, Phys. Rev. Lett. 87, 068103 (2001). 10.1103/PhysRevLett.87. 068103
78. M. Sturrock, A. Hellander, M. Matzavinos, and M. Chaplain, J. R. Soc. Interface 10, 20120988 (2013). 10.1098/rsif.2012.0988
79. D. Fange and J. Elf, PLoS Comput. Biol. 2 (6), e80 (2006). 10.1371/journal.pcbi.0020080
80. C. Gardiner, K. McNeil, D. Walls, and I. Matheson, J. Stat. Phys. 14, 307 (1976). 10.1007/BF01030197
81. See, for example, D. Gillespie and E. Seitaridou, Simple Brownian Diffusion (Oxford University Press, 2012), Chap. 5. That chapters revised Sec. 5.6 (which can be downloaded from the books webpage on the publishers website) shows that below a certain value of h, no propensity function can give a physically accurate modeling of diffusive transfers of molecules between voxels.
82. Equation also follows by applying the centered finite difference method. in numerical analysis to the diffusion Eq..
83. T. Kurtz, J. Appl. Probab. 7, 49 (1970). 10.2307/3212147
84. S. Isaacson and C. Peskin, SIAM J. Sci. Comput. 28, 47 (2006). 10.1137/040605060
85. S. Isaacson, D. McQueen, and C. Peskin, Proc. Nat. Acad. Sci. U.S.A. 108, 3815 (2011). 10.1073/pnas.1018821108
86. S. Engblom, L. Ferm, A. Hellander, and P. Lötstedt, SIAM J. Sci. Comput. 31, 1774 (2009). 10.1137/080721388
87. B. Drawert, S. Engblom, and A. Hellander, BMC Syst. Biol. 6, 76 (2012). 10.1186/1752-0509-6-76
88. H.-W. Kang, L. Zheng, and H. Othmer, J. Math. Biol. 65, 1017 (2012). 10.1007/s00285-011-0469-6
89. D. Rossinelli, B. Bayati, and P. Koumoutsakos, Chem. Phys. Lett. 451, 136 (2008). 10.1016/j.cplett.2007.11.055
90. K. Iyengar, L. Harris, and P. Clancy, J. Chem. Phys. 132, 094101 (2010). 10.1063/1.3310808
91. T. Marquez-Lago and K. Burrage, J. Chem. Phys. 127, 104101 (2007). 10.1063/1.2771548
92. S. Lampoudi, D. Gillespie, and L. Petzold, J. Chem. Phys. 130, 094104 (2009). 10.1063/1.3074302
93. W. Koh and K. Blackwell, J. Chem. Phys. 134, 154103 (2011). 10.1063/1.3572335
94. B. Drawert, M. Lawson, L. Petzold, and M. Khammash, J. Chem. Phys. 132, 074101 (2010). 10.1063/1.3310809
95. L. Ferm, A. Hellander, and P. Lötstedt, J. Comput. Phys. 229, 343 (2010). 10.1016/j.jcp.2009.09.030
96. D. Gillespie, S. Lampoudi, and L. Petzold, J. Chem. Phys. 126, 034302 (2007). 10.1063/1.2424461
97. S. Lampoudi, D. Gillespie, and L. Petzold, J. Comput. Phys. 228, 3656 (2009). 10.1016/j.jcp.2009.02.002
98. R. Grima, J. Chem. Phys. 132, 185102 (2010). 10.1063/1.3427244
99. J. Zon, S. van Zon, and P. Rein ten Wolde, Phys. Rev. Lett. 94, 128103 (2005). 10.1103/PhysRevLett.94.128103
100. R. Erban and J. Chapman, Phys. Biol. 6, 046001 (2009). 10.1088/1478-3975/6/4/046001
101. D. Fange, O. G. Berg, P. Sjöberg, and J. Elf, Proc. Natl. Acad. Sci. U.S.A. 107, 19820 (2010). 10.1073/pnas.1006565107
102. S. Isaacson, SIAM J. Appl. Math. 70, 77 (2009). 10.1137/070705039
103. S. Isaacson and D. Isaacson, Phys. Rev. E 80, 066106 (2009). 10.1103/PhysRevE.80.066106
104. S. Hellander, A. Hellander, and L. Petzold, Phys. Rev. E 85, 042901 (2012). 10.1103/PhysRevE.85.042901
105. S. Isaacson, A convergent reaction-diffusion master equation., preprint arXiv:1211.6772v1 (2012).
106. M. Doi, J. Phys. A: Math. Gen. 9, 1465 (1976); 10.1088/0305-4470/9/9/008
107. M. Doi, J. Phys. A: Math. Gen. 9, 1479 (1976). 10.1088/0305-4470/9/9/009
108. A. Hellander, S. Hellander, and P. Lötstedt, Multiscale Model. Simul. 10, 585 (2012). 10.1137/110832148
109. M. Flegg, S. Chapman, and R. Erban, J Roy. Soc. Interface 9, 859 (2012). 10.1098/rsif.2011.0574
110. K. Sanft, S. Wu, M. Roh, J. Fu, R. Lim, and L. Petzold, Bioinformatics 27, 2457 (2011). 10.1093/bioinformatics/btr401