Dear all,
I'm not 100% sure if this question is best directed at the r-list, or a mailing list concerned with Bayesian analysis, so please accept my apologies if another audience may be more appropriate.
I have been using the rjags package to run Jags models with multiple chains and store the results in a Coda based mcmc list. For instance, having created a jags model and done initial adapting and updating, I run the following command:
coda_odp_gini_only <- coda.samples(odp_gini_only,
variable.names=c("beta","sigma2.u2", "deviance"),
n.iter=itercount, thin=thincount)
This create an object with 4 separate chains as requested in my initial call to Jags which created the object "odp_gini_only". I then use the Coda package to look at the results stored in "coda_odp_gini_only". Prior to running any analysis with Coda, I used the command,
coda.options(combine.plots=TRUE, combine.stats=TRUE)
to ask for results that combine the four separate chains. Sure enough, if I enter "summary(coda_odp_gini_only)", I am given a single set of output combining the four chains. However, if I enter "HPDinterval(coda_odp_gini_only)" I receive 4 sets of HPD figures, one for each chain. Is it possible to combine the four chains together to receive a single set of HPD estimates?
In a similar vein, is it possible to use the Coda object to estimate the proportion of a given parameters distribution which is above (or below) a given value, for instance, the proportion of the distribution of beta[1] greater than zero? Again, in doing so, is it possible to combine the results of the four chains into a single estimate?
Kind regards, and many thanks in advance for any advice anyone can offer me.
Paul
Dr Paul Norris
Lecturer in Social Policy
University of Edinburgh
The University of Edinburgh is a charitable body, registered in Scotland, with registration number SC005336.