Incorporating extrinsic noise into the stochastic simulation of biochemical reactions: A comparison of approaches

The chemical Langevin equation for biochemical systems in dynamic environments

Modelling and simulation of complex biochemical reaction networks form cornerstones of modern biophysics. Many of the approaches developed so far capture temporal fluctuations due to the inherent stochasticity of the biophysical processes, referred to as intrinsic noise. Stochastic fluctuations, however, predominantly stem from the interplay of the network with many other - and mostly unknown - fluctuating processes, as well as with various random signals arising from the extracellular world; these sources contribute extrinsic noise. Here we provide a computational simulation method to probe the stochastic dynamics of biochemical systems subject to both intrinsic and extrinsic noise. We develop an extrinsic chemical Langevin equation-a physically motivated extension of the chemical Langevin equation- to model intrinsically noisy reaction networks embedded in a stochastically fluctuating environment. The extrinsic CLE is a continuous approximation to the Chemical Master Equation (CME) with time-varying propensities. In our approach, noise is incorporated at the level of the CME, and can account for the full dynamics of the exogenous noise process, irrespective of timescales and their mismatches. We show that our method accurately captures the first two moments of the stationary probability density when compared with exact stochastic simulation methods, while reducing the computational runtime by several orders of magnitude. Our approach provides a method that is practical, computationally efficient and physically accurate to study systems that are simultaneously subject to a variety of noise sources.

A mathematical modelling framework for the regulation of intra-cellular OCT4 in human pluripotent stem cells

Human pluripotent stem cells (hPSCs) have the potential to differentiate into all cell types, a property known as pluripotency. A deeper understanding of how pluripotency is regulated is required to assist in controlling pluripotency and differentiation trajectories experimentally. Mathematical modelling provides a non-invasive tool through which to explore, characterise and replicate the regulation of pluripotency and the consequences on cell fate. Here we use experimental data of the expression of the pluripotency transcription factor OCT4 in a growing hPSC colony to develop and evaluate mathematical models for temporal pluripotency regulation. We consider fractional Brownian motion and the stochastic logistic equation and explore the effects of both additive and multiplicative noise. We illustrate the use of time-dependent carrying capacities and the introduction of Allee effects to the stochastic logistic equation to describe cell differentiation. We conclude both methods adequately capture the decline in OCT4 upon differentiation, but the Allee effect model has the advantage of allowing differentiation to occur stochastically in a sub-set of cells. This mathematical framework for describing intra-cellular OCT4 regulation can be extended to other transcription factors and developed into predictive models.

Stochastic simulation algorithms for computational systems biology: Exact, approximate, and hybrid methods

Nowadays, mathematical modeling is playing a key role in many different research fields. In the context of system biology, mathematical models and their associated computer simulations constitute essential tools of investigation. Among the others, they provide a way to systematically analyze systems perturbations, develop hypotheses to guide the design of new experimental tests, and ultimately assess the suitability of specific molecules as novel therapeutic targets. To these purposes, stochastic simulation algorithms (SSAs) have been introduced for numerically simulating the time evolution of a well-stirred chemically reacting system by taking proper account of the randomness inherent in such a system. In this work, we review the main SSAs that have been introduced in the context of exact, approximate, and hybrid stochastic simulation. Specifically, we will introduce the direct method (DM), the first reaction method (FRM), the next reaction method (NRM) and the rejection-based SSA (RSSA) in the area of exact stochastic simulation. We will then present the τ-leaping method and the chemical Langevin method in the area of approximate stochastic simulation and an implementation of the hybrid RSSA (HRSSA) in the context of hybrid stochastic-deterministic simulation. Finally, we will consider the model of the sphingolipid metabolism to provide an example of application of SSA to computational system biology by exemplifying how different simulation strategies may unveil different insights into the investigated biological phenomenon. This article is categorized under: Models of Systems Properties and Processes > Mechanistic Models Analytical and Computational Methods > Computational Methods.

Uniformization techniques for stochastic simulation of chemical reaction networks

This work considers the method of uniformization for continuous-time Markov chains in the context of chemical reaction networks. Previous work in the literature has shown that uniformization can be beneficial in the context of time-inhomogeneous models, such as chemical reaction networks incorporating extrinsic noise. This paper lays focus on the understanding of uniformization from the viewpoint of sample paths of chemical reaction networks. In particular, an efficient pathwise stochastic simulation algorithm for time-homogeneous models is presented which is complexity-wise equal to Gillespie's direct method. This new approach therefore enlarges the class of problems for which the uniformization approach forms a computationally attractive choice. Furthermore, as a new application of the uniformization method, we provide a novel variance reduction method for (raw) moment estimators of chemical reaction networks based upon the combination of stratification and uniformization.

Efficient finite-difference method for computing sensitivities of biochemical reactions

Sensitivity analysis of biochemical reactions aims at quantifying the dependence of the reaction dynamics on the reaction rates. The computation of the parameter sensitivities, however, poses many computational challenges when taking stochastic noise into account. This paper proposes a new finite-difference method for efficiently computing sensitivities of biochemical reactions. We employ propensity bounds of reactions to couple the simulation of the nominal and perturbed processes. The exactness of the simulation is preserved by applying the rejection-based mechanism. For each simulation step, the nominal and perturbed processes under our coupling strategy are synchronized and often jump together, increasing their positive correlation and hence reducing the variance of the estimator. The distinctive feature of our approach in comparison with existing coupling approaches is that it only needs to maintain a single data structure storing propensity bounds of reactions during the simulation of the nominal and perturbed processes. Our approach allows to compute sensitivities of many reaction rates simultaneously. Moreover, the data structure does not require to be updated frequently, hence improving the computational cost. This feature is especially useful when applied to large reaction networks. We benchmark our method on biological reaction models to prove its applicability and efficiency.

A Critical Comparison of Rejection-Based Algorithms for Simulation of Large Biochemical Reaction Networks

The rejection-based simulation technique has been applying to improve the computational efficiency of the stochastic simulation algorithm (SSA) in simulating large reaction networks, which are required for a thorough understanding of biological systems. We compare two recently proposed simulation methods, namely the composition-rejection algorithm (SSA-CR) and the rejection-based SSA (RSSA), aiming for this purpose. We discuss the right interpretation of the rejection-based technique used in these algorithms in order to make an informed choice when dealing with different aspects of biochemical networks. We provide the theoretical analysis as well as the detailed runtime comparison of these algorithms on concrete biological models. We highlight important factors that are omitted in previous analysis of these algorithms. The numerical comparison shows that for reaction networks where the search cost is expensive then SSA-CR is more efficient, and for reaction networks where the update cost is dominant, often the case in practice, then RSSA should be the choice.