Bug in mcmcsamp
On 8/12/07, Jarrod Hadfield <j.hadfield at ed.ac.uk> wrote:
Hi,
I'm using lme4 0.99875-6 and I think there may be a bug in mcmcsamp. When I fit 2 random effects, the posterior variances always seem to be approximately the same irrespective of the true variances (or the REML estimates of the true variances). For example:
A<-gl(100,10)
B<-gl(100,1,1000)
y<-rnorm(100,0,sqrt(1))[A]+rnorm(100,0,sqrt(3))[B]+rnorm
(1000,0,sqrt(3))
model1<-lmer(y~(1|A)+(1|B))
model1MCMC<-mcmcsamp(model1,n=5000,trans=FALSE)
plot(model1MCMC)
Thanks for reporting the problem, Jarrod. I tried your script with the development version of lme4 and got the expected results, which I enclose Could you try the mcmcsamp.R script in lme4 0.99875-6 and tell me if the results are similar? -------------- next part -------------- A non-text attachment was scrubbed... Name: Hadfield.pdf Type: application/pdf Size: 338640 bytes Desc: not available URL: <https://stat.ethz.ch/pipermail/r-sig-mixed-models/attachments/20070812/dd518120/attachment.pdf> -------------- next part -------------- ## bug report by Jarrod Hadfield (2007-08-12) library(lme4) library(coda) set.seed(1) A<-gl(100,10) B<-gl(100,1,1000) y<-rnorm(100,0,sqrt(1))[A]+rnorm(100,0,sqrt(3))[B]+rnorm(1000,0,sqrt(3)) model1<-lmer(y~(1|A)+(1|B)) model1MCMC<-mcmcsamp(model1,n=5000,trans=FALSE) head(as.data.frame(model1MCMC), n = 20) HPDinterval(model1MCMC) pdf(file = "Hadfield.pdf", height = 8, width = 6) xyplot(model1MCMC, strip = FALSE, strip.left = TRUE, scales = list(x=list(axs='i')), type = c("g", "l")) densityplot(model1MCMC, plot.points = FALSE) dev.off()