0

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.

  • 1
    Thanks, I made a mistake while copying and pasting. This is the full code. – Claudio Laudani Jul 23 '22 at 00:38
  • 1
    `w[i,1:k-1]` looks dangerous to me. Did you mean, perhaps, `w[i,1:(k-1)]`? Your code is syntactically correct. We know that because it runs. The problem is with your data. Or the link between your data and your model. You’ve not given us sufficient context to investigate either of those. – Limey Jul 23 '22 at 06:26
  • @Limey The data are those written down to create the matrixes (actually this is a shorter version with only 8 rows per matrix, but it doesn't work either way); how can I give you more context? I'll try with the suggestion you gave me – Claudio Laudani Jul 23 '22 at 07:38
  • @Limey IT WORKED! But now it give me the error "Error in jags.model(model.file, data = data, inits = init.values, n.chains = n.chains, : RUNTIME ERROR: Incorrect number of parameters in function rank" I removed the ranking part of the code and works perfectly, do you know why is that? – Claudio Laudani Jul 23 '22 at 07:46

1 Answers1

1

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!