gillespie {smfsb}

Simulate a sample path from a stochastic kinetic model described by a stochastic Petri net

Description

This function simulates a single realisation from a discrete stochastic kinetic model described by a stochastic Petri net (SPN).

Usage

gillespie(N, n, ...)

Arguments

N

An R list with named components representing a stochastic Petri net (SPN). Should contain N$M, a vector representing the initial marking of the net, N$Pre, a matrix representing the LHS stoichiometries, N$Post, a matrix representing the RHS stoichiometries, and N$h, a function representing the rates of the reaction processes. N$h should have first argument x, a vector representing the current state of the system, and second argument t, a scalar representing the current simulation time (in the typical time-homogeneous case, N$h will ignore this argument). N$h may posses additional arguments, representing reaction rates, for example.

n

An integer representing the number of events to simulate, excluding the initial state, N$M.

...

Additional arguments (such as reactions rates) will be passed into the function N$h.

Value

A list with first component t, a vector of length n containing event times and second component x, a matrix with n+1 rows containing the state of the system. The ith row of x contains the state of the system prior to the ith event.

See Also

simpleEuler, rdiff, discretise, gillespied, StepGillespie

Examples

# load the LV model
data(spnModels)
# simulate a realisation of the process and plot it
out = gillespie(LV,10000)
op = par(mfrow=c(2,2))
plot(stepfun(out$t,out$x[,1]),pch="")
plot(stepfun(out$t,out$x[,2]),pch="")
plot(out$x,type="l")

# use the "discretise" function to map it to an R "ts" object
plot(discretise(out,dt=0.01),plot.type="single",lty=c(1,2))
par(op)