rjagsr2jags

Index out of range taking subset error with JAGS while performing network metanalysis


I was trying to develop a model to run a network metanalysis with jags and tried to use a previous code I've already used with WinBUGS (and worked properly).

However, when I try to use it with JAGS it gives me a coding error, but it still works perfectly fine with WINBugs.

This is the code:

   model{                            
for(i in 1:ns){                     
  w[i,1] <- 0    
  delta[i,1] <- 0           
  mu[i] <- log(p[i,1])
  p[i,1] ~ dunif(0,1)     
  for (k in 1:na[i]) {             
    r[i,k] ~ dbin(p[i,k],n[i,k]) 
   }
  for (k in 2:na[i]) {             
    log(p[i,k]) <- mu[i] + min(delta[i,k], -log(p[i,1]))
    delta[i,k] ~ dnorm(md[i,k],taud[i,k])
    md[i,k] <-  d[t[i,k]] - d[t[i,1]] + sw[i,k]
    taud[i,k] <- tau *2*(k-1)/k
    w[i,k] <- (delta[i,k] - d[t[i,k]] + d[t[i,1]])
    sw[i,k] <- sum(w[i,1:k-1])/(k-1)
   }
 }   

d[1]<-0       
for (k in 2:nt){  d[k] ~ dnorm(0,.1) }
sd ~ dunif(0,2)    
tau <- pow(sd,-2)  

for (c in 1:(nt-1)) {  
  for (k in (c+1):nt)  { 
      lnRR[c,k]<- d[c] - d[k]
      RR[c,k]<- exp(d[c]-d[k])
      lnRR[k,c]<- d[k] - d[c]
      RR[k,c]<- exp(d[k]-d[c])
   }  
 }

for (k in 1:nt) { 
  rk[k] <- rank(d[],k)  
  best[k] <- equals(rk[k],1)    

  for (h in 1:nt){  prob[h,k] <- equals(rk[k],h)  }
 }
}            

                   

And here are the imput I give:

   col_1<- c(1007 , 105, 1912,  467,  577 , 709 , 142,  100)
    
    col_2 <-c(1012 , 109, 1913 , 460 , 586 , 714,  153 , 103)
    
    n<-matrix(data=c(col_1, col_2),      # matrix of patients in each arm per trial
              nrow=length(col_1), ncol=2)
    
    col_1 <- c(79 , 0  ,9 , 7, 10, 31,  0,  5)
    
    col_2 <- c( 95 , 0 ,13, 17,  9 ,31 , 0 , 3)
    
    r<-matrix(data=c(col_1, col_2),      # matrix of events in each arm per trial
              nrow=length(col_1), ncol=2)
    
    col_1 <- c( 12,  2, 13,  4,  4,  5,  6 , 8)
    
    col_2 <- c( 95 , 0 ,13, 17,  9 ,31 , 0 , 3)
    
    t<-matrix(data=c(col_1, col_2),       # matrix of treatments (go from 1 to 14 in the database)
              nrow=length(col_1), ncol=2)

na<-rep(2, dim(n))

ns<-dim(n)[1]

nt<-14

info <- list ("n", "r", "t", "na", "ns", "nt")

parameters <- c("RR","lnRR", "rk",  "best")

library(R2jags)

jags(data= info, 
     parameters.to.save = parameters, 
     model.file = model.file, 
     n.chains = 3, 
     n.thin = nt, 
     n.iter = 3000,
     n.burnin = 500)

And this is the error it gives to me:

Error in jags.model(model.file, data = data, inits = init.values, n.chains = n.chains,  : 
  RUNTIME ERROR:
Compilation error on line 16.
Index out of range taking subset of  w

I think there are some differences with the WINBugs code that are still not clear to me;

Thanks in advance for the help.


Solution

  • Answer here:

    Thanks to the suggestions of Limey I corrected the code and then, after some research, I also corrected the ranking code.

    Here is the final version of the JAGS code:

    model{                            
    for(i in 1:ns){                     
      w[i,1] <- 0    
      delta[i,1] <- 0           
      mu[i] <- log(p[i,1])
      p[i,1] ~ dunif(0,1)     
      for (k in 1:na[i]) {             
        r[i,k] ~ dbin(p[i,k],n[i,k]) 
       }
      for (k in 2:na[i]) {             
        log(p[i,k]) <- mu[i] + min(delta[i,k], -log(p[i,1]))
        delta[i,k] ~ dnorm(md[i,k],taud[i,k])
        md[i,k] <-  d[t[i,k]] - d[t[i,1]] + sw[i,k]
        taud[i,k] <- tau *2*(k-1)/k
        w[i,k] <- (delta[i,k] - d[t[i,k]] + d[t[i,1]])
        sw[i,k] <- sum(w[i,1:(k-1)])/(k-1)
       }
     }   
    
    d[1]<-0       
    for (k in 2:nt){  d[k] ~ dnorm(0,.1) }
    sd ~ dunif(0,2)    
    tau <- pow(sd,-2)  
    
    for (c in 1:(nt-1)) {  
      for (k in (c+1):nt)  { 
          lnRR[c,k]<- d[c] - d[k]
          RR[c,k]<- exp(d[c]-d[k])
          lnRR[k,c]<- d[k] - d[c]
          RR[k,c]<- exp(d[k]-d[c])
       }  
     }
    
    rk_bad <- rank(d[])  
        rk_good <- rank(-d[])
        for(k in 1:nt){
            rk[k] <- rk_good[k]*equals(good_event, 1) + rk_bad[k] * equals(good_event, 0)
        }
    
        for(k in 1:nt){
            best[k] <- equals(rk[k], 1)
        }
    }   
    

    Since I'm doing a metanalysis of "bad events" (death or disease) I gave the information of good_event<-0

    Thanks all for the advices!