我看了一个类似的问题，但它对我没有帮助。这是我的JAGS代码（包含在文本文件 simulazione_vivek.txt 中）：
model {
###############################################
######### Analisi per dati ISTAT ##########
###############################################
for (i in 1:ntot) {
bmi[i]~dnorm(mubmi[regione[i]],sdbmi[regione[i]])
fumo[i]~dbin(pfumo[regione[i]],1)
grassi[i]~dbin(pgrassi[regione[i]],1)
reddito[i]~dweib(vreddito[regione[i]],lambdareddito[regione[i]])
}
for (i in 1:popita) {
xreddito[i]~dweib(vreddito[reg[i]],lambdareddito[reg[i]])
xbmi[i]~dnorm(mubmi[reg[i]],sdbmi[reg[i]])
xfumo[i]~dbin(pfumo[reg[i]],1)
xgrassi[i]~dbin(pgrassi[reg[i]],1)
casi[i]~dbin(pcasi[i],1)
logit(pcasi[i])<-alpha+betar*xreddito[i]+betaf*xfumo[i]+betab*xbmi[i]+betag*xgrassi[i]
}
for (i in 1:nreg)
{
mubmi[i]~dnorm(0,0.001)
sdbmi[i]~dgamma(0.001,0.001)
vreddito[i]~dgamma(0.001,0.001)
lambdareddito[i]~dgamma(0.001,0.001)
pfumo[i]~dunif(0,1)
pfisica[i]~dunif(0,1)
pgrassi[i]~dunif(0,1)
}
alpha~dnorm(0,.001)
betar~dnorm(0,.001)
betaf~dnorm(0,.001)
betab~dnorm(0,.001)
betafi~dnorm(0,.001)
betag~dnorm(0,.001)
}
这是我的R代码：
library(rjags)
# simulation data from a survey with exposure but NOT response
nreg=10 #number of regions
personreg=50 #number of person per region
mureddito=rnorm(nreg,30000,2000) #simulation of values for income for each region
mubmi=rnorm(nreg,22,1) #simulation of values for bmi for each region
pfumo=runif(nreg,.3,.5)
pgrassi=runif(nreg,.1,.3)
reddito=rnorm(personreg,mureddito[1],100)
regione=rep(1,personreg)
bmi=rnorm(personreg,mubmi[1],.5)
fumo=rbinom(personreg,1,pfumo[1])
grassi=rbinom(personreg,1,pgrassi[1])
for (i in 2:nreg)
{
reddito=c(reddito,rnorm(personreg,mureddito[i],100))
regione=c(regione,rep(i,personreg))
bmi=c(bmi,rnorm(personreg,mubmi[i],.5))
fumo=c(fumo,rbinom(personreg,1,pfumo[i]))
grassi=c(grassi,rbinom(personreg,1,pfumo[i]))
}
ntot=length(reddito)
# simulation aggregated data
# breast #number of breast cancer for each region
# pop population
breast<-floor(runif(nreg,20,40))
pop<-floor(runif(nreg,1000,2000))
# data manipulation for the program
casi<-c(rep(1,breast[1]),rep(0,pop[1]-breast[1]))
reg<-c(rep(1,pop[1]))
for (i in 2:nreg)
{
casi<-c(casi,rep(1,breast[i]),rep(0,pop[i]-breast[i]))
reg<-c(reg,rep(i,pop[i]))
}
popita=length(casi)
nupdate=10
niter=10
input=list("ntot"=ntot,"regione"=regione,"nreg"=nreg,
"reddito"=reddito, "bmi"=bmi,"grassi"=grassi,
"fumo"=fumo,'casi'=casi,"casireg"=reg)
init=list("mubmi"=rep(0,nreg),"sdbmi"=rep(0.01,nreg),
"pfumo"=rep(0.5,nreg),"pgrassi"=rep(0.5,nreg),
"pfisica"=rep(0.5,nreg),
"vreddito"=rep(1,nreg),"lambdareddito"=rep(1,nreg))
m <- jags.model("simulazione_vivek.txt", data=input, init)
update(m, nupdate)
但我收到此错误消息：
Error in jags.model("simulazione_vivek.txt", data = input, init) :
RUNTIME ERROR:
Compilation error on line 16.
Cannot evaluate upper index of counter i
In addition: Warning message:
In jags.model("simulazione_vivek.txt", data = input, init) :
Unused variable "casireg" in data
> update(m, nupdate)
Error in update(m, nupdate) : object 'm' not found
请帮忙。