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(); }