Re: [R] ANOVA and lmer

From: Douglas Bates <bates_at_stat.wisc.edu>
Date: Fri, 8 Feb 2008 08:58:38 -0600

On Feb 7, 2008 3:43 PM, Eric Imbert <eric.imbert_at_univ-montp2.fr> wrote:

> I am analyzing from a very simple experiment.
> I have measured plants of two different colours (yellow and purple) in 9
> different populations.

> So, I have two different factors : a fixed effect (Colour with two
> levels) and a random one (Population with 9 levels).

> I first analyzed the data with the aov function
> LargS is the variable

> aov(formula = LargS ~ Col + Error(Col/Pop))

> Terms:
> Col
> Sum of Squares 3.440351
> Deg. of Freedom 1

> Estimated effects are balanced

> Stratum 2: Col:Pop

> Terms:
> Residuals
> Sum of Squares 3017.112
> Deg. of Freedom 16

> Residual standard error: 13.73206

> Stratum 3: Within

> Terms:
> Residuals
> Sum of Squares 3347.385
> Deg. of Freedom 302

> To test for the interaction Col*Pop, I used the following F-ratio =
> (3017/16)/(3347/302) = 188. Highly significant !

> Now, let's go to the analysis performed by lmer - First I do the linear
> model without the Col*Pop interaction :
> m3=lmer(LargS ~ Col + (1 | Pop)

> And next with the interaction : m2=lmer(LargS ~ Col + (Col | Pop))

I don't think this model provides a comparison that is directly comparable with the test you describe above.

There is a model that is intermediate to your m2 and m3 models in terms of the number of parameters. (Actually when I see these results I realize that there is a bug in anova for lmer models in that the count of the parameters 1 less than it should be. I forgot to include the residual variance as a parameter to be estimated. However, this will not affect comparisons here because it is only the differences in the counts that matter.)

You model m2 is shown as having 5 parameters, which should be 6 (2 fixed effects, 3 variance components for the random effects, 1 residual variance estimate), while m2 is shown as having 3 parameters (should be 4: 2 f.e. + 1 r.e. var + 1 resid var). The model

LargS ~ Col + (1|Pop) + (1|Col:Pop)

which allows for random effects for each population and for each Col:Pop combination, has 5 parameters (currently shown as 4). The random effects for Pop are independent of each other and with a constant variance while the random effects for the Col:Pop combinations are independent of each other and of the Pop random effects and with another constant variance. Thus there are 2 f.e. + 2 r.e. var + 1 resid. var.

I think it is important to plot the data first and decide which of these models is indicated by the data. It sounds like the structure of your experiment is similar to the structure of the "barley" data set in the lattice package. Deepayan Sarkar's forthcoming book "Lattice: Multivariate Data Visualization with R" (to be published by Springer and due to ship next month) describes very effective ways of plotting such data. You can get a preview at his web site for the book http://lmdvr.r-forge.r-project.org

> Comparing both models : anova(m2,m3) :
>
> Df AIC BIC logLik Chisq Chi Df Pr(>Chisq)
> m3 3 1710.67 1721.97 -852.33
> m2 5 1714.59 1733.43 -852.30 0.0746 2 0.9634
>
> => Conclusion : the interaction Col*Pop is not significant !
>
> I guess I am missing something.
> Who can help ?
>
> Eric
>
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help_at_r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



R-help_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Fri 08 Feb 2008 - 15:04:12 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Fri 08 Feb 2008 - 16:30:13 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive