Jags RUNTIME ERROR: ошибка компиляции в строке 17. Попытка переопределить дельту узла [15,1,1] - PullRequest
Новая
       45

Jags RUNTIME ERROR: ошибка компиляции в строке 17. Попытка переопределить дельту узла [15,1,1]

0 голосов
/ 08 октября 2018

У меня проблема с компиляцией модели Jags с использованием следующего кода.Кажется, что параметр delta пытается быть переопределен (строка 17: delta [i, si [i, k], hh] ~ dnorm (md [i, si [i, k], hh], taud [i,si [i, k], hh])) в начале третьего цикла в модели, но я не уверен, почему?

Модель: 

        model
{
    for (hh in 1:nrepeat) {
        for (i in 1:ns) {
            w[i, 1, hh] <- 0.00000E+00
            j[i, 1, hh] <- 0.00000E+00
            delta[i, bi[i, hh], hh] <- 0.00000E+00
            mu[i, hh] ~ dnorm(0.00000E+00, 1.00000E-05)
            beta0[i, hh] ~ dnorm(0.00000E+00, 1.00000E-05)
            b1[i, 1, hh] <- 0.00000E+00
            for (k in 1:na[i]) {
                index[i, k, hh] <- split[i, hh] * (equals(t[i, 
                  k], pair[hh, 1]) + equals(t[i, k], pair[hh, 
                  2]))
            }
            for (k in 2:na[i]) {
                delta[i, si[i, k], hh] ~ dnorm(md[i, si[i, k], 
                  hh], taud[i, si[i, k], hh])
                md[i, si[i, k], hh] <- (d[si[i, k], hh] - d[bi[i, 
                  hh], hh] + sw[i, k, hh] + b1[i, k, hh]) * (1 - 
                  index[i, m[i, k], hh]) + direct[hh] * index[i, 
                  m[i, k], hh]
                j[i, k, hh] <- k - (equals(1, split[i, hh]) * 
                  step(k - 3))
                taud[i, si[i, k], hh] <- tau[hh] * 2 * (j[i, 
                  k, hh] - 1)/j[i, k, hh]
                w[i, k, hh] <- (delta[i, si[i, k], hh] - ((d[si[i, 
                  k], hh] - d[bi[i, hh], hh]) + b1[i, k, hh])) * 
                  (1 - index[i, k, hh])
                sw[i, k, hh] <- sum(w[i, 1:(k - 1), hh])/(j[i, 
                  k, hh] - 1)
                b1[i, k, hh] <- (beta[t[i, k], hh] * (x[i] - 
                  mx))
            }
        }
        for (l in 1:np) {
            r[l] ~ dbin(p[l], n[l])
            logit(p[l]) <- mu[s[l], hh] + beta0[s[l], hh] * (x[l] - 
                mx) + delta[s[l], tipd[l], hh] + (deltab[l, hh] * 
                (1 - equals(tipd[l], bi[s[l], hh])))
            rhat[l] <- p[l] * n[l]
            dev[l] <- 2 * (r[l] * (log(r[l]) - log(rhat[l])) + 
                (n[l] - r[l]) * (log(n[l] - r[l]) - log(n[l] - 
                  rhat[l])))
            index2[l, hh] <- split[s[l], hh] * (equals(tipd[l], 
                pair[hh, 1]) + equals(tipd[l], pair[hh, 2]))
            deltab[l, hh] <- (beta[tipd[l], hh] - beta[bi[s[l], 
                hh], hh]) * (x[l] - mx) * (1 - index2[l, hh]) + 
                directbeta[hh] * (x[l] - mx) * (index2[l, hh])
        }
        totresdev[hh] <- sum(dev[])
        direct[hh] ~ dnorm(0.00000E+00, 1.00000E-05)
        directbeta[hh] ~ dnorm(0.00000E+00, 1.00000E-05)
        d[1, hh] <- 0.00000E+00
        beta[1, hh] <- 0.00000E+00
        sd[hh] ~ dunif(0.00000E+00, 10)
        tau[hh] <- pow(sd[hh], -2)
        tausq[hh] <- sd[hh] * sd[hh]
        for (k in 2:nt) {
            d[k, hh] ~ dnorm(0.00000E+00, 1.00000E-05)
            beta[k, hh] ~ dnorm(0.00000E+00, 1.00000E-05)
        }
        for (k in 1:nt) {
            for (v in 1:nz) {
                dz[v, k, hh] <- d[k, hh] - (beta[k, hh]) * (mx - 
                  z[v])
            }
        }
        for (c in 1:(nt - 1)) {
            for (k in (c + 1):nt) {
                betas[c, k, hh] <- beta[k, hh] - beta[c, hh]
                lor[c, k, hh] <- (d[k, hh] - d[c, hh])
                for (v in 1:nz) {
                  lorz[v, c, k] <- (dz[v, k, hh] - dz[v, c, hh])
                }
            }
        }
        for (v in 1:nz) {
            directz[v, hh] <- direct[hh] - (directbeta[hh]) * 
                (mx - z[v])
            directorz[v, hh] <- exp(directz[v, hh])
        }
        for (v in 1:nz) {
            diff[v, hh] <- directz[v, hh] - lorz[v, pair[hh, 
                1], pair[hh, 2]]
            prob[v, hh] <- step(diff[v, hh])
        }
    }
}

Файл данных можетполучить здесь: файл данных и модели

Ввод для JAGS: 

rm(list = ls())
library(coda)
library(rjags)
library(netmeta);

#LOAD FUNCTIONS TO SHAPE DATA
#CHECK IF PAIR(X,Y) IN ROW I OF DATA AND GIVE BASELINE FOR DATA ROW I
PairXY <- function(treat, pair)
{
  N <- nrow(treat)
  out <- cbind(split=rep(0,N), b=rep(0,N))
  for (i in 1:N) {
    pos <- match(pair, treat[i,], nomatch=0) # lenght = length(pair) = 2
    out[i,1] <- ifelse(prod(pos)>0, 1, 0) # 1 if pair in line i, 0 o.w.
    out[i,2] <- ifelse(prod(pos)==0, 1, pos[1])
  }
  out
}
# GIVES NA-1 INDEXES TO SWEEP NON-BASELINE ARMS ONLY
NonbaseSweep <- function(index, na)
{
  N <- NROW(na)
  C <- max(na)
  out <- matrix(nrow=N, ncol=C)
  for (i in 1:N) {
    for (k in 2:na[i]) {
      out[i,k] <- k - (index[i,"b"] >= k)


    }
  }
  out
}
# BUILDS MATRIX WITH NON-BASELINE TREATMENTS
Sweeptreat <- function(treat, m)
{
  N <- NROW(treat)
  C <- NCOL(m)
  out <- matrix(nrow=N, ncol=C)
  for (i in 1:N) {
    for (k in 2:C) {
      out[i,k] <- treat[i,m[i,k]]
    }
  }
  out
}
## BUILDS VECTOR WITH BASELINE TREATMENTS
Basetreat <- function(treat, b)
{
  N <- nrow(treat)
  out <- rep(0,N)
  for (i in 1:N) {
    out[i] <- treat[i,b[i]]
  }
  out
}

# Data
datasrc <- read.csv('OR_Responders.csv', header=TRUE)
datasrc[] <- lapply(datasrc, function(x){
  x[is.nan(x)] <- NA
  x
}) 
datasrc<-subset(datasrc,datasrc$XSev1!="NA")

mx<-mean(c(datasrc$XSev1,datasrc$XSev2,datasrc$XSev3),na.rm = TRUE); 
dat1<-datasrc
dat1$study<-1:length(dat1$t1)
dat2 <- pairwise(list(t1, t2, t3),event=list(r1,r2,r3),n=list(n1, n2, n3),data=dat1, studlab=study)

na=dat1$na #NUMBER OF ARMS IN EACH STUDY
t=cbind(dat1$t1,dat1$t2,dat1$t3, deparse.level = 0) #TREATMENT NUMBER
s=dat2$studlab #STUDY NUMBER
#y=dat2$TE #OUTCOME
r=dat2$event1+dat2$event2
n=dat2$n1+dat2$n2
tipd=dat2$treat2 #TREATMENT NUMBER
x=dat2$XSevag #COVARIATE VALUES
#x=dat2$age/12 #COVARIATE VALUES
ns=max(s) #NMUBER OF TRIALS
nt=22 #NUMBER OF TREATMENTS
np=length(x) #NUMBER OF PATIENTS
mx=mean(x) #AVERAGE COVARIATE VALUE
z=c(10,20,30,40,50) #CHOSEN COVARIATE VALUES AT WHICH TREATMENT EFFECTS ARE REQUIRED TO BE ESTIMATED.
nz=length(z) #NUMBER OF CHOSEN COVARIATE VALUES


#CHOOSE NODE TO SPLIT
parameters = c("prob","diff");model.file1 = "ModelC1.txt"
contre <- read.csv('Treat.csv', header=TRUE)

pair<-contre[,2:3];pair<-data.matrix(pair, rownames.force = NA)
split <- array(0, dim=c(ns,length(contre[,1])))
bi <- array(0, dim=c(ns,length(contre[,1])))


for (i in 1:length(contre[,1])) {
  pairx <- c(contre[i,2],contre[i,3]);
  checkPair <- PairXY(t, pairx);split[,i]<-checkPair[,"split"]
  bix <- Basetreat(t, checkPair[,"b"]);bi[,i] = bix;
  m <- NonbaseSweep(checkPair, na);#m[,i] = mtt;
  si <- Sweeptreat(t,m);#si[,i] = six;
}

# Defining some MCMC parameters for JAGS
nchains  <- 3; # How Many Chains?
nadapt<-100
nburnin  <- 100; # How Many Burn-in Samples?
nsamples <- 100;  # How Many Recorded Samples?
nthin <- 10;

datastruct<-list(z=z, nz=nz,ns=ns,nt=nt,mx=mx,s=s,na=na,nrepeat=length(pair[,1]),split = split,m=m,bi=bi,si=si,pair=pair,x=x,tipd=tipd,t=t,r=r,n=n,np=np)
mod1 <- jags.model(file =model.file1, data=datastruct, n.chains=nchains, n.adapt=nadapt);
# update(mod1, nburnin);
# sim<-coda.samples(mod1, parameters, nsamples, thin = nthin, na.rm=TRUE);
# simsum<-summary(sim);

Сообщение об ошибке: 

Compiling model graph
   Resolving undeclared variables
   Allocating nodes
Deleting model

Error in jags.model(file = model.file1, data = datastruct, n.chains = nchains,  : 
  RUNTIME ERROR:
Compilation error on line 17.
Attempt to redefine node delta[15,1,1]

Любая подсказка будетбыть очень ценным.Спасибо Лисо

...