The key is the `keep.order' argument to terms.formula:
> options(contrasts=c("contr.helmert", "contr.poly"))
> library(MASS)
> aov(yield ~ block + N * P + K, npk)
Call:
aov(formula = yield ~ block + N * P + K, data = npk)
Terms:
block N P K N:P Residuals
Sum of Squares 343.2950 189.2817 8.4017 95.2017 21.2817 218.9033
Deg. of Freedom 5 1 1 1 1 14
Residual standard error: 3.954232
Estimated effects are balanced> aov(terms(yield ~ block + N * P + K, keep.order=TRUE), npk)
Call:
aov(formula = terms(yield ~ block + N * P + K, keep.order = TRUE),
data = npk)
Terms:
block N P N:P K Residuals
Sum of Squares 343.2950 189.2817 8.4017 21.2817 95.2017 218.9033
Deg. of Freedom 5 1 1 1 1 14
Residual standard error: 3.954232
Estimated effects are balanced
Now, I wonder how I knew that ... oh yes, it is in the vr.glm demo.
On Thu, 18 Dec 2003, Pascal A. Niklaus wrote:
> Is there a way to prevent the re-ordering of factors by aov? I do have a
> three-way interaction that I do want to fit before a two-way interaction
> (different factors, so they are not nested), but R moves the two-way
> interaction to the front. I know it generally makes sense to fit the
> two-way interactions first, but in this case I think I know what I'm
doing.
>
> + hv
> + hv:rep
> + hv:soil
> + hv:spdiv
> + hv:spdivnom
> + hv:fgdiv
> + hv:fgdivnom
> + hv:soil:spdiv
> + hv:soil:spdivnom
> + hv:soil:fgdiv
> + hv:soil:fgdivnom
> + hv:mx #<-- that one here
> + hv:soil:mx
>
> mx is a random factor, and I do want to semi-manually test the three-way
> interactions against "hv:mx". soil is treated as fixed effect.
>
> I know I can fit the first part, then save the residuals, and then fit
> the last two terms, but it would be really nice not to have to do this
> in two steps.
--
Brian D. Ripley, ripley at stats.ox.ac.uk
Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel: +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UK Fax: +44 1865 272595