# Why do the results of a MANOVA change when the order of the predictor variables is changed?

So, for example, using the Iris data and treating iris species as the predictor variable and sepal length, sepal width, petal length, and petal width as the dependent variables we get MANOVA output that looks like this:

set.seed(2)
# Creating a matrix of the 4 dependent variables (DVs)
Y <- as.matrix(iris[,c(1:4)])

# MANOVA looking at the effect of species on DVs
summary(manova(Y ~ iris$Species)) # Df Pillai approx F num Df den Df Pr(>F) # iris$Species    2   1.1919    53.466    8       290       < 2.2e-16 ***


That seems to make sense. Species has a significant effect on our DVs. Now what if we add another predictor variable (a random one which we shouldn’t expect to have an effect on the DVs)?

# Creating a random dummy variable to be used as a predictor variable
iris$random.dummy <- sample(x = c(0,1), size = 150, replace = TRUE) # MANOVA looking at the effect of species + our random dummy on DVs summary(manova(Y ~ iris$Species + iris$random.dummy)) # Df Pillai approx F num Df den Df Pr(>F) # iris$Species        2   1.19339   53.263    8       288     <2e-16 ***
# iris$random.dummy 1 0.03784 1.406 4 143 0.2349  That also seems to make sense. Species is significant still, but our random dummy variable is not. Now what if we simply switch the order of those variables? # Switching the order of our two predictor variables in the formula summary(manova(Y ~ iris$random.dummy + iris$Species)) # Df Pillai approx F num Df den Df Pr(>F) # iris$random.dummy   1   0.13031   5.357     4       143     0.0004764 ***
# iris$Species 2 1.19526 53.470 8 288 < 2.2e-16 ***  Now, the Pillai’s trace and approximate F-values change and our random dummy variable has become significant. So my questions are these. Why do the results of a MANOVA change when the order of the predictor variables is changed? and What does this mean for those of us trying to use and interpret a MANOVA? • The answer here stats.stackexchange.com/questions/11127 is very relevant and explains large parts of this conundrum. Does not explain it completely though: I don't understand why the outcome of summary(manova(Y ~ iris$random.dummy + iris$Species)) differs from summary(manova(Y ~ iris$random.dummy)). Dec 12, 2015 at 22:21
• Update: it's because the error term is not the same. I added [anova] tag to your question, because this issue is not specific to MANOVA. Dec 13, 2015 at 1:12
• Note that there isn’t balance in the species and random factors (table(iris$random.dummy, iris$Species)). ANOVA and MANOVA for unbalanced data will always be problematic. If you instead generate random data that are balanced across species (c(sample(rep(0:1, each=25)), sample(rep(0:1, each=25)), sample(rep(0:1, each=25)))), you will not have the problem of inconsistent results for the random term from the three models (including one without Species). (Because of the error term, the results will not be identical, but they will be very similar, and usually not statistically significant.) Dec 13, 2015 at 10:28
• Just to clarify my earlier comment: With balanced data, you will get identical results for the effect of random.dummy regardless of the order of variables, but you will not get identical results (though very similar) if you exclude the Species variable. For unbalanced data, expect (possibly) very different results for random.dummy for the three models. This demonstrates how difficult it is to correctly interpret the results of (M)ANOVA models for unbalanced data. Dec 13, 2015 at 10:35
• Now that I have had the chance to take a closer look at the link amoeba shared, I see the second issue has to do with how the restricted vs. unrestricted models are set up in base R. The car package MANOVA function gave me the same output regardless of order. Thanks all for the assistance! Dec 14, 2015 at 18:22