I created a population model that at points, was doing exactly what I wanted it to do. However, as I began to continue rerunning the code to yield continued correct results, to verify I had done the code correctly, the data would change in an incorrect way. Essentially, the data would space out and the altogether disappear into all zeros. This occurred without me making any changes to the code.
I tried commenting out some of my "if statement" parameters, changing the dataset from multiplying matrix "d1" and "df2" to "d1" and "pop_model", and then changed the plot from df2 to pop_model. I even moved up my df2 variable to before my arguments, thinking this was the problem.
projection_matrix <- matrix(c(0,4.4,8.8,0.2,0,0,0,.4,.4),nrow=3,ncol=3,byrow=T)
projection_matrix
disease_matrix_1 <- matrix(c(0,.22,.44,.05,0,0,0,.05,.05),nrow=3,ncol=3,byrow=T)
disease_matrix_1
d1 <- disease_matrix_1
disease_matrix_2 <- matrix(c(0,.66,1.32,.15,0,0,0,.15,.15),nrow=3,ncol=3,byrow=T)
disease_matrix_2
d2 <- disease_matrix_2
disease_matrix_3 <- matrix(c(0,.22,.44,.05,0,0,0,.05,.05),nrow=3,ncol=3,byrow=T)
disease_matrix_3
d3 <- disease_matrix_3
disease_matrix_4 <- matrix(c(0,.22,.44,.05,0,0,0,.05,.05), nrow=3, ncol=3, byrow=T)
disease_matrix_4
d4 <- disease_matrix_4
disease_matrix_5 <- matrix(c(0,.22,.44,.05,0,0,0,.05,.05), nrow=3, ncol=3, byrow=T)
disease_matrix_5
d5 <- disease_matrix_5
disease_matrix_6 <- matrix(c(0,.22,.44,.05,0,0,0,.05,.05), nrow=3, ncol=3, byrow=T)
disease_matrix_6
d6 <- disease_matrix_6
disease_matrix_7 <- matrix(c(0,.22,.44,.05,0,0,0,.05,.05), nrow=3, ncol=3, byrow=T)
disease_matrix_7
d7 <- disease_matrix_7
disease_matrix_8 <- matrix(c(0,.22,.44,.05,0,0,0,.05,.05), nrow=3, ncol=3, byrow=T)
disease_matrix_8
d8 <- disease_matrix_8
disease_matrix_9 <- matrix(c(0,.22,.44,.05,0,0,0,.05,.05), nrow=3, ncol=3, byrow=T)
disease_matrix_9
d9 <- disease_matrix_9
disease_matrix_10 <- matrix(c(0,.22,.44,.05,0,0,0,.05,.05), nrow=3, ncol=3, byrow=T)
disease_matrix_10
d10 <- disease_matrix_10
abundance_yr0 <- c(0,0,12)
abundance_yr0
Year1 <- projection_matrix %*% abundance_yr0 # matrix multiplication!
Year1
Year2 <- projection_matrix %*% Year1
Year2
t=200; N0=abundance_yr0; x=projection_matrix
-------------------------------------------------------
pop_model <- matrix(0,nrow=nrow(x),ncol=t+1)
pop_model[,1] <- N0
df2 <- replace(pop_model, pop_model > 2.00e+9, sample(2.00e+9:2.5e+9, 10000, replace=T))
for(i in 2:(t+1)){
if(i == 135){
pop_model[,i] <- d1 %*% df2[,i-1]
}
else{
pop_model[,i] <- x %*% df2[,i-1]}
if(i == 155){
pop_model[,i] <- d2 %*% df2[,i-1]
}
}
if(i == 160){
pop_model[,i] <- d3 %*% df2[,i-1]
}
if(i == 161){
pop_model[,i] <- d4 %*% df2[,i-1]
}
if(i == 162){
pop_model[,i] <- d5 %*% df2[,i-1]
}
if(i == 163){
pop_model[,i] <- d6 %*% df2[,i-1]
}
if(i == 164){
pop_model[,i] <- d7 %*% df2[,i-1]
}
if(i == 165){
pop_model[,i] <- d8 %*% df2[,i-1]
}
if(i == 166){
pop_model[,i] <- d9 %*% df2[,i-1]
}
if(i == 167){
pop_model[,i] <- d10 %*% df2[,i-1]
}
}
plot(1,1,pch="",ylim=c(0,3.00e+9),xlim=c(0,t+1),xlab="Years",ylab="Abundance",xaxt="n")
cols <- topo.colors(ncol(x))
max(pop_model)
for(s in 1:ncol(x)){
points(df2[s,],col=cols[s],type="l",lwd=2)
}
axis(1,at=seq(1,t+1),labels = seq(0,t))
legend("topleft",col=cols,lwd=rep(2,ncol(x)),legend=c("Infant", "Adolescent", "Adult"))
When run correctly it should produce a population model graph with a large population crash around 155-165 depending if you run the whole code or not.
When running something like df2[156] or df2[3,156] the numbers should be significantly lower than the rest of the graph, almost near 0.
