lmer model and mcmcsamp
Hi,
I have had similar problems - I came to the conclusion that a) I had
to make some transformation such as log something/divide by sigma etc
b) the "locally-uniform" priors are i) not uniform or ii) uniform but
in the wrong locality or c) mcmcsamp is broken.
It would be good to get some clarity; mcmcsamp even seems to behave
oddly for the example data set:
fm1 <- lmer(Reaction ~ Days + (1|Subject), sleepstudy)
fm1.mcmc<-mcmcsamp(fm1, n=1000)
summary(fm1)
Linear mixed model fit by REML
Formula: Reaction ~ Days + (1 | Subject)
Data: sleepstudy
AIC BIC logLik deviance REMLdev
1794 1807 -893.2 1794 1786
Random effects:
Groups Name Variance Std.Dev.
Subject (Intercept) 1378.18 37.124
Residual 960.46 30.991
Number of obs: 180, groups: Subject, 18
Fixed effects:
Estimate Std. Error t value
(Intercept) 251.4051 9.7459 25.80
Days 10.4673 0.8042 13.02
HPDinterval(VarCorr(fm1.mcmc,type="varcov"),prob=0.95)
lower upper
[1,] 307.8815 868.9589
[2,] 862.6868 1363.5824
attr(,"Probability")
[1] 0.95
Cheers,
Jarrod
On 11 Aug 2011, at 17:36, Emma Jones wrote:
I posted this message back in February and got no response....but am still puzzled. Can anyone help? Thanks Emma Hi, I am new to using the package lme4 and wondered if anyone can help me? I have fitted a model (note the data is unbalanced)- see below with the results:
/ a=lmer(RW~1+(1|Year),chron,na.action=na.omit)
/>/ a
/Linear mixed model fit by REML
Formula: RW ~ 1 + (1 | Year)
Data: chron
AIC BIC logLik deviance REMLdev
8771 8790 -4383 8761 8765
Random effects:
Groups Name Variance Std.Dev.
Year (Intercept) 0.48050 0.69318
Residual 0.52704 0.72598
Number of obs: 3719, groups: Year, 239
Fixed effects:
Estimate Std. Error t value
(Intercept) -0.02793 0.04803 -0.5815
This I am happy with and I know is correctly fitted- my problems don't
lie here.
What I would ideally like now is to produce 95% credible intervals on
the variance components. I have found some code on a help sheet online
and followed this...
samp=mcmcsamp(a,n=10000)
#gives a credible interval for mcmc
HPDinterval(VarCorr(samp,type="varcov"),prob=0.95)
lower upper
[1,] 0.1958900 0.2678522
[2,] 0.5317379 0.5853902
attr(,"Probability")
0.95
As it can be seen, both of my point estimates lie outside the credible
interval.
I am confused by this, can anyone advise??
Many thanks in advance,
Emma
[[alternative HTML version deleted]]
_______________________________________________ R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
The University of Edinburgh is a charitable body, registered in Scotland, with registration number SC005336.