^{1,a),b)}, Paul Kirk

^{1,a)}and Michael P. H. Stumpf

^{1,b)}

### Abstract

Moment approximation methods are gaining increasing attention for their use in the approximation of the stochastic kinetics of chemical reaction systems. In this paper we derive a general moment expansion method for any type of propensities and which allows expansion up to any number of moments. For some chemical reaction systems, more than two moments are necessary to describe the dynamic properties of the system, which the linear noise approximation is unable to provide. Moreover, also for systems for which the mean does not have a strong dependence on higher order moments, moment approximation methods give information about higher order moments of the underlying probability distribution. We demonstrate the method using a dimerisation reaction, Michaelis-Menten kinetics and a model of an oscillating p53 system. We show that for the dimerisation reaction and Michaelis-Menten enzyme kinetics system higher order moments have limited influence on the estimation of the mean, while for the p53 system, the solution for the mean can require several moments to converge to the average obtained from many stochastic simulations. We also find that agreement between lower order moments does not guarantee that higher moments will agree. Compared to stochastic simulations, our approach is numerically highly efficient at capturing the behaviour of stochastic systems in terms of the average and higher moments, and we provide expressions for the computational cost for different system sizes and orders of approximation. We show how the moment expansion method can be employed to efficiently quantify parameter sensitivity. Finally we investigate the effects of using too few moments on parameter estimation, and provide guidance on how to estimate if the distribution can be accurately approximated using only a few moments.

A.A. and M.P.H.S. acknowledge support from the Biotechnology and Biological Sciences Research Council, UK (Grant No. BB/G020434/1). P.K and M.P.H.S. acknowledge support from a Human Frontier Science Project (HFSP) grant (Grant No. RGP0061/2011). M.P.H.S. is a Royal Society Wolfson Research Merit Award holder.

I. INTRODUCTION

II. MOMENT EXPANSION METHOD

III. RESULTS

A. Dimerisation

B. Michaelis-Menten enzyme kinetics

C. P53 system

D. Parameter sensitivity estimation

E. Simple heuristics for moment expansions

F. Computational complexity

IV. CONCLUSION

### Key Topics

- Stochastic processes
- 15.0
- Enzyme kinetics
- 10.0
- Probability theory
- 10.0
- Chemical reactions
- 9.0
- Reaction kinetics modeling
- 6.0

##### C12

## Figures

Study of dimerization system, initial values x = [301, 0], parameters k = [1.66 × 10−3, 0.2]. (a) Single SSA realisation. (b) Average of 100 000 SSA simulations. (c) and (d) Histogram of SSA runs (grey bars) and probability density of normal distribution (blue line) calculated from mean and variance of SSA runs corresponding to points c and d in (b). (e) Mean for both variables, calculated using SSA, the moment approximation using 1 moment (deterministic) and two central moments (2m). (f) Variance of x 1 calculated using SSA, two central moments (2m), three central moments with Gaussian closure (3m), and four central moments (4m). (g) Third central moment calculated using SSA and moment approximation method, fourth central moment calculated using SSA, and moment approximation method.

Study of dimerization system, initial values x = [301, 0], parameters k = [1.66 × 10−3, 0.2]. (a) Single SSA realisation. (b) Average of 100 000 SSA simulations. (c) and (d) Histogram of SSA runs (grey bars) and probability density of normal distribution (blue line) calculated from mean and variance of SSA runs corresponding to points c and d in (b). (e) Mean for both variables, calculated using SSA, the moment approximation using 1 moment (deterministic) and two central moments (2m). (f) Variance of x 1 calculated using SSA, two central moments (2m), three central moments with Gaussian closure (3m), and four central moments (4m). (g) Third central moment calculated using SSA and moment approximation method, fourth central moment calculated using SSA, and moment approximation method.

Study of Michaelis-Menten kinetics, with parameters , and initial conditions, S(0) = 301, P(0) = 0, and E(0) = 120. (a) Single SSA realisation. Trajectories calculated using (b) moment approximation including only the mean (deterministic). (c) and (d) Variance of S and covariance between S and P calculated using SSA and approximation using 2 moments. (e) Skewness of S calculated using SSA and approximation with 3 moments. (f) Kurtosis calculated using SSA and approximation up to 4 and 6 moments.

Study of Michaelis-Menten kinetics, with parameters , and initial conditions, S(0) = 301, P(0) = 0, and E(0) = 120. (a) Single SSA realisation. Trajectories calculated using (b) moment approximation including only the mean (deterministic). (c) and (d) Variance of S and covariance between S and P calculated using SSA and approximation using 2 moments. (e) Skewness of S calculated using SSA and approximation with 3 moments. (f) Kurtosis calculated using SSA and approximation up to 4 and 6 moments.

Study of p53 model, parameter set q 1 = [90, 0.002, 1.7, 1.1, 0.93, 0.96, 0.01]. (a) Single SSA realisation. (b) Average of 100 000 SSA simulations. Trajectories calculated using (c) moment approximation including only the mean (deterministic), (d) linear noise approximation, (e) mean and the variance, (f) up to three central moments, (g) five central moments, and (h) six central moments. (i) and (j) Cumulative difference between mean trajectory calculated with SSA and trajectories calculated using (i) 2 moments and (j) 6 moments.

Study of p53 model, parameter set q 1 = [90, 0.002, 1.7, 1.1, 0.93, 0.96, 0.01]. (a) Single SSA realisation. (b) Average of 100 000 SSA simulations. Trajectories calculated using (c) moment approximation including only the mean (deterministic), (d) linear noise approximation, (e) mean and the variance, (f) up to three central moments, (g) five central moments, and (h) six central moments. (i) and (j) Cumulative difference between mean trajectory calculated with SSA and trajectories calculated using (i) 2 moments and (j) 6 moments.

Analysis of distribution for p53 model. (a) Variance calculated based on SSA runs. (b) Variance calculated with SSA, LNA, and moment approximation method. (c) Skewness calculated based on SSA runs (blue line) and skewness for normal distribution (cyan dashed line). (d) and (f) Histograms calculated based on SSA for points d, e, and f in Figure 4(c) , and probability density function of normal distribution calculated using mean and variance based on SSA (cyan line).

Analysis of distribution for p53 model. (a) Variance calculated based on SSA runs. (b) Variance calculated with SSA, LNA, and moment approximation method. (c) Skewness calculated based on SSA runs (blue line) and skewness for normal distribution (cyan dashed line). (d) and (f) Histograms calculated based on SSA for points d, e, and f in Figure 4(c) , and probability density function of normal distribution calculated using mean and variance based on SSA (cyan line).

Contour plots for the p53 system of the distance d between SSA trajectories and trajectories calculated with the moment expansion method. (a) and (b) Contour for varying k 1 and k 3, calculated using expansion up to 2 and 6 moments, respectively. (c) and (d) Contour for varying k 5 and k 6, calculated using expansion up to 2 and 6 moments, respectively. Red dots indicate parameter values used for SSA simulations.

Contour plots for the p53 system of the distance d between SSA trajectories and trajectories calculated with the moment expansion method. (a) and (b) Contour for varying k 1 and k 3, calculated using expansion up to 2 and 6 moments, respectively. (c) and (d) Contour for varying k 5 and k 6, calculated using expansion up to 2 and 6 moments, respectively. Red dots indicate parameter values used for SSA simulations.

Assessing the sensitivity to parameter k 1 of the dimerisation system of Sec. III A . Initial conditions and parameters are the same as previously, except that we additionally consider a perturbed value of k 1 = 1.1 × 1.66 × 10−3. (a) Average of 100 000 SSA simulations for the original and perturbed k 1 values. (b) As (a), except the mean is estimated using the proposed moment expansion approximation with 6 moments. (c) The sensitivity coefficients ∂μ1(t, [k 1, k 2])/∂k 1 and ∂μ2(t, [k 1, k 2])/∂k 2 (where μ1 and μ2 represent the means of x 1 and x 2, respectively) estimated from the SSA and 6m output via finite different approximations. (d), (e), and (f) As (a), (b), and (c), but for the second central moment. (g), (h), and (i) As (a), (b), and (c), but for the third central moment.

Assessing the sensitivity to parameter k 1 of the dimerisation system of Sec. III A . Initial conditions and parameters are the same as previously, except that we additionally consider a perturbed value of k 1 = 1.1 × 1.66 × 10−3. (a) Average of 100 000 SSA simulations for the original and perturbed k 1 values. (b) As (a), except the mean is estimated using the proposed moment expansion approximation with 6 moments. (c) The sensitivity coefficients ∂μ1(t, [k 1, k 2])/∂k 1 and ∂μ2(t, [k 1, k 2])/∂k 2 (where μ1 and μ2 represent the means of x 1 and x 2, respectively) estimated from the SSA and 6m output via finite different approximations. (d), (e), and (f) As (a), (b), and (c), but for the second central moment. (g), (h), and (i) As (a), (b), and (c), but for the third central moment.

Study of the deviation from the mean (x − μ). (a) Deterministic mean and single SSA trajectory for dimerisation system. (b) Deterministic mean and 9 points taken from different SSA trajectories for dimerisation system. (c) Deterministic mean and single SSA trajectory for p53 system. (d) Deterministic mean and 9 points taken from different SSA trajectories for p53 system.

Study of the deviation from the mean (x − μ). (a) Deterministic mean and single SSA trajectory for dimerisation system. (b) Deterministic mean and 9 points taken from different SSA trajectories for dimerisation system. (c) Deterministic mean and single SSA trajectory for p53 system. (d) Deterministic mean and 9 points taken from different SSA trajectories for p53 system.

## Tables

Error between mean, second, and third central moment calculated with SSA and approximation methods for the dimerization system.

Error between mean, second, and third central moment calculated with SSA and approximation methods for the dimerization system.

Number of potentially nonzero central moment terms to include in moment approximation for different numbers of variables (columns) and number of central moments (rows).

Number of potentially nonzero central moment terms to include in moment approximation for different numbers of variables (columns) and number of central moments (rows).

Article metrics loading...

Full text loading...

Commenting has been disabled for this content