Stochastic simulation is a widely used method for estimating quantities in models of chemical reaction networks where uncertainty plays a crucial role. However, reducing the statistical uncertainty of the corresponding estimators requires the generation of a large number of simulation runs, which is computationally expensive. To reduce the number of necessary runs, we propose a variance reduction technique based on control variates. We exploit constraints on the statistical moments of the stochastic process to reduce the estimators' variances. We develop an algorithm that selects appropriate control variates in an on-line fashion and demonstrate the efficiency of our approach on several case studies.
On the other end of the spectrum there are mean-field approximations, which model the mean densities faithfully in the system size limit @cite_15 . In between there are techniques such as moment closure @cite_14 , that not only consider the mean, but also the variance and other higher order moments. These methods depend on ad-hoc approximations of higher order moments to close the ODE system given by the moment equations. Yet another class of methods approximate molecular counts continuously and approximate the dynamics in such a continuous space, e.g. the system size expansion @cite_0 and the chemical Langevin equation @cite_26 .
@cite_0 @cite_15 @cite_14 @cite_26
2517280042 2150055571 2786311729 2051813288
Stochastic fluctuations of molecule numbers are ubiquitous in biological systems. Important examples include gene expression and enzymatic processes in living cells. Such systems are typically modelled as chemical reaction networks whose dynamics are governed by the Chemical Master Equation. Despite its simple structure, no analytic solutions to the Chemical Master Equation are known for most systems. Moreover, stochastic simulations are computationally expensive, making systematic analysis and statistical inference a challenging task. Consequently, significant effort has been spent in recent decades on the development of efficient approximation and inference methods. This article gives an introduction to basic modelling concepts as well as an overview of state of the art methods. First, we motivate and introduce deterministic and stochastic methods for modelling chemical networks, and give an overview of simulation and exact solution methods. Next, we discuss several approximation methods, including the chemical Langevin equation, the system size expansion, moment closure approximations, time-scale separation approximations and hybrid methods. We discuss their various properties and review recent advances and remaining challenges for these methods. We present a comparison of several of these methods by means of a numerical case study and highlight some of their respective advantages and disadvantages. Finally, we discuss the problem of inference from experimental data in the Bayesian framework and review recent methods developed the literature. In summary, this review gives a self-contained introduction to modelling, approximations and inference methods for stochastic chemical kinetics. In the stochastic formulation of chemical reactions, the dynamics of the first M-order moments of the species populations generally do not form a closed system of differential equations, in the sense that the time-derivatives of first M-order moments generally depend on moments of order higher than M. However, for analysis purposes, these dynamics are often made to be closed by approximating the needed derivatives of the first M-order moments by nonlinear functions of the same moments. These functions are called the moment closure functions. Recent results have introduced the technique of derivative-matching, where the moment closure functions are obtained by first assuming that they exhibit a certain separable form, and then matching time derivatives of the exact (not closed) moment equations with that of the approximate (closed) equations for some initial time and set of initial conditions. However, for multi-species reactions these results have been restricted to second order truncations, i.e, M = 2. This paper extends these results by providing explicit formulas to construct moment closure functions for any arbitrary order of truncation M. We show that with increasing M the closed moment equations provide more accurate approximations to the exact moment equations. Striking features of these moment closure functions are that they are independent of the reaction parameters (reaction rates and stoichiometry) and moreover the dependence of higher-order moment on lower order ones is consistent with the population being jointly lognormally distributed. To illustrate the applicability of our results we consider a simple bi-molecular reaction. Moment estimates from a third order truncation are compared with estimates obtained from a large number of Monte Carlo simulations Applying the method of moments to the chemical master equation appearing in stochastic chemical kinetics often leads to the so-called closure problem. Recently, several authors showed that this problem can be partially overcome using moment-based semidefinite programs (SDPs). In particular, they showed that moment-based SDPs can be used to calculate rigorous bounds on various descriptions of the stochastic chemical kinetic system’s stationary distribution(s)—for example, mean molecular counts, variances in these counts, and so on. In this paper, we show that these ideas can be extended to the corresponding dynamic problem, calculating time-varying bounds on the same descriptions. In order to solve the master equation by a systematic approximation method, an expansion in powers of some parameter is needed. The appropriate parameter is the reciprocal size of the system, defined as the ratio of intensive and extensive variables. The lowest approximation yields the phenomenological law for the approach to equilibrium. The next approximation determines the mean square of the fluctuations about the phenomenological behavior. In equilibrium this approximation has the form of a linear Fokker–Planck equation. The higher approximations describe the effect of the non-linearity on the fluctuations, in particular on their spectral density. The method is applied to three examples: density fluctuations, Alkemade's diode, and Rayleigh's piston. The relation to the expansion recently given by Siegel is also discussed.