GillespieSSA test. Single species

GillespieSSA test. Single species

h qin

December 13, 2015

Example from Pineda-Krch08 J of Statistical Software.

rm(list=ls())

library(GillespieSSA)

parms = c(b=2, d=1, K=1000)

x0 = c(N=500)

nu = matrix(c(+1,-1), ncol=2)

a = c("b*N", "(b+(b-d)*N/K)*N")

tf = 10

method = "D"

simName = "Logistic growth"

parms

##    b    d    K 

##    2    1 1000

nu 

##      [,1] [,2]

## [1,]    1   -1

a

## [1] "b*N"             "(b+(b-d)*N/K)*N"

out = ssa(x0, a, nu, parms, tf, method, simName, verbose = TRUE, consoleInterval = 1)

## Running D method with console output every 1 time step

## Start wall time: 2015-12-13 11:55:55...

## t=0 : 500

## (0.204s) t=1.000299 : 370

## (0.343s) t=2.000827 : 256

## (0.439s) t=3.001757 : 214

## (0.539s) t=4.001413 : 194

## (0.624s) t=5.003206 : 164

## (0.694s) t=6.001344 : 128

## (0.762s) t=7.000501 : 138

## (0.825s) t=8.000797 : 120

## (0.884s) t=9.00031 : 107

## t=10.0004 : 83

## tf: 10.0004

## TerminationStatus: finalTime

## Duration: 0.95 seconds

## Method: D

## Nr of steps: 8381

## Mean step size: 0.001193222+/-0.001501032

## End wall time: 2015-12-13 11:55:56

## --------------------

str(out)

## List of 3

##  $ data : num [1:8383, 1:2] 0 0.000869 0.003032 0.003286 0.003868 ...

##   ..- attr(*, "dimnames")=List of 2

##   .. ..$ : chr [1:8383] "timeSeries" "" "" "" ...

##   .. ..$ : chr [1:2] "" "N"

##  $ stats:List of 8

##   ..$ startWallime      : chr "2015-12-13 11:55:55"

##   ..$ endWallTime       : chr "2015-12-13 11:55:56"

##   ..$ elapsedWallTime   : Named num 0.95

##   .. ..- attr(*, "names")= chr "elapsed"

##   ..$ terminationStatus : chr "finalTime"

##   ..$ nSteps            : int 8381

##   ..$ meanStepSize      : num 0.00119

##   ..$ sdStepSize        : num 0.0015

##   ..$ nSuspendedTauLeaps: num 0

##  $ args :List of 18

##   ..$ x0                 : Named num 500

##   .. ..- attr(*, "names")= chr "N"

##   ..$ a                  : chr [1:2] "b*N" "(b+(b-d)*N/K)*N"

##   ..$ nu                 : num [1, 1:2] 1 -1

##   ..$ parms              : Named num [1:3] 2 1 1000

##   .. ..- attr(*, "names")= chr [1:3] "b" "d" "K"

##   ..$ tf                 : num 10

##   ..$ method             : chr "D"

##   ..$ tau                : num 0.3

##   ..$ f                  : num 10

##   ..$ epsilon            : num 0.03

##   ..$ nc                 : num 10

##   ..$ hor                : num NaN

##   ..$ dtf                : num 10

##   ..$ nd                 : num 100

##   ..$ ignoreNegativeState: logi TRUE

##   ..$ consoleInterval    : num 1

##   ..$ censusInterval     : num 0

##   ..$ verbose            : logi TRUE

##   ..$ simName            : chr "Logistic growth"