[R] define F-ratio computations with aov

Galanidis Alexandros agal at env.aegean.gr
Thu Mar 25 15:49:19 CET 2010


Many thanks to both of you!

I used the summary(obj.aov)$Mean way and I have the correct values now, although slightly different from the analogous SPSS computations. 

I'll try splplot and Error methods in a calmer period, but I have the feeling that in such cases aov is not the best method anyways and one should proceed to lme4 package.

aleXis

________________________________________
Από: RICHARD M. HEIBERGER [rmh(at)temple.edu]
Αποστολή: Τετάρτη, 17 Μαρτίου 2010 4:21 μμ
Προς: Galanidis Alexandros
Κοιν.: R Maillist
Θέμα: Re: [R] define F-ratio computations with aov

Please look at the Error term in ?aov

On Wed, Mar 17, 2010 at 4:25 AM, Galanidis Alexandros  wrote:
Greetings to all,

This is my model: aov.fit<-aov(Y~A+B+C+D+E+A:C+A:E)

In summary(aov.fit) all F values are comptuted by eg MS(A)/MS(Residuals). This is not correct (or what I want), except for F(B) and F(A:E). I suppose P values are not correct either.

Is it possible with aov to define the way F computations will be done? I 'd like them to be like this: F(A)=MS(A)/MS(E), F(C)=MS(C)/MS(E), F(D)=MS(D)/MS(E), F(E)=MS(E)/MS(A:E), F(A:C)=MS(A:C)/MS(A:E)

For your example, it looks like

mydata.aov <- aov(Y ~ A+B+C+D + Error(E/A), data=mydata)

might be what you need.  It isn't possible to be sure of the correct statement
for your example without seeing the actual treatment plan.

The standard split-plot design is specified with

splplot <- data.frame(y=rnorm(72),
                      blocks=factor(rep(1:6, each=12)),
                      plots=factor(rep(rep(1:3, each=4), 6)),
                      subplots=factor(rep(1:4, 18)),
                      A=factor(rep(c(3,1,2, 3,1,2, 2,3,1, 3,2,1, 2,1,3, 1,2,3), each=4)),
                      B=factor(
                        c(4, 3, 2, 1, 1, 2, 4, 3, 1, 2, 3, 4,
                          3, 1, 2, 4, 4, 1, 2, 3, 2, 1, 3, 4,
                          2, 3, 4, 1, 4, 2, 3, 1, 1, 4, 2, 3,
                          3, 4, 1, 2, 1, 3, 4, 2, 2, 3, 4, 1,
                          4, 1, 3, 2, 3, 4, 1, 2, 3, 4, 2, 1,
                          3, 1, 4, 2, 4, 3, 1, 2, 1, 2, 3, 4)))

splplot.aov <- aov(y ~ A*B +  Error(blocks/plots/subplots), data=splplot)
summary(splplot.aov)
> summary(splplot.aov)
Error: blocks
          Df Sum Sq Mean Sq F value Pr(>F)
Residuals  5 2.0231 0.40463
Error: blocks:plots
          Df  Sum Sq Mean Sq F value Pr(>F)
A          2  0.2876  0.1438  0.1135 0.8938
Residuals 10 12.6646  1.2665
Error: blocks:plots:subplots
          Df Sum Sq Mean Sq F value Pr(>F)
B          3  2.814 0.93797  0.8987 0.4493
A:B        6  3.358 0.55964  0.5362 0.7778
Residuals 45 46.965 1.04367
>

As you see, the correct F tests are automatically determined and displayed.

Rich

________________________________________
Από: Michael Rennie [mdrennie(at)gmail.com]
Αποστολή: Τετάρτη, 17 Μαρτίου 2010 2:32 μμ
Προς: Galanidis Alexandros
Κοιν.: R Maillist
Θέμα: Re: [R] define F-ratio computations with aov

Howdy,

In the past, I've just run the ANOVA as normal, and then just grabbed
the appropriate MS for the estimation of F ratios. Eg, this will get you
the MS in your anova object:

summary(obj.aov)[[1]][3]

or

summary(obj.aov)$Mean

And if you want a specific MS,

summary(obj.aov)[[1]][[1,3]]

or

summary(obj.aov)[[1]]$Mean[1]


Then you can just put whichever MS over whichever other MS, estimate
your F-ratios, with something like:

Ffact<- summary(obj.aov)[[1]]$Mean[1]/summary(obj.aov)[[1]]$Mean[3]

estimate the p-values with:

pFfact<-1-pf(Ffact, summary(obj.aov)[[1]]$Df[1],
summary(obj.aov)[[1]]$Df[3])

  and you're off to the races.

You can also specify error strata in the aov() model, but then all you
get is the MS and you have to estimate your F-ratios anyway (though the
indexing is a little different). E.g., if you had a nested anova, you
could specify it as:

ex.aov<-aov(Fixed ~ Nested + Error(Nested/Fixed))

At least this way, the summary() doesn't give you the wrong F-ratios, so
you aren't temped to interpret them incorrectly (as you would in the
previous example).

HTH,

Mike


More information about the R-help mailing list