# [R] Gelman-Rubin convergence diagnostics via coda package

fbielejec fbielejec at gmail.com
Thu Mar 17 18:25:14 CET 2011

```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,

update(jags, 1000)

out <- coda.samples(model = jags,
variable.names =  c('beta'),
n.iter = 25000,
thin = 5)

gelman.plot(out)

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

```