A multilevel adaptive reaction-splitting simulation method for stochastic reaction networks

TL;DR

This paper introduces a multilevel adaptive simulation method for stochastic reaction networks that efficiently estimates system observables by adaptively partitioning reactions into high and low activity groups, reducing computational cost.

Contribution

The paper presents a novel multilevel Monte Carlo approach with adaptive reaction partitioning and a new variance reduction technique for stochastic reaction network simulation.

Findings

Achieves error bounds with complexity O(TOL^{-2})

Substantially outperforms standard SSA and previous hybrid methods

Uses a control variate based on Kurtz's stochastic time change representation

Abstract

Stochastic modeling of reaction networks is a framework used to describe the time evolution of many natural and artificial systems, including, biochemical reactive systems at the molecular level, viral kinetics, the spread of epidemic diseases, and wireless communication networks, among many other examples. In this work, we present a novel multilevel Monte Carlo method for kinetic simulation of stochastic reaction networks that is specifically designed for systems in which the set of reaction channels can be adaptively partitioned into two subsets characterized by either "high" or "low" activity. Adaptive in this context means that the partition evolves in time according to the states visited by the stochastic paths of the system. To estimate expected values of observables of the system at a prescribed final time, our method bounds the global computational error to be below a prescribed tolerance, $TOL$, within a given confidence level. This is achieved with a computational complexity of order $O(TOL^{-2})$, the same as with an exact method, but with a smaller constant. We also present a novel control variate technique based on the stochastic time change representation by Kurtz, which may dramatically reduce the variance of the coarsest level at a negligible computational cost. Our numerical examples show substantial gains with respect to the standard Stochastic Simulation Algorithm (SSA) by Gillespie and also our previous hybrid Chernoff tau-leap method.