  Stochasticity plays a fundamental role in various biochemical processes, such
as cell regulatory networks and enzyme cascades. Isothermal, well-mixed systems
can be modelled as Markov processes, typically simulated using the Gillespie
Stochastic Simulation Algorithm (SSA). While easy to implement and exact, the
computational cost of using the Gillespie SSA to simulate such systems can
become prohibitive as the frequency of reaction events increases. This has
motivated numerous coarse-grained schemes, where the "fast" reactions are
approximated either using Langevin dynamics or deterministically. While such
approaches provide a good approximation when all reactants are abundant, the
approximation breaks down when one or more species exist only in small
concentrations and the fluctuations arising from the discrete nature of the
reactions becomes significant. This is particularly problematic when using such
methods to compute statistics of extinction times for chemical species, as well
as simulating non-equilibrium systems such as cell-cycle models in which a
single species can cycle between abundance and scarcity. In this paper, a
hybrid jump-diffusion model for simulating well- mixed stochastic kinetics is
derived. It acts as a bridge between the Gillespie SSA and the chemical
Langevin equation. For low reactant reactions the underlying behaviour is
purely discrete, while purely diffusive when the concentrations of all species
is large, with the two different behaviours coexisting in the intermediate
region. A bound on the weak error in the classical large volume scaling limit
is obtained, and three different numerical discretizations of the
jump-diffusion model are described. The benefits of such a formalism are
illustrated using computational examples.
