Gillespie algorithm (source code)

= Gillespie algorithm
{c}
{wiki}

= Stochastic simulation algorithm
{c}
{synonym}

The <Gillespie algorithm> samples statistically exact trajectories of a well-mixed <chemical reaction network> modeled by a <continuous-time Markov chain>. At copy-number state $x$, each <reaction propensity function> $a_r(x)$ is a nonnegative event rate and each <stoichiometric vector> $\nu_r$ specifies the update. Assume the rates have no explicit time dependence. Put $a_0(x)=\sum_r a_r(x)$. If $a_0(x)=0$, no further <chemical reaction> can occur: the state is an <absorbing state>.

For $a_0(x)>0$, sample a waiting time $\tau$ from the <exponential distribution> with rate $a_0(x)$, independently choose channel $r$ with <probability> $a_r(x)/a_0(x)$, and update $(t,x)\leftarrow(t+\tau,x+\nu_r)$. Recompute the <reaction propensity functions> and repeat. The no-event <probability> is $S(\tau)=e^{-a_0\tau}$ because $S'=-a_0S$; multiplying by $a_r$ gives the joint waiting-time/channel <probability density function> $a_re^{-a_0\tau}$. This factorization explains the two draws.<Approximate accelerated stochastic simulation of chemically reacting systems by Daniel T. Gillespie>[Section II, equations (5)–(9), review of the exact stochastic simulation algorithm]{ref}

For $A\to B$ with <reaction propensity function> $a(x)=kx_A$, the waiting time has mean $1/(kx_A)$, and $x_A=0$ is an <absorbing state>. Exactness concerns the specified <chemical master equation>; spatial mixing errors and explicitly time-varying rates require a different model or sampling rule.