Model Simulation - Part I
2. A simple approximation of the distribution of Xt
From Unit 3 we know that P(Xt+1 = x) can be obtained as a sum of all the possible paths X1 = xk,1 , … , Xt = xk,t, Xt+1 = x, where xk,j is a vector of 0s and 1s. Let assume that the initial metapopulation state is certain, that is π assigns the unit mass to a given vector (see Remark 3.2, Unit 3). Then, what the simulation function “sim1” does is just to provide one of such a paths, the matrix “sim”, thanks to the transition probabilities “E” and “C” (the j-th column of “sim” is xk,j ). Therefore, it can be intuitively understood that by repeating this sampling procedure many times to obtain several “likely” paths and then summing up over all the simulations we can obtain an approximation of P(Xt+1 = x).
Example 7.3. Basic IFM, whitout rescue effect.
As far as the basic IFM is concerned, simulations are obtained by applying the function “sim1” of Example 7.1, where “C” is computed a part based on “pm”.
sim1<-function(o){
o<- o > 0
if(any(o)){
cond<-ifelse(o,E,C)
o<-ifelse(runif(length(o))<cond,!o,o)
}
as.numeric(o)
}
# To carry out ns ≥ 2 simulations we repeat the procedure of Example # 7.1 and we sum all together the ns matrices sim to obtain the matrix # simT. Then, the [i,j]-th cell of simT says how many of the ns # simulations assign 1 to the patch i. All the simulations start from # p1.
ns=10000
nt=10
simT<-matrix(0,nrow=length(p1),ncol=nt+1)
for(k in 1:ns) {
sim<-matrix(0,nrow=length(p1),ncol=nt+1)
sim[,1]<-p1
for(t in 1:nt) sim[ ,t+1]<-sim1(sim[ ,t])
simT=simT+sim
}
simT[1:13,]
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11]
[1,] 10000 2374 7763 4049 6561 4943 5977 5271 5654 5453 5646
[2,] 0 9996 3316 7811 4854 6745 5464 6380 5713 6253 5853
[3,] 10000 2429 8197 3806 7093 4638 6489 5083 6175 5363 5898
[4,] 0 9998 2860 7964 4416 6847 5206 6296 5550 6097 5675
[5,] 0 9984 3044 7848 4518 6775 5262 6303 5597 6041 5813
[6,] 0 10000 3414 7691 4840 6811 5497 6358 5852 6072 5911
[7,] 0 7682 4392 5753 5212 5404 5335 5353 5371 5307 5458
[8,] 10000 2879 7947 4244 6973 5004 6442 5385 6179 5631 5970
[9,] 10000 3249 7817 4760 6756 5512 6289 5725 6087 5896 6016
[10,] 10000 3467 7741 4984 6774 5599 6368 5939 6138 5996 6105
[11,] 10000 3261 7780 4713 6845 5389 6385 5688 6125 5864 6075
[12,] 10000 2845 7783 4387 6684 5048 6271 5444 5961 5585 5945
[13,] 0 9990 3371 7765 4833 6810 5417 6375 5783 6268 5839
simT[1:13,]/ns
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11]
[1,] 1 0.2374 0.7763 0.4049 0.6561 0.4943 0.5977 0.5271 0.5654 0.5453 0.5646
[2,] 0 0.9996 0.3316 0.7811 0.4854 0.6745 0.5464 0.6380 0.5713 0.6253 0.5853
[3,] 1 0.2429 0.8197 0.3806 0.7093 0.4638 0.6489 0.5083 0.6175 0.5363 0.5898
[4,] 0 0.9998 0.2860 0.7964 0.4416 0.6847 0.5206 0.6296 0.5550 0.6097 0.5675
[5,] 0 0.9984 0.3044 0.7848 0.4518 0.6775 0.5262 0.6303 0.5597 0.6041 0.5813
[6,] 0 1.0000 0.3414 0.7691 0.4840 0.6811 0.5497 0.6358 0.5852 0.6072 0.5911
[7,] 0 0.7682 0.4392 0.5753 0.5212 0.5404 0.5335 0.5353 0.5371 0.5307 0.5458
[8,] 1 0.2879 0.7947 0.4244 0.6973 0.5004 0.6442 0.5385 0.6179 0.5631 0.5970
[9,] 1 0.3249 0.7817 0.4760 0.6756 0.5512 0.6289 0.5725 0.6087 0.5896 0.6016
[10,] 1 0.3467 0.7741 0.4984 0.6774 0.5599 0.6368 0.5939 0.6138 0.5996 0.6105
[11,] 1 0.3261 0.7780 0.4713 0.6845 0.5389 0.6385 0.5688 0.6125 0.5864 0.6075
[12,] 1 0.2845 0.7783 0.4387 0.6684 0.5048 0.6271 0.5444 0.5961 0.5585 0.5945
[13,] 0 0.9990 0.3371 0.7765 0.4833 0.6810 0.5417 0.6375 0.5783 0.6268 0.5839