我正在尝试使用状态空间模型来估计人口统计数据(繁殖力、存活率、人口增长、人口规模)。我们有 4 个不同的年龄状态。
# J0 = number of individuals 0-1
# surv1 = survivorship from 0-1
# J1 = number of individuals 0-1
# surv2 = survivorship from 1-2
# J2= = number of individuals 0-1
# surv3 = survivorship from 2-3
# J3= number of individuals 0-1
# survad = survivorship >3 "adult")
# Data as vectors (Talek clan from 1988-2013)
# X0 = individuals 0-1 in years
# X1 = individuals 1-2 in years
# X2 = individuals 2-3 in years
# X3 = individuals 3+ in years
# Total = group size
X0 <- c(7, 9, 4, 8, 9, 5, 8, 5, 7, 5, 5, 8, 10, 3, 5, 7, 2, 6, 6, 11, 14, 12, 15, 9, 10)
X1 <- c( 4, 4, 3, 4, 8, 5, 2, 4, 3, 4, 4, 5, 3, 7, 0, 5, 6, 3, 3, 5, 10, 12, 10, 13, 8)
X2 <- c(3, 2, 3, 3, 3, 8, 4, 1, 1, 2, 2, 4, 2, 2, 5, 0, 5, 5, 4, 3, 3, 10, 12, 7, 10)
X3 <- c(18, 16, 13, 16, 29, 29, 26, 22, 21, 18, 16, 15, 16, 15, 11, 14, 9, 12, 16, 18, 21, 23, 33, 32, 31)
Total <- c(32, 31, 23, 31, 49, 47, 40, 32, 32, 29, 27, 32, 31, 27, 21, 26, 22, 26, 29, 37, 48, 57, 70, 61, 59)
这是错误代码:
sink(file = "HyenaIPM_all.txt")
cat("
model {
# Specify the priors for all parameters in the model
N.est[1] ~ dnorm(50, tau.proc)T(0,) # Initial abundance
mean.lambda ~ dunif(0, 5)
sigma.proc ~ dunif(0, 50)
tau.proc <- pow(sigma.proc, -2)
for (t in 1:TT) {
fec[t] ~ dunif(0, 5) # per capita fecundidty
surv1[t] ~ dunif(0, 1) # survivorship from 0-1
surv2[t] ~ dunif(0, 1) # survivorship from 1-2
surv3[t] ~ dunif(0, 1) # survivorship from 2-3
survad[t] ~ dunif(0, 1) # adult survivorship
}
# Estimate fecundity and survivorship
for (t in 2:TT) {
# Fecundity
J0[t+1] ~ dpois(survad[t]*fec[t])
J0[t+1] <- J3[t] * fec[t]
# Survivorship
J1[t+1] ~ dbin(surv1[t], J0[t])
J1[t+1] <- J0[t]*surv1[t]
J2[t+1] ~ dbin(surv2[t], J1[t])
J2[t+1] <- J1[t]*surv2[t]
J3[t+1] ~ dbin(surv3[t], J2[t-1])
J3[t+1] <- J2[t]*surv3[t] + J3[t]*survad[t]
A[t+1] ~ dbin(survad[t], A[t])
A[t+1] <- J3[t]*surv3[t] + A[t]*survad[t]
# Lambda
lambda[t+1] ~ dnorm(mean.lambda, tau.proc)
N.est[t+1] <- N.est[t]*lambda[t]
}
# Population size
for (t in 1:TT){
N[t] ~ dpois(N.est[t])
}
}
", fill = T)
sink()
# Parameters monitored
sp.params <- c("fec", "surv1", "surv2", "surv3", "survad", "lambda")
# MCMC settings
ni <- 200
nt <- 10
nb <- 100
nc <- 3
# Initial values
sp.inits <- function()list(mean.lambda = runif(1, 0, 1))
#Load all the data
sp.data <- list(N = Total, TT = length(Total), J0 = X0, J1 = X1, J2 = X2, J3 = X3)
library(R2jags)
hyena_model <- jags(sp.data, sp.inits, sp.params, "HyenaIPM_all.txt", n.chains = nc, n.thin = nt, n.iter = ni, n.burnin = nb)
不幸的是,当我运行代码时出现以下错误。
Error in jags.model(model.file, data = data, inits = init.values, n.chains = n.chains, :
RUNTIME ERROR:
Index out of range for node J0
有人对我们为什么会收到此错误有任何建议吗?不知道为什么 J0 的分布是错误的。