[R-sig-ME] "mixed" MANOVAs
peter dalgaard
pdalgd at gmail.com
Sat Jun 30 21:22:48 CEST 2012
On Jun 29, 2012, at 19:08 , Andrés Egea wrote:
> Dear all:
>
>
> I need to perform a MANOVA including both fixed and random factors. I know to run a traditional MANOVA (i.e., all factors are fixed) and even a MANCOVA (it seems that manova() function admits covariates), but not how to include random factors. Somewhere I read that MANOVA does not admit random factors, although I have no clue on the reason. Is it that true?
Well, according to ?manova,
‘manova’ does not support multistratum analysis of variance, so
the formula should not include an ‘Error’ term.
There is no specific reason that it cannot work (for balanced error designs, but that's true for aov as well) other than the fact that noone has implemented it. (It's what some sarcastically call a SMOP -- Simple Matter of Programming...)
And true enough, if you try:
> example(aov)
> npk$foo <- rnorm(24)
> m <- manova(cbind(yield,foo) ~ N*P*K + Error(block), npk)
Warning message:
In aov(cbind(yield, foo) ~ N * P * K + Error(block), npk) :
Error() model is singular
or:
> m <- aov(cbind(yield,foo) ~ N*P*K + Error(block), npk)
Warning message:
In aov(cbind(yield, foo) ~ N * P * K + Error(block), npk) :
Error() model is singular
(The two do basically the same thing, manova() does little more than adding an extra class to the result.)
The warning is a bit obscure; it probably just means that something inside aov() didn't expect a matrix response.
HOWEVER, it turns out that aov has done all the right calculations!!! The result (in m) is composed of a list of objects of class "mlm", for which anova.mlm (which didn't exist at the time manova() was written) knows how do the right thing:
> lapply(m[-1],anova)
$block
Analysis of Variance Table
Df Pillai approx F num Df den Df Pr(>F)
N:P:K 1 0.38534 0.94038 2 3 0.4819
Residuals 4
$Within
Analysis of Variance Table
Df Pillai approx F num Df den Df Pr(>F)
N 1 0.58359 7.7083 2 11 0.008079 **
P 1 0.04395 0.2528 2 11 0.781005
K 1 0.44757 4.4560 2 11 0.038241 *
N:P 1 0.21078 1.4689 2 11 0.272010
N:K 1 0.15182 0.9845 2 11 0.404269
P:K 1 0.30928 2.4627 2 11 0.130665
Residuals 12
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
And lo and behold:
> m <- manova(cbind(yield,foo) ~ N*P*K + Error(block), npk)
Warning message:
In aov(cbind(yield, foo) ~ N * P * K + Error(block), npk) :
Error() model is singular
> summary(m)
Error: block
Df Pillai approx F num Df den Df Pr(>F)
N:P:K 1 0.38534 0.94038 2 3 0.4819
Residuals 4
Error: Within
Df Pillai approx F num Df den Df Pr(>F)
N 1 0.58359 7.7083 2 11 0.008079 **
P 1 0.04395 0.2528 2 11 0.781005
K 1 0.44757 4.4560 2 11 0.038241 *
N:P 1 0.21078 1.4689 2 11 0.272010
N:K 1 0.15182 0.9845 2 11 0.404269
P:K 1 0.30928 2.4627 2 11 0.130665
Residuals 12
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
So it seems that things are pretty darn close to actually just working.
One caveat: I am pretty sure that the above results are more than superficially sane, but I haven't actually checked. For good measure, if you have a textbook MANOVA example, perhaps you should try that first and see if you can reproduce its test statistics.
-pd
> Anyway, could you tell a way (either the suitable package or even the script) to run MANOVAs with both fixed and random factors, please? My database is balanced.
>
> Any assisstance I can get will be greatly appreciated. Thank you very much in advance for your help.
>
> Looking forward to hearing from you, sincerely,
> Andrés Egea
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
--
Peter Dalgaard, Professor,
Center for Statistics, Copenhagen Business School
Solbjerg Plads 3, 2000 Frederiksberg, Denmark
Phone: (+45)38153501
Email: pd.mes at cbs.dk Priv: PDalgd at gmail.com
More information about the R-sig-mixed-models
mailing list