Youyi Fong
2012-05-21
Hi, I am wondering what is the variance of MNMetropolis proposal distribution
in JAGS. To compare, MCMCglmm uses the posterior samples up to any time point
to guide the proposal (page 117 of the MCMCglmm course note). In my
experience, it seems to do a better job than JAGS for the following Poisson
GLMM model:
model {
for (i in 1:n) {
dat ~ dpois(exp(b[dat,1] + b[dat,2]dat + b[dat,3]dat + b[dat,4]dat +
b[dat,5]dat) )
}
for (i in 1:m) {
b ~ dmnorm(theta, Tau) # b does not work
}
theta ~ dmnorm (mu0, inverse(Lambda0))
Tau ~ dwish(R0, 7)
}
Thanks,
Youyi
Youyi Fong
2012-05-22
In a flurry of activities, I might have compared orange with apple. In a
properly-run, very preliminary test, JAGS and MCMCglmm appear to perform
similarly for this example.
Martyn Plummer
2012-05-23
Well that's good, but you should reallly see what happens when you load the
glm module. You should get the auxiliary mixture variables sampler AuxMix.
Youyi Fong
2012-05-23
I loaded glm module and didn't see much difference Here is how I did it.
list.modules()
"basemod" "bugs" "glm"
jmodel = jags.model(file= "jags_models/tumor_glmm.txt", data=jags.data,
inits=jags.inits)
Compiling model graph
Resolving undeclared variables
Allocating nodes
Graph Size: 5230
Initializing model
|++++++++++++++++++++++++++++++++++++++++++++++++++| 100%
list.samplers(jmodel)
$MNormMetropolis
"b"
Is AuxMix a sampler that is supposed to get used after glm is loaded? I cannot
find more information about it in the jags manual.
Youyi Fong
2012-05-24
After centering each covariate, the mixing behavior greatly improves. By the
way, not sure why, but in the first post, the 'random effect line of the model
didn't come out right. Here it is again
model {
for (i in 1:n) {
dat ~ dpois(exp(b[dat,1] + b[dat,2:5] %*% dat) )
}
for (i in 1:m) {
b ~ dmnorm(theta, Tau) # random effects
}
theta ~ dmnorm (mu0, inverse(Lambda0))
Tau ~ dwish(R0, 7)
}
Youyi Fong
2012-05-24
The model still didn't come out right. Here it is in quotes. Hopefully that
will fix it.
model {
for (i in 1:n) {
dat ~ dpois(exp(b[dat,1] + b[dat,2:5] %*% dat) )
}
for (i in 1:m) {
b ~ dmnorm(theta, Tau)
}
theta ~ dmnorm (mu0, inverse(Lambda0))
Tau ~ dwish(R0, 7)
}
Youyi Fong
2012-05-24
Still not right, some indexes are missing in the likelihood like and random
effect line. Try it one more time using the code tag.
model { for (i in 1:n) { dat[i,2] ~ dpois(exp(b[dat[i,1],1] + b[dat[i,1],2:5] %*% dat[i,3:6]) ) } for (i in 1:m) { b[i,1:5] ~ dmnorm(theta, Tau) } theta ~ dmnorm (mu0, inverse(Lambda0)) Tau ~ dwish(R0, 7) }