The regime-conversion method: a hybrid technique for simulating well-mixed chemical reaction networks

Abstract

There exist several methods for simulating biological and physical systems as represented by chemical reaction networks. Systems with low numbers of particles are frequently modeled as discrete-state Markov jump processes and are typically simulated via a stochastic simulation algorithm (SSA). An SSA, while accurate, is often unsuitable for systems with large numbers of individuals, and can become prohibitively expensive with increasing reaction frequency. Large systems are often modeled deterministically using ordinary differential equations, sacrificing accuracy and stochasticity for computational efficiency and analytical tractability. In this paper, we present a novel hybrid technique for the accurate and efficient simulation of large chemical reaction networks. This technique, which we name the regime-conversion method, couples a discrete-state Markov jump process to a system of ordinary differential equations by simulating a reaction network using both techniques simultaneously. Individual molecules in the network are represented by exactly one regime at any given time, and may switch their governing regime depending on particle density. In this manner, we model high copy-number species using the cheaper continuum method and low copy-number species using the more expensive, discrete-state stochastic method to preserve the impact of stochastic fluctuations at low copy number. The motivation, as with similar methods, is to retain the advantages while mitigating the shortfalls of each method. We demonstrate the performance and accuracy of our method for several test problems that exhibit varying degrees of inter-connectivity and complexity by comparing averaged trajectories obtained from both our method and from exact stochastic simulation.