[R-sig-ME] Side effects from mcmcsamp()

Juan Pedro Steibel steibelj at msu.edu
Fri May 15 22:16:53 CEST 2009


Hello,
I noticed this a wile ago. One thing that occurred to me was to copy to 
another object before running mcmcsamp, example:
obmix0<-obmix
where obmix is the lmer object.
then do: mcmcsamp(obmix)....
For some reason obmix0 gets modified too.

It's been a while since I checked this so it may have been fixed.

Anyways, a way around this is to extract all the info needed from your 
lmer object to other variables before running mcmcsamp, that's what I do 
now and it works fine.
Thanks
JP



John Maindonald wrote:
> Dear Douglas -
> The following is at the very least a serious trap!  The output
> from VarCorr() and print() changes remarkably after running
> mcmcsamp() with the lmer object as argument.  The estimate
> of sigma is always (for these data) high, but roughly in the
> middle of the range of values that come out of ant111b.samp at sigma
>
> I have also been able to reproduce this behaviour
> (a) under lme4_0.999375-28 [below, I use 30]
> (b) under Mac OSX.  In fact, I thought initially that this
> was an OSX problem!
>
> While mcmcsamp() is in view, it would be very helpful
> to have a progress report on where it is at.
>
> Many thanks
> John Maindonald.
>
> I run the following code:
>
> library(DAAG)
> library(lme4)
> ant111b.lmer <- lmer(harvwt ~ 1 + (1 | site), data=ant111b)
> ant111b.lmer
> VarCorr(ant111b.lmer)
> ant111b.samp <- mcmcsamp(ant111b.lmer, n=1000)
>
> ## Now, extract the VarCorr and summary correlations
> VarCorr(ant111b.lmer)
> summary(ant111b.lmer)
>
> > ant111b.lmer <- lmer(harvwt ~ 1 + (1 | site), data=ant111b)
> > ant111b.lmer
> Linear mixed model fit by REML
> Formula: harvwt ~ 1 + (1 | site)
>   Data: ant111b
>   AIC   BIC logLik deviance REMLdev
> 100.4 104.8 -47.21    95.08   94.42
> Random effects:
> Groups   Name        Variance Std.Dev.
> site     (Intercept) 2.36773  1.53874
> Residual             0.57754  0.75996
> Number of obs: 32, groups: site, 8
>
> Fixed effects:
>            Estimate Std. Error t value
> (Intercept)   4.2917     0.5603   7.659
> > library(DAAG)
> > library(lme4)
> > ant111b.lmer <- lmer(harvwt ~ 1 + (1 | site), data=ant111b)
> > ant111b.lmer
> Linear mixed model fit by REML
> Formula: harvwt ~ 1 + (1 | site)
>   Data: ant111b
>   AIC   BIC logLik deviance REMLdev
> 100.4 104.8 -47.21    95.08   94.42
> Random effects:
> Groups   Name        Variance Std.Dev.
> site     (Intercept) 2.36773  1.53874
> Residual             0.57754  0.75996
> Number of obs: 32, groups: site, 8
>
> Fixed effects:
>            Estimate Std. Error t value
> (Intercept)   4.2917     0.5603   7.659
> > ant111b.samp <- mcmcsamp(ant111b.lmer, n=1000)
> > VarCorr(ant111b.lmer)
> $site
>            (Intercept)
> (Intercept)    8.363436
> attr(,"stddev")
> (Intercept)
>   2.891961
> attr(,"correlation")
>            (Intercept)
> (Intercept)           1
>
> attr(,"sc")
> sigmaREML
> 1.428289
> > summary(ant111b.lmer)
> Linear mixed model fit by REML
> Formula: harvwt ~ 1 + (1 | site)
>   Data: ant111b
>   AIC   BIC logLik deviance REMLdev
> 112.1 116.5 -53.04    105.7   106.1
> Random effects:
> Groups   Name        Variance Std.Dev.
> site     (Intercept) 8.3634   2.8920
> Residual             2.0400   1.4283
> Number of obs: 32, groups: site, 8
>
> Fixed effects:
>            Estimate Std. Error t value
> (Intercept)   4.2917     0.4171   10.29
>
> > sessionInfo()
> R version 2.9.0 (2009-04-17)
> i386-pc-mingw32
>
> locale:
> LC_COLLATE=English_United States.1252;LC_CTYPE=English_United 
> States.1252;LC_MONETARY=English_United 
> States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252
>
> attached base packages:
> [1] stats     graphics  grDevices utils     datasets  methods   base
>
> other attached packages:
> [1] lme4_0.999375-30   Matrix_0.999375-26 lattice_0.17-22
> [4] DAAG_0.99-1        MASS_7.2-47
>
> loaded via a namespace (and not attached):
> [1] grid_2.9.0
>
>
>
> John Maindonald             email: john.maindonald at anu.edu.au
> phone : +61 2 (6125)3473    fax  : +61 2(6125)5549
> Centre for Mathematics & Its Applications, Room 1194,
> John Dedman Mathematical Sciences Building (Building 27)
> Australian National University, Canberra ACT 0200.
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>


-- 
=============================
Juan Pedro Steibel

Assistant Professor
Statistical Genetics and Genomics

Department of Animal Science & 
Department of Fisheries and Wildlife

Michigan State University
1205-I Anthony Hall
East Lansing, MI
48824 USA 

Phone: 1-517-353-5102
E-mail: steibelj at msu.edu




More information about the R-sig-mixed-models mailing list