Accelerated Algorithms for Stochastic Simulation of Chemically Reacting Systems
Stochastic models are widely used in the simulation of biochemical systems at a cellular level. For well mixed models, the system state can be represented by the population of each species. The probabilities for the system to be in each state are governed by the Chemical Master Equation (CME), which is generally a huge ordinary differential equation (ODE) system. The cost of solving the CME directly is generally prohibitive, due to its huge size.
The Stochastic Simulation Algorithm (SSA) provides a kinetic Monte Carlo approach to obtain the solution to the CME. It does this by simulating every reaction event in the system. A great many stochastic realizations must be performed, to obtain accurate probabilities for the states. The SSA can generate a highly accurate result, however the computation of many SSA realizations may be expensive if there are many reaction events. Tau-leaping is an approximate algorithm that can speed up the simulation for many systems. It advances the system with a selected stepsize. In each step, it directly samples the number of reaction events in each reaction channel, which yields a faster simulation than SSA. The error in tau-leaping is controlled by selecting the stepsize properly.
We have developed a new, accelerated tau-leaping algorithm for discrete stochastic simulation that make use of the fact that exact (time-dependent) solutions are known for some of the most common reaction motifs (subgraphs of the network of chemical species and reactants). This idea can be extended to spatial stochastic simulation, by treating the diffusion network as a special motif for which there is an exact time dependent solution. We describe the well-mixed and spatial stochastic time dependent solution algorithms, along with numerical experiments illustrating their effectiveness.