R2jags运行报索引越界错误且无法查看模型结果如何解决？

错误原因 #

• 观测变量名大小写不匹配：JAGS模型里的似然部分用的是Y[i,j]，但你传给jags.data的观测矩阵命名是小写的y，JAGS找不到匹配的Y变量，索引时就会触发越界错误。
• 数据读取格式错误：读取Thoropa.txt时如果没有跳过表头行，会把第一行的PCR1~PCR12字符串也读成数据，导致矩阵混入非数值，或行列数计算错误。
• 初始值函数参数错误：你写的runif(0,1)第一个参数是生成随机数的个数，填0会生成空值，正确写法为runif(1, 0, 1)，初始值格式错误也会导致运行报错。
• 多余的先验声明：你在R环境里写的psi ~ dunif(0, 1)这类代码是JAGS语法，不能直接在R里运行，直接删除即可。
• 初始脚本找不到结果的原因：你把jags返回的模型对象命名为Thoropa_model，但提取结果时用的是model$summary，变量名不匹配，替换为Thoropa_model$summary即可读取结果。

修复后可运行代码 #

# 加载依赖包
library(R2jags)

# 第一步：读取检测矩阵，跳过表头，转换为数值矩阵
Thoropa <- read.table("Thoropa.txt", header = TRUE, sep = "\t")
y <- as.matrix(Thoropa)

# 第二步：写入模型文件
sink("Thoropa2.txt")
cat("model {
    # Priors
    psi ~ dunif(0,1)
    p11 ~ dunif(0,1)
    p10 ~ dunif(0,p10_max)
    
    # Likelihood 
    for (i in 1:S){
    z[i] ~ dbern(psi)
    p[i] <- z[i]*p11 + (1-z[i])*p10
    for (j in 1:K){
    Y[i,j] ~ dbern(p[i])
    }
    }
    } ",fill=TRUE)
sink()

# 第三步：准备模型输入参数，变量名和模型内保持一致
S <- nrow(y)
K <- ncol(y)
p10_max <- 0.05
jags.data <- list(Y = y, S = S, K = K, p10_max = p10_max)

# 第四步：修正初始值函数，新增z的初始值提升收敛速度
jags.inits <- function(){
  list(
    psi = runif(1, 0.05, 0.95),
    p11 = runif(1, p10_max, 1),
    p10 = runif(1, 0, p10_max),
    z = ifelse(rowSums(y) > 0, 1, 0)
  )
}

jags.params <- c("psi", "p11", "p10")

# 运行模型
Thoropa_model <- jags.parallel(
  data = jags.data, 
  inits = jags.inits, 
  parameters.to.save = jags.params, 
  model.file = "Thoropa2.txt", 
  n.chains = 4, 
  n.thin = 10, 
  n.iter = 100000, 
  n.burnin = 50000, 
  jags.seed = 333
)

模型结果查看方法 #

运行完成后使用以下命令查看结果：

# 打印完整汇总结果，包含均值、分位数、Rhat收敛诊断值（Rhat<1.01说明收敛良好）
print(Thoropa_model)
# 单独查看参数的汇总矩阵
Thoropa_model$summary
# 查看特定参数的中位数
Thoropa_model$summary["psi", "50%"]
Thoropa_model$summary["p11", "50%"]
Thoropa_model$summary["p10", "50%"]
# 绘制参数后验分布和链轨迹图，检查收敛情况
plot(Thoropa_model)

内容的提问来源于stack exchange，提问作者Paula Cabral Eterovick