二项分布年度N值的JAGS估计：报错排查与修正

问题背景 #

我需要用贝叶斯方法从二项分布中估计N值，手头有40年的计数数据和估计的成功率，目标是得到每年的N值估计。作为JAGS新手，我尝试给n添加索引n[i]来获取年度估计，结果报错："Index out of range taking subset of n"。我已经能跑通生成单一N值的模型（移除似然中n的索引即可），但不确定当前思路是否正确。

初始错误模型与测试数据 #

sink("file.jags")
cat("
    model {

    ## Likelihood    
    for (i in 1:nyear) {
      
      x[i] ~ dbin(theta, n[i])
     
    }

    ## Priors    
    lambda ~ dgamma(0.005, 0.005)    
    theta ~ dbeta(22, 37)
    mu <- lambda/theta
    n ~ dpois(mu)
    
}    
", fill = TRUE)

sink()

# Initial values
inits <- function() { 
  list(
    n = sample(1:100
               , size = 1) 
    , lambda = runif(1, 1, 100)
  )
}

# Parameters monitored
params <- c('theta'
            , 'n'
)

# MCMC settings
ni <- 2000
nb <- 200
nc <- 3

# Bundle data
jags.data <- list(n = c(round(runif(40, 1, 100)))
                  , nyear = 40)

# Fit model
out <- jagsUI(data = jags.data
              , inits = inits
              , parameters.to.save = params
              , model.file = "file.jags"
              , n.chains = nc
              , n.iter = ni
              , n.burnin = nb
              , verbose = TRUE)

# Summarize posteriors
print(out$summary, dig = 2)

修正后的模型 #

根据建议，将n的先验放入循环中生成年度先验，同时从初始值列表中移除n，解决了"setParameters"维度不匹配的问题，模型成功运行。此外将lambda的先验从gamma分布改为均匀分布，测试后两种先验结果无显著差异，最终保留均匀分布先验。

"model {

    # Likelihood    
    for (i in 1:nyear) {
      
      reported_morts[i] ~ dbin(report_rate, total_morts[i])
      
      # Annual prior for total mortality. Necessary to get annual estimates.
      total_morts[i] ~ dpois(mu)
    }

    # Priors
    lambda ~ dunif(1, 100)
    report_rate ~ dbeta(22, 37)
    mu <- round(lambda/report_rate)
    
}"  

# Initial values. Removed "n = sample(1:100, size = 1)".
inits <- function() { 
  list(
    lambda = runif(1, 1, 100)
  )
}  

内容的提问来源于stack exchange，提问作者Matt