r - R2OpenBugs 错误消息：预期集合运算符 c 错误 pos 2174

Question

以下代码在 OpenBugs 中运行良好，但是当我在 R2OpenBugs 中运行它时，我收到一条错误消息“预期集合运算符 c 错误 pos 2174”。感谢您帮助修复错误。谢谢。

R2OpenBugsCode

# Load the R2OpenBUGS library
library(coda)
library(boot)
library(R2OpenBUGS)
library(mcmcplots)
# workingdirectory
setwd("C:/Users/)
sink("model3.txt")        
cat("# Binomial likelihood, logit link, subgroup
# Fixed effects model with one covariate
model{                          # *** PROGRAM STARTS
for(i in 1:ns){                 # LOOP THROUGH STUDIES
    mu[i] ~ dnorm(0,.0001)      # vague priors for all trial baselines
    for (k in 1:na[i])  {       # LOOP THROUGH ARMS
        r[i,k] ~ dbin(p[i,k],n[i,k])    # binomial likelihood
# model for linear predictor, covariate effect relative to treat in arm 1
        logit(p[i,k]) <- mu[i] + d[t[i,k]] - d[t[i,1]]
                        + (beta[t[i,k]]-beta[t[i,1]]) * x[i]
# expected value of the numerators 
        rhat[i,k] <- p[i,k] * n[i,k]
#Deviance contribution
        dev[i,k] <- 2 * (r[i,k] * (log(r[i,k])-log(rhat[i,k]))
             +  (n[i,k]-r[i,k]) * (log(n[i,k]-r[i,k]) - log(n[i,k]-rhat[i,k])))
      }
# summed residual deviance contribution for this trial
    resdev[i] <- sum(dev[i,1:na[i]])
     }   
totresdev <- sum(resdev[])      # Total Residual Deviance
d[1] <- 0    # treatment effect is zero for reference treatment
beta[1] <- 0 # covariate effect is zero for reference treatment
for (k in 2:nt){  # LOOP THROUGH TREATMENTS
    d[k] ~ dnorm(0,.0001) # vague priors for treatment effects
    beta[k] <- B # common covariate effect
  }
B ~ dnorm(0,.0001) # vague prior for covariate effect
# treatment effect when covariate = z[j]
for (k in 1:nt){  # LOOP THROUGH TREATMENTS
    for (j in 1:nz) { dz[j,k] <- d[k] + (beta[k]-beta[1])*z[j] }
  }
# pairwise ORs and LORs for all possible pair-wise comparisons, if nt>2
for (c in 1:(nt-1)) {  
    for (k in (c+1):nt)  { 
# when covariate is zero
        or[c,k] <- exp(d[k] - d[c])
        lor[c,k] <- (d[k]-d[c])
# at covariate=z[j]
        for (j in 1:nz) {
            orz[j,c,k] <- exp(dz[j,k] - dz[j,c])
            lorz[j,c,k] <- (dz[j,k]-dz[j,c])
          }
     }  
 }
# Provide estimates of treatment effects T[k] on the natural (probability) scale
# Given a Mean Effect, meanA, for 'standard' treatment A, 
# with precision (1/variance) precA and covariate value z[j]
#A ~ dnorm(meanA,precA)
#for (k in 1:nt) { 
#    for (j in 1:nz){
#        logit(T[j,k]) <- A + d[k] + (beta[k]-beta[1]) * z[j]  
#      }
#  }
}                                                     # *** PROGRAM ENDS
", fill = TRUE)
sink() 
cat(
"t[,1]  t[,2]   na[]    r[,1]   n[,1]   r[,2]   n[,2]   x[] #   ID  name
1   2   2   256 2223    182 2221    1   #   1   4S
1   2   2   4   125 1   129 1   #   2   Bestehorn
1   2   2   0   52  1   94  1   #   3   Brown
1   2   2   2   166 2   165 1   #   4   CCAIT
1   2   2   77  3301    80  3304    0   #   5   Downs
1   2   2   3   1663    33  6582    0   #   6   EXCEL
1   2   2   8   459 1   460 1   #   7   Furberg
1   2   2   3   155 3   145 1   #   8   Haskell
1   2   2   0   42  1   83  1   #   9   Jones
1   2   2   4   223 3   224 0   #   10  KAPS
1   2   2   633 4520    498 4512    1   #   12  LIPID
1   2   2   1   124 2   123 1   #   13  MARS
1   2   2   11  188 4   193 1   #   14  MAAS
1   2   2   5   78  4   79  1   #   15  PLAC 1
1   2   2   6   202 4   206 1   #   16  PLAC 2
1   2   2   3   532 0   530 0   #   17  PMSGCRP
1   2   2   4   178 2   187 1   #   18  Riegger
1   2   2   1   201 3   203 1   #   19  Weintraub
1   2   2   135 3293    106 3305    0   #   20  Wscotland
", file="mydata.txt")
mydata=read.table("mydata.txt", header = TRUE,col.names=c("t1","t2","na","r1","n1","r2","n2","x"),sep="")
t<-cbind(mydata$t1,mydata$t2)
t<- as.matrix(t)
na<-as.vector(mydata$na)
r<-cbind(mydata$r1,mydata$r2)
r<-as.matrix(r)
n<-cbind(mydata$n1,mydata$n2)
n<-as.matrix(n)
x<-as.vector(mydata$x)
ns=19
nt=2
z=c(1)
nz=1
N = nrow(mydata)
mydata3 = list("N","t","na","r","n","x", "ns", "nt", "nz", "z")
# Parameters to be monitored (= to estimate)
params = c( "d", "B")
inits <- list(
inits1=list(d=c( NA, 0), mu=c(0,0,0,0,0,   0,0,0,0,0,   0,0,0,0,0,   0,0,0,0), B=0),
inits2=list(d=c( NA, -1), mu=c(-3,-3,3,-3,3,    -3,3,-3,3,-3,    -3,-3,3,3,-3,   3,-3,-3,3), B=-1),
inits3=list(d=c( NA, 2), mu=c(-3,5,-1,-3,7,    -3,-4,-3,-3,0,    5,0,-2,-5,1,   -2,5,3,0), B=1.5))
# MCMC settings
nc <- 3      #number of MCMC chains to run
ni <- 50000  #number of samples for each chain     
nb <- 20000   #number of samples to discard as burn-in
nt <- 1      #thinning rate, increase this to reduce autocorrelation
output <- bugs(data=mydata3, inits=inits, parameters.to.save=params, model.file="model3.txt", 
               n.chains=nc, n.iter=ni, n.burnin=nb, n.thin=nt, debug=TRUE, DIC=TRUE,
               working.directory=getwd())

我得到的日志文件错误是这样的：

model is syntactically correct
**expected the collection operator c error pos 2174**
variable ns is not defined
model must have been compiled but not updated to be able to change RN generator
BugsCmds:NoCompileInits
BugsCmds:NoCompileInits
BugsCmds:NoCompileInits
model must be compiled before generating initial values
model must be initialized before updating
model must be initialized before monitors used
model must be initialized before monitors used
model must be initialized before monitors used
model must be initialized before DIC can be monitored
model must be initialized before updating
model must be initialized before monitors used
DIC monitor not set

标签： ropenbugs