Am I correct to understand that the order in which variables are specified in a multifactorial ANOVA makes a difference but that the order does not matter when doing a multiple linear regression?

So assuming an outcome such as measured blood loss y and two categorical variables

  1. adenoidectomy method a,
  2. tonsillectomy method b.

The model y~a+b is different to the model y~b+a (or so my implementation in R seems to indicate).

Am I correct to understand that the term here is that ANOVA is a hierarchical model since it first attributes as much variance as it can to the first factor before trying to attribute residual variance to the second factor?

In the example above the hierarchy makes sense because I always do the adenoidectomy first before doing the tonsillectomy but what would happen if one had two variables with no inherent order?


2 Answers 2


This question evidently came from a study with an unbalanced two-way design, analyzed in R with the aov() function; this page provides a more recent and detailed example of this issue.

The general answer to this question, as to so many, is: "It depends." Here it depends on whether the design is balanced and, if not, which flavor of ANOVA is chosen.

First, it depends on whether the design is balanced. In the best of all possible worlds, with equal numbers of cases in all cells of a factorial design, there would be no difference due to the order of entering the factors into the model, regardless of how ANOVA is performed.* The cases at hand, evidently from a retrospective clinical cohort, seem to be from a real world where such balance was not found. So the order might matter.

Second, it depends on how the ANOVA is performed, which is a somewhat contentious issue. The types of ANOVA for unbalanced designs differ in the order of evaluating main effects and interactions. Evaluating interactions is fundamental to two-way and higher-order ANOVA, so there are disputes over the best way to proceed. See this Cross Validated page for one explanation and discussion. See the Details and the Warning for the Anova() (with a capital "A") function in the manual for the car package for a different view.

The order of factors does matter in unbalanced designs under the default aov() in R, which uses what are called type-I tests. These are sequential attributions of variance to factors in the order of entry into the model, as the present question envisioned. The order does not matter with the type-II or type-III tests provided by the Anova() function in the car package in R. These alternatives, however, have their own potential disadvantages noted in the above links.

Finally, consider the relation to multiple linear regression as with lm() in R, which is essentially the same type of model if you include interaction terms. The order of entry of variables in lm() does not matter in terms of regression coefficients and p-values reported by summary(lm()), in which a k-level categorical factor is coded as (k-1) binary dummy variables and a regression coefficient is reported for each dummy.

It is, however, possible to wrap the lm() output with anova() (lower-case "a," from the R stats package) or Anova() to summarize the influence of each factor over all of its levels, as one expects in classical ANOVA. Then the ordering of factors will matter with anova() as for aov(), and will not matter with Anova(). Similarly, the disputes over which type of ANOVA to use would return. So it's not safe to assume order-independence of factor entry with all downstream uses of lm() models.

*Having equal numbers of observations in all cells is sufficient but, as I understand it, not necessary for the order of factors to be irrelevant. Less demanding types of balance may allow for order-independence.

  • $\begingroup$ Indeed yes, that observational data was unbalanced, very unbalanced. $\endgroup$
    – Farrel
    Aug 6, 2018 at 0:41
  • $\begingroup$ Hopefully this comment still gets an answer here: You say that, under a balanced study design, the SS estimate will never be order dependent, irregardless of the type of anova test (typeI, II, III) chosen. I am not sure if I understand this. using the 'anova' function in R (which uses type I tests) on a linear model based on data that is balanced, surely the feature order matters, no? $\endgroup$
    – PejoPhylo
    Dec 19, 2019 at 13:20
  • 1
    $\begingroup$ @PejoPhylo when the data are balanced then you can have what is called an orthogonal design. With an orthogonal design there is one unique way to partition the sums of squares among the treatments and their interactions, so the order of entry of treatments will not matter with respect to estimates of effects and their p-values. This page provides a mathematical explanation. This is not immediately obvious; the question I just linked was asked by a member of this site with one of the highest reputations. Unbalanced data can destroy orthogonality. $\endgroup$
    – EdM
    Dec 19, 2019 at 16:27
  • $\begingroup$ Thanks a lot for your answer @EdM $\endgroup$
    – PejoPhylo
    Dec 20, 2019 at 7:30
  • 1
    $\begingroup$ @stan I haven't worked through the calculations, but I suspect that comes from the bias introduced by omitting the intercept and how that plays out in an unbalanced design. You almost certainly should not be omitting an intercept. Search for "regression" "omit" "intercept" on this site. $\endgroup$
    – EdM
    Feb 4, 2021 at 13:14

The term hierarchical model refers to the structure between the factors. For example, a multi-center study is hierarchical: You have the patients nested within the hospitals treating them. Each hospital treats patients with placebo and verum, but recieving each of them in either hospital A or B is slightly different due to some common effect of the hospital governing on all their patients (might even be an interaction effect with the experimental agent). So it's called hierarchical effect.

Now your ectomy methods may be hierarchical: Is it plausible that a certain tonsillectomy method is slightly different (in itself, not yet in the effect, because that's what your are going to estimate and test) depending on the adenoidectomy method used prior on the same patient? If yes, you should specify it in your model.

Your observation that y~a+b may be different from y~b+a indicates that there is something wrong. Additive effects commute, so there should not be a difference (apart from small numerical differences). It is neither plausible nor desired that the effect of the surgery methods may depend on the order in which the statistician later specifies the effects. So you probably picked the wrong approach to feed R with the data.

  • 1
    $\begingroup$ I am not sure I follow the last paragraph. In unbalanced factorial ANOVA the p-values for each factor computed via Type I (sequential) sum of squares will certainly dependent on the order of the factors. I believe this is the whole point of the question. $\endgroup$
    – amoeba
    May 12, 2016 at 13:54
  • $\begingroup$ I'm not sure if @Farrel got Type I SS. I remember I once observed SAS to output different Type III SS due to some disparate sorting in the data set and model statement. Maybe this can happen with R, too? $\endgroup$ May 12, 2016 at 14:25
  • 2
    $\begingroup$ I cannot know for a fact and he might not remember himself given that the Q was asked five years ago. But I think this is by far the most parsimonious interpretation of his words "The model y~a+b is different to the model y~b+a (or so my implementation in R seems to indicate)", in particular given the fact that aov command in R uses Type I SS by default. When I offered the bounty, I expected to get an answer explaining the issues behind unbalanced anova design, differences between Type I/II/III SS, and some comments on whether linear regression does or doesn't have the same issues. $\endgroup$
    – amoeba
    May 12, 2016 at 14:31
  • 1
    $\begingroup$ No. Design matrix is singular in anova even if it is balanced, when there is no difference between SS I/II/III. The SS I/II/III are different only in the unbalanced case because the factors become non-orthogonal (unlike in the balanced case). In my understanding, this corresponds to a linear regression with correlated predictors, which is a very common situation. My answer is that the same issue occurs in regression too, it's just that it is standard to compute a p-value of one predictor after accounting for the effects of all other predictors; this corresponds to Type III SS in anova. $\endgroup$
    – amoeba
    May 12, 2016 at 14:47
  • 1
    $\begingroup$ Such questions about variable order in ANOVA keep coming in, like this one migrated from Stack Overflow yesterday. I think it's safe to assume that this 5-year-old question was similarly based on aov rather than lm, and it would be helpful to have an answer to this question of the type that @amoeba indicated in the comment from May 12, 14:31. $\endgroup$
    – EdM
    May 15, 2016 at 18:40

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.