User:Carl Boettiger/Notebook/Stochastic Population Dynamics/2010/05/21

From OpenWetWare
Jump to navigationJump to search
Owwnotebook icon.png Stochastic Population Dynamics Report.pngMain project page
Resultset previous.pngPrevious entry      Next entryResultset next.png

Puzzled over the two-step jump process results

The variance dynamics are slightly different if the continuous time Markov process allows events to contain two steps (such as a transition of an individual from one age class to another, or one patch to another, etc) as well as one step birth/death) transitions. [math]\displaystyle{ V\tau_{ij}(X/V) }[/math] to represent the rate of transitions from state i to state j. We change variables for our state space from the discrete count state variable X by dividing by the system size (volume) V, to obtain the continuous, mean-field approximation of the density p = X/V. Then:

[math]\displaystyle{ \begin{align} & \frac{d}{dt} E(p) = b(p)-d(p) + \sum_{j \neq i} \tau_{ji}(p) - \tau_{ij}(p) =: f(p)\\ & \frac{d}{dt} M = \textrm{diag}( b(p) + d(p) ) - T + Df(p) M + M Df(p)^T \end{align} }[/math]

where T is a symmetric matrix with off-diagonal elements i,j [math]\displaystyle{ \tau_{ij}(p) + \tau{ji}(p) }[/math] and diagonal elements chosen so that rows sum to zero. This approximation is implemented in noise_approx.R

Example: Crowley model

(Note this is the original Crowley model where the c coefficients are equal. In the original model this represented competition for space by overgrowth and hence would more naturally be expressed as a two step process).

[math]\displaystyle{ \begin{align} \dot x &= b_1 x (K - x - y)/K - d_1 x + c x y/K \\ \dot y &= b_2 y (K - x - y)/K - d_1 y - c x y/K \end{align} }[/math]

Interpreting the cxy terms as transitions rather than independent birth and death events results in the same mean-field dynamics but differing fluctuations. For instance, the covariance dynamics are:

[math]\displaystyle{ \frac{d}{dt} \textrm{Cov}(x,y) = \left(b_1 (1-2x-y)-d_1 +c y + b_2(1-x-2y)-d_2 -c x \right)\textrm{Cov}(x,y) - (b_1 x +c x )\sigma_y^2 -(b_2 y -c y)\sigma_x^2 - cxy }[/math]

Whereas interpreted as single-step birth and death the cxy term at the end is absent.

The variance dynamics are similarly altered.

Numerical implementation

This individual-based simulation is implemented in metapop.c and callable from R interface from ind_based_models.R as metapop_ibm(). The file metapop.R specifies the dynamics of birth, death, transition matrix and the Jacobian so that linear_noise_approx() function from above can be used to calculate the moments. This file also runs the crolwey_ibm model (from crowley.c, interface in the ind_based.models.R) and the linear noise approximation to it from crowley.R to compare how the two-step process changes the noise dynamics.

This comparison can be run simply by loading the package and running metapop_example() function.

For some reason the individual based-simulation of the two-step process seems to resemble that of the one-step process more, while the two-step linear noise approximation alternately underestimates (Fig 1) or overestimates (Fig 2) the variance. Looking for bugs in the implementation at the moment.

Fig 1. Parameters: b_1 = .2, b_2 = .1, d_1 = d_2 = c = .1
Fig. 2. Parameters: b_1 = .2, b_2 = .1, d_1 = d_2 = .1, c = .15