[R] Gelman-Rubin convergence diagnostics via coda package

[fbielejec]

Dear,

I'm trying to run diagnostics on MCMC analysis (fitting a log-linear
model to rates data). I'm getting an error message when trying
Gelman-Rubin shrink factor plot:

>gelman.plot(out)
Error in chol.default(W) : 
  the leading minor of order 2 is not positive definite

I take it that somewhere, somehow a matrix is singular, but how can
that be remedied? 

My code:

library(rjags)

#------------JAGS MODEL3 CODE-------------#
data {
   for (i in 1:N) {
     y[i] <- z2[i] - z1[i]
   }
}

model {

  for (i in 1 : N) {
  
      mean[i] <- (beta[i] * (n[i] - z1[i]) * (n[i] - z2[i]) ) / n[i]
            
      y[i] ~ dpois( mean[i] )
      log(beta[i]) <-  a1 * type[i] + a0
}

a1 ~ dnorm(0, 1.0E-4)
a0 ~ dnorm(0, 1.0E-4)

}
#----------------------------------# 

N = 12
n = c(63, 33, 18, 16, 36, 3, 33, 101, 64, 30, 10, 36)
z1 = c(33, 8, 14, 4, 2, 8, 24, 76, 35, 14, 6, 22)
z2 = c(37, 8, 16, 6, 3, 10, 28, 91, 56, 17, 8, 28)
type = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1)

init0 = list(a0 = 0.0, a1 = 0.0, .RNG.seed = 1234, .RNG.name =
"base::Super-Duper") init1 = list(a0 = 1.0, a1 = 1.0, .RNG.seed =
1234, .RNG.name = "base::Super-Duper")

jags <- jags.model('model3.bug',
                   data = list('N' = N,
                               'n' = n,
                               'z1' = z1,
                               'z2' = z2,
                               'type' = type),
                   inits = list(init0, init1),            
                   n.chains = 2,
                   n.adapt = 100)
      
update(jags, 1000)
      
out <- coda.samples(model = jags,
                    variable.names =  c('beta'),
                    n.iter = 25000,
                    thin = 5)

gelman.plot(out)

-- 
while(!succeed) { try(); }