ssa.exact {adaptivetau}  R Documentation 
Simulates the trajectory of a continuoustime Markov process (aka. the Gillespie simulation algorithm).
ssa.exact(init.values, transitions, rateFunc, params, tf)
init.values 
Vector of initial values for all variables. This should be a simple onedimensional vector of real numbers. IMPORTANT: variables must never take negative values. 
transitions 
Twodimensional matrix of integers specifying how each state variable (rows) should be changed for a given transition (columns). Generally this will be a sparse matrix of primarily 1s and 1s. 
rateFunc 
Rfunction that returns instantaneous transition rates for each
transition in the form a realvalued onedimensional vector with
length equal to the number of transitions. The order of these
rates must match the order in which transitions were specified in
the

params 
any R variable to be passed asis to rateFunc, presumably specifying useful parameters. 
tf 
final time of simulation. Due to the approximate nature of the tauleaping algorithm, the simulation may overshoot this time somewhat. 
This function is supplied as a bonus with the adaptivetau package, since the C++ function that underlies this (exact) stochastic simulation algorithm is used in the (approximate) adaptive tau stochastic simulation as well.
The initial values, transition matrix, and transition rates must all be designed such that variables are always nonnegative. The algorithm relies on this invariant.
Consider calling enableJIT(1)
before running
ssa.exact. In most circumstances, this should yield some speedup by
bytecode compiling the rate function.
Twodimensional matrix with rows for every timepoint at which the rateFunc was evaluated and columns giving the value for every state variable at that time. The first column specifies the time.
Philip Johnson
This function is a bonus the comes along with the approximate (but
faster) ssa.adaptivetau
. For systems with sparse
transition matrices, see helper function ssa.maketrans
.
## LotkaVolterra example lvrates < function(x, params, t) { with(params, { return(c(preygrowth*x["prey"], ## prey growth rate x["prey"]*x["pred"]*eat, ## prey death / predator growth rate x["pred"]*preddeath)) ## predator death rate }) } params=list(preygrowth=10, eat=0.01, preddeath=10); r=ssa.exact(c(prey = 1000, pred = 500), matrix(c(1,0, 2,1, 0,1), nrow=2), lvrates, params, tf=2) matplot(r[,"time"], r[,c("prey","pred")], type='l', xlab='Time', ylab='Counts') legend("topleft", legend=c("prey", "predator"), lty=1:2, col=1:2)