The effect of the number of replicates in different cells on the results of ANOVA

I have written 3way-anova in c++. I have 3 factors, lets say a,b and c and my aim to check the strength of all possible interactions and main effects.The result of my code is the same as MATLAB when I use type I sum of squares.
But when I change the data so that the number of replicates is so high in some cells and too low in others: I don't get the same results as MATLAB. (to be precise, only SSt, SSe and SSa are the same as MATLAB) . My question is that is it possible that since I have a large difference in the number of replicates, I should use type III sum of squares? Or is there is special way that matlab treats the data in such cases so its results differs from mine?

-
did you implement the test in Matlab yourself, or are you checking against code from the statistics toolbox? At risk of adding a second gold standard, have you checked it against R's implementation? –  shabbychef May 24 '11 at 22:41
@shabbychef I implemented the code in c++. and the results is the same as matlab when I have between 1 to 5 replicates per level. but when this range changes to 1 to 100, some of my SS values become higher in comparison to corresponding SS values of MATLAB. and since each SS value is computed based on the others, some of them become negative! It does not happen in matlab. So I know that something is wrong in my code. –  Pegah May 25 '11 at 5:55

I don't have Matlab but from what I've read in the on-line help for N-way analysis of variance it's not clear to me whether Matlab would automatically adapt the type (1--3) depending on your design. My best guess is that yes you got different results because the tests were not designed in the same way.

Generally, with an imbalanced design it is recommended to use Type III sum of squares (SS), where each term is tested after all other (the difference with Type II sum of squares is only apparent when an interaction term is present), while with an incomplete design it might be interesting to compare Type III and Type IV SS. Note that the use of type III vs. Type II in the case of unbalanced data is subject to discussion in the literature.

(The following is based on a French tutorial that I cannot found anymore on the original website. Here is a personal copy, and here is another paper that discussed the different ways to compute SS in factorial ANOVAs: Which Sums of Squares Are Best In Unbalanced Analysis of Variance?)

The difference between Type I/II and Type III (also called Yates's weighted squares of means) lies in the model that serves as a reference model when computing SS, and whether factors are treated in the order they enter the model or not. Let's say we have two factors, A and B, and their interaction A*B, and a model like y ~ A + B + A:B (Wilkinson's notation).

With Type I SS, we first compute SS associated to A, then B, and finally A*B. Those SS are computed as the difference in residual SS (RSS) between the largest model omitting the term of interest and the smallest one including it.

For Type II and III, SS are computed in a sequantial manner, starting with those associated to A*B, then B, and finally A. For A*B, it is simply the difference between the RSS in the full model and the RSS in the model without interaction. The SS associated to B is computed as the difference between RSS for a model where B is omitted and a model where B is included (reference model); with Type III SS, the reference model is the full model (A+B+A*B), whereas for Type I and II SS, it is the additive model (A+B). This explains why Type II and III will be identical when no interaction is present in the full model. However, to obtain the first SS, we need to use dummy variables to code the levels of the factor, or more precisely difference between those dummy-coded levels (which also means that the reference level considered for a given factor matters; e.g., SAS consider the last level, whereas R consider the first one, in a lexicographic order). To compute SS for the A term, we follow the same idea: we consider the difference between the RSS for the model A+B+A*B and that for the reduced model B+A*B (A omitted), in case of Type III SS; with Type II SS, we consider A+B vs. B.

Note that in a complete balanced design, all SS will be equal. Moreover, with Type I SS, the sum of all SS will equal that of the full model, whatever the order of the terms in the model is. (This is not true for Type II and Type III SS.)

A detailed and concrete overview of the different methods is available in one of Howell's handout: Computing Type I, Type II, and Type III Sums of Squares directly using the general linear model. That might help you check your code. You can also use R with the car package, by John Fox who dicussed the use of incremental sum of squares in his textbook, Applied Regression Analysis, Linear Models, and Related Methods (Sage Publications, 1997, § 8.2.4--8.2.6). An example of use can be found on Daniel Wollschläger website.

Finally, the following paper offers a good discussion on the use of Type III SS (§ 5.1):

Venables, W.N. (2000). Exegeses on Linear Models. Paper presented to the S-PLUS User’s Conference Washington, DC, 8-9th October, 1998.

(See also this R-help thread, references therein, and the following post Anova – Type I/II/III SS explained.)

-
Good references! I think this one is also worth reading: Blair & Higgins, "Test of hypotheses for unbalanced factorial designs under various regression/coding method combinations" (1977) epm.sagepub.com/content/38/3/621.full.pdf+html Like Venables, they stress that type III model comparisons give different results for different coding methods. For correct type III SS, one needs a sum-to-zero coding scheme like effect coding or Helmert. Dummy coding (default in R) gives the wrong results. –  caracal May 25 '11 at 8:58
@caracal (+1) Thanks for the reference. I'm pretty sure David Howell also discussed this issue at length in his textbook, but I cannot retrieve it actually (should learn to clean up my desk one day :-) –  chl May 25 '11 at 9:01
Dear chl, thank you so much for the very good references. I was trying to find the exact mathematical formula for each of these three types. Now I am using the formulas in Karpinski page 11 and 12. When I use these formulas, my SSb and SSc becomes higher than they should be (in comparison to MATLAB), and since SSab, SSbc, SSac, and SSabc are calculated based on them, they even become negative sometimes. Is it possible that a high number of replicate has such an effect? –  Pegah May 27 '11 at 11:34