Tau-leaping


In probability theory, tau-leaping, or τ-leaping, is an approximate method for the simulation of a stochastic system. It is based on the Gillespie algorithm, performing all reactions for an interval of length tau before updating the propensity functions. By updating the rates less often this sometimes allows for more efficient simulation and thus the consideration of larger systems.
Many variants of the basic algorithm have been considered.

Algorithm

The algorithm is analogous to the Euler method for deterministic systems, but instead of making a fixed change
the change is
where is a Poisson distributed random variable with mean.
Given a state with events occurring at rate and with state change vectors, the method is as follows:
  1. Initialise the model with initial conditions.
  2. Calculate the event rates.
  3. Choose a time step. This may be fixed, or by some algorithm dependent on the various event rates.
  4. For each event generate, which is the number of times each event occurs during the time interval.
  5. Update the state by
  6. :
  7. :where is the change on state variable due to event. At this point it may be necessary to check that no populations have reached unrealistic values.
  8. Repeat from Step 2 onwards until some desired condition is met.

Algorithm for efficient step size selection

This algorithm is described by Cao et al. The idea is to bound the relative change in each event rate by a specified tolerance . This is achieved by bounding the relative change in each state variable by, where depends on the rate that changes the most for a given change in. Typically is equal the highest order event rate, but this may be more complex in different situations.
This algorithm typically requires computing auxiliary values, and should only require reusing previously calculated values. An important factor in this is that since is an integer value, there is a minimum value by which it can change, preventing the relative change in being bounded by 0, which would result in also tending to 0.
  1. For each state variable, calculate the auxiliary values
  2. :
  3. :
  4. For each state variable, determine the highest order event in which it is involved, and obtain
  5. Calculate time step as
  6. :
This computed is then used in Step 3 of the leaping algorithm.