Re: [R] nested anova and multiple comparisons

From: Kingsford Jones <kingsfordjones_at_gmail.com>
Date: Sat, 26 Apr 2008 14:52:48 -0700

Hi Stephen,

On Sat, Apr 26, 2008 at 10:29 AM, Stephen Cole <swbcole_at_gmail.com> wrote: ...<snip>...
>
> I have managed to accomplish my first two goals by analyzing the data
> as a 3 level nested Anova with the following code
>
> mod1 <- aov(Count~Region/Location/Site, data=data)
>
> This allows me to get the MS for a Anova table. However R does not
> compute the correct F-statistics (because it uses the residual MS for
> the denominator in all calculations which it shouldnt) so i manually
> computed the F-stats and variance components by hand.
>

R's just doing what it was asked it to do. The fact that it used the residual MS isn't a surprise given your code.

> >From reading the help guide i learned about and tried using the
> Error(Region/Location/Site) command but then i can only get MS and no
> F-statistics and still hand to compute the rest by hand.

Yes, if you want to use Method-of-Moments estimation to solve for expected mean squares for the variance components and calculate F-stats w/ the 'correct' MS, then then specifying error strata is reasonable for balanced data (assuming you've met model assumptions). However, I believe most statisticians these days are more inclined to use (Restricted) Maximum Likelihood estimation (e.g. using lme or lmer) because of the added fliexibility it provides (also it won't produce negative variance estimates when the mean squared error is larger than the mean square between the groups of interest)

As for why you are only getting MS and no F-stats, it's hard to say without a reproducible example (see the posting guide). In my experience this will occur when there are insufficient degrees of freedom for calculating an F-stat at a given level.

>
> My problem now is that i would liek to use TukeyHSD for multiple
> comparisons. Howeber since R is unable to compute the correct F
> statistics in this nested design i also think it is using the wrong MS
> and df in calculating tukeys q. for example when i use
>

Again, it's not that R is unable to, it's that you asked R not to.

> TukeyHSD(mod1, "Region") i will get values however i do not think
> they have been calculated correctly.
>
> Furthermore when i use the Error(Region/Location/Site) term i can then
> no longer use TukeyHSD as i get the error message that there is no
> applicable use of tukey on this object.

Looking at the methods available for TukeyHSD shows that it will work for an object of class "aov"

> methods(TukeyHSD)

[1] TukeyHSD.aov

But ?aov states:

    Value:

         An object of class 'c("aov", "lm")' or for multiple responses of
         class 'c("maov", "aov", "mlm", "lm")' or for multiple error strata
         of class '"aovlist"'.  There are 'print' and 'summary' methods
         available for these.

So your model with error strata is of class "aovlist" not "aov".

>
> i am just wondering if there is any way to use Multiple comparisons
> with R in a nested design like mine.
>

I'm not sure but the package 'multcomp' might be of help.

> I have thought about using lme or lmer but if i understand them right
> with a balanced design i should be able to get the same result using
> aov
>

Even with balanced data, you don't get the exact same thing with aov. As mentioned above lme and lmer use different estimation methods. One of the advantages of using (RE)ML is a lot more flexibility in how you specifiy the structure of random effects covariance matrices, and (at least with lme) you have a lot of flexibility in how you structure the residual covariance matrix as well. This is particularly useful when you are not meeting assumptions of constant variance and/or independence of observations (which seems to be the rule rather than the exception with ecological data with a spatial component, such as you appear to have).

The theory behind analyses you want to do are quite complex with many potential pitfalls. If at all possible, I highly recommend finding the help of someone who is an expert in these types of models (known as mixed, hierarchical, multilevel, random effects, variance components, nested ANOVA, etc....).

best,

Kingsford Jones

> Thanks
>
> ______________________________________________
> 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 Sat 26 Apr 2008 - 21:56:43 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 Sat 26 Apr 2008 - 22:30:32 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