Recall that the alternative hypothesis is usually what we suspect to be true and hope to conclude. Typically we are looking to find differences among at least one pair of our treatment means. Because of this, the null hypothesis (the opposite of the alternative) states that there are no differences among the treatment group means.

The idea behind ANOVA methods is to compare different sources of variability: between sample variability and within sample variability. To test the Null hypothesis, traditionally written as \(H_0 \colon \mu_1 = \mu_2 = ⋯ = \mu_T\), we need to compute a test (F) statistic that compares the between sample variability to within sample variability.

To understand the computation of this statistic it is helpful to look at the ANOVA table. The table below is an example of a blank (no entries yet) ANOVA table.

ANOVA

To define the elements of the table and fill in these quantities, let’s return to our example data (Lesson1 Data) for the hypothetical greenhouse experiment:

Notation

Each observation in the dataset can be referenced by two indicator subscripts, i and j as \(Y_{ij}\)

For those of you not familiar with this notation, we use Y to indicate that it is a response variable. The subscript i refers to the \(i^{th}\) level of the treatment; our example has 4 treatments so i will take on the values 1, 2, 3, and 4. The subscript j refers to the \(j^{th}\) observation; our example has 6 observations for each treatment so j takes the values 1, 2, 3, 4, 5, and 6. It is important to note that the \(j^{th}\)observation is occurring within the \(i^{th}\) treatment level. For example, it can be seen in the table below that \(Y_{4,2} = 27.5\).

We now can define the various means explicitly using these subscripts. The overall or Grand Mean is given by

\(\text{Grand Mean }=\bar{Y}_{..}\)

where the dots indicate that the quantity has been averaged over that subscript. For the Grand Mean, we have averaged over all j observations in all i treatment levels. Alternatively, the treatment means are given by

\(\text{Treatment Mean }=\bar{Y}_{i.}\)

Indicating that we have averaged over the j observations in each of the i treatment levels.

The means can be found in the output from the summary procedure generated in SAS, seen below. These and other coding details will be discussed in Lesson 3.

In the output, we see the column heading _TYPE_. The summary procedure in SAS calculates all possible means when specified, thus the _TYPE_ indicates what mean is being computed. _TYPE_ = 0 is the Grand Mean, and we can see this from the number of observations (given by _FREQ_) of 24. Further, each of the treatment means is listed as _TYPE_ = 1, and we confirm that 6 replications were made for each treatment level (remember that j took on values 1 through 6). Note that SAS automatically has ordered the treatment levels alphabetically.

In this example, the grand mean and treatment means are all we need to compute the quantities for the ANOVA table, continued in the next section.

When working with ANOVA, we start with the total variability in the response variable and divide or "partition" it into different parts: the variability due to our treatment (i.e. the between sample variability) and the leftover or residual variability (i.e. the within sample variability). The variability that is due to our treatment we hope is significantly large whereas the variability in the response that is leftover can be thought of as the nuisance, "error", or "residual" variability.

To imagine this further, think about the data storage capacity of a computer. If you have 8GB of storage total, you can ask your computer to show the types of files that are occupying the storage. We can visualize the response variability similarly, seen below. The ANOVA model is (in a very elementary fashion) going to compare the variability due to the treatment to the variability left over.

From elementary statistics, when we think of computing a variance of a random variable (say X), we use the expression:

\(\text{variance }=\dfrac{\sum(X_i-\bar{X})^2}{N-1}=\dfrac{SS}{df}\)

The numerator of this expression is referred to as the Sum of Squares, or Sum of Squared deviations from the mean, or simply SS. (If you don't recognize this then we suggest you sharpen your introductory statistics skills!) The denominator is the degrees of freedom, (N - 1), or df.

The ANOVA table is set up to generate quantities analogous to the simple variance calculation above. In our greenhouse experiment example:

1. We start by considering the TOTAL variability in the response variable. This is done by calculating the SS_Total

   \(\text{Total SS }=\sum_{i}\sum_{j}(Y_{ij}-\bar{Y}_{..})^2\)

   = 312.47.

   The degrees of freedom for the Total SS is N - 1 = 24 – 1 = 23, where N is the total sample size.

2. Our next step determines how much of the variability in Y is accounted for by our treatment. We now calculate SS_Treatment or SS_Trt:

   \(\text{Treatment SS }=\sum_{i}n_i(\bar{Y_{i.}}-\bar{Y}_{..})^2\)

   Note: The sum of squares for the treatment is the deviation of the group mean from the grand mean. So in some sense, we are "aggregating" all of the responses from that group and representing the "group effect" as the group mean.

   For our example:

   \begin{aligned}\text{Treatment SS }= 6*(21.0-26.1667)^2 + 6*(28.6-26.1667)^2 \\ +6*(25.8667-26.1667)^2 + 6*(29.2-26.1667)^2 = 251.44\end{aligned}

   Note that in this case we have equal numbers of observations (6) per treatment level, and it is, therefore, a balanced ANOVA.

3. Finally, we need to determine how much variability is "leftover". This is the Error or Residual sums of squares. SS_Error or SSE is calculated by subtraction:

   \(\text{Error SS }=\sum_{i}\sum_{j}(Y_{ij}-\bar{Y}_{i.})^2 = \text{ Total SS - Treatment SS}\)

   \(SS_{Error} = 312.47 – 251.44 = \mathbf{61.033}\)

   Note here that the "leftover" is really the deviation of any score from its group mean.

We can now fill in the table:

ANOVA

We have T treatment levels and so we use T - 1 for the df for the treatment. In our example, there are 4 treatment levels (the control and the 3 fertilizers) so T = 4 and T - 1 = 4 - 1 or 3. Finally, we obtain the error df by subtraction as we did with the SS.

The Mean Squares (MS) can now be calculated as:

\(MS_{Trt}=\dfrac{SS_{Trt}}{df_{Trt}}=\dfrac{251.44}{3}=83.813\)

and

\(MS_{Error}=\dfrac{SS_{Error}}{df_{Error}}=\dfrac{61.033}{20}=3.052\)

It is important to note that \(MS_{\text{Error}}\) is commonly referred to as \(MSE\). Also notice we do not need to calculate the \(MS_{\text{Total}}\).

ANOVA

Finally, we can compute the F statistic for our ANOVA. Conceptually we are comparing the ratio of the variability due to our treatment (remember we expect this to be relatively large) to the variability leftover, or due to error (and of course, since this is an error we want this to be small). Following this logic, we expect our F to be a large number. If we go back and think about the computer storage space we can picture most of the storage space taken up by our treatment, and less of it taken up by error. In our example, the F is calculated as:

\(F=\dfrac{MS_{Trt}}{MS_{Error}}=\dfrac{83.813}{3.052}=27.46\)

ANOVA

So how do we know if the F is "large enough" to conclude we have a significant amount of variability due to our treatment? We look up the critical value of F and compare it to the value we calculated. Specifically the critical F is \(F_\alpha = F_{(0.05, 3, 20)} = 3.10\). The critical value can be found using tables or technology. An example using SAS is seen below. 

data Fvalue;                                                                                                                               
                   q=finv(0.95, 3, 20);                                                                                                                   
                   put q=;                                                                                                                           
                run;
                proc print data=work.Fvalue;
                 run;

Most F tables actually index this value as \(1 - \alpha = .95\)

The \(F_{\text{calculated}}\) > \(F_\alpha\) so we reject \(H_0\) in favor of the alternative \(H_A\). The p-value (which we don't typically calculate by hand) is the area under the curve to the right of the \(F_{\text{calculated}}\) and is the way the process is reported in statistical software. Note that in the unlikely event that the \(F_{\text{calculated}}\) is exactly equal to the \(F_{\alpha}\) then the \(\text{p-value} = \alpha\). As the calculated F statistic increases beyond the \(F_{\alpha}\) and we go further into the rejection region, the area under the curve (hence the p-value) gets smaller and smaller. This leads us to the analogous decision rule: If the p-value is \(<\alpha\) then we reject \(H_0\).

If (and only if) we reject the null hypothesis, we then conclude at least one group is different from one other (importantly we do NOT conclude that all the groups are different).

If we reject the null, then we want to know WHICH group, or groups, are different. In our example we are not satisfied knowing at least one treatment level is different; we want to know where the difference is and the nature of the difference. To answer this question, we can follow up the ANOVA with a mean comparison procedure to find out which means differ from each other and which ones do not.

You might think we could not bother with the ANOVA and proceed with a series of t-tests to compare the groups. While that is intuitively simple, it creates inflation of the type I error. How does this inflation of type I error happen? For a single test, the probability of type I error is

\(\alpha=1-(.95) \)

The probability of committing a type I error for two simultaneous tests follows from the Multiplication Rule for independent events in probability. Recall that for two independent events A and B, the probability of A and B both occurring is P(A and B) = P(A) * P(B). So for two tests, we have

\(\alpha = 1 - ( (.95)*(.95) ) = 0.0975\)

which is now larger than the \(\alpha\) that we originally set. For our example, we have 6 comparisons, so

\(\alpha = 1 - (.95^6) = 0.2649\)

which is a much larger (inflated) probability of committing a type I error than we originally set (0.05). The multiple comparison procedures are designed to compensate for the type I error inflation (although each does so in a slightly different way).

There are several multiple comparison procedures that can be employed, but we will start with the one most commonly used, the Tukey procedure. In the Tukey procedure, we compute a "yardstick" value ( \(w\)) based on the \(MS_{\text{Error}}\) and the number of means being compared. If any two means differ by more than the Tukey \(w\) value, then we conclude they are significantly different.

1. Step 1: Compute Tukey’s \(w\) value

   \(w=q_{\alpha(p, df_{Error})}\cdot s_{\bar{Y}}\)

   where \(q_\alpha\) is obtained from a Table of Tukey \(q\) values,

   Tukey q Values Table  

   df for Error Term	\(\alpha\)	p = Number of Treatments
   		2	3	4	5	6	7	3	9	10
   5	0.05 0.01	3.64 5.70	4.6 6.98	5.22 7.80	5.67 8.42	6.03 8.91	6.33 9.32	6.58 9.67	6.80 9.97	6.99 10.24
   6	0.05 0.01	3.46 5.24	4.34 6.33	4.90 7.03	5.30 7.56	5.63 7.97	5.90 8.32	6.12 8.61	6.32 8.87	6.49 9.10
   7	0.05 0.01	3.34 4.95	4.16 5.92	4.68 6.54	5.06 7.01	5.36 7.37	5.61 7.68	5.82 7.94	6.00 8.17	6.16 8.37
   8	0.05 0.01	3.26 4.75	4.04 5.64	4.53 6.20	4.89 6.62	5.17 6.96	5.40 7.24	5.60 7.47	5.77 7.68	5.92 7.86
   9	0.05 0.01	3.20 4.60	3.95 5.43	4.41 5.96	4.76 6.35	5.02 6.66	5.24 6.91	5.43 7.13	5.59 7.33	5.74 7.49
   10	0.05 0.01	3.15 4.48	3.88 5.27	4.33 5.77	4.65 6.14	4.91 6.43	5.12 6.67	5.30 6.87	5.46 7.05	5.60 7.21
   11	0.05 0.01	3.11 4.39	3.82 5.15	4.26 5.62	4.57 5.97	4.82 6.25	5.03 6.48	5.20 6.67	5.35 6.84	5.49 6.99
   12	0.05 0.01	3.08 4.32	3.77 5.05	4.20 5.50	4.51 5.84	4.75 6.10	4.95 6.32	5.12 6.51	5.27 6.67	5.39 6.81
   13	0.05 0.01	3.06 4.26	3.73 4.96	4.15 5.40	4.45 5.73	4.69 5.98	4.88 6.19	5.05 6.37	5.19 6.53	5.32 6.67
   14	0.05 0.01	3.03 4.21	3.70 4.89	4.11 5.32	4.41 5.63	4.64 5.88	4.83 6.08	4.99 6.26	5.13 6.41	5.25 6.54
   15	0.05 0.01	3.01 4.17	3.67 4.84	4.08 5.25	4.37 5.56	4.59 5.80	4.78 5.99	4.94 6.16	5.08 6.31	5.20 6.44
   16	0.05 0.01	3.00 4.13	3.65 4.79	4.05 5.19	4.33 5.49	4.56 5.72	4.74 5.92	4.90 6.08	5.03 6.22	5.15 6.35
   17	0.05 0.01	2.98 4.10	3.63 4.74	4.02 5.14	4.30 5.43	4.52 5.66	4.70 5.85	4.86 6.01	4.99 6.15	5.11 6.27
   18	0.05 0.01	2.97 4.07	3.61 4.70	4.00 5.09	4.28 5.38	4.49 5.60	4.67 5.79	4.82 5.94	4.96 6.08	5.07 6.20
   19	0.05 0.01	2.96 4.05	3.59 4.67	3.98 5.05	4.25 5.33	4.47 5.55	4.65 5.73	4.79 5.89	4.92 6.02	5.04 6.14
   20	0.05 0.01	2.95 4.02	3.58 4.64	3.96 5.02	4.23 5.29	4.45 5.51	4.62 5.69	4.77 5.84	4.90 5.97	5.01 6.09
   24	0.05 0.01	2.92 3.96	3.53 4.55	3.90 4.91	4.17 5.17	4.37 5.37	4.54 5.54	4.68 5.69	4.81 5.81	4.92 5.92
   30	0.05 0.01	2.89 3.89	3.49 4.45	3.84 4.80	4.10 5.05	4.30 5.24	4.46 5.40	4.60 5.54	4.72 5.65	4.83 5.76
   40	0.05 0.01	2.86 3.82	3.44 4.37	3.79 4.70	4.04 4.93	4.23 5.11	4.39 5.27	4.52 5.39	4.63 5.50	4.74 5.60

2. and  \(p\) = the number of treatment levels
   \(s_\bar{Y}\) = standard error of a treatment mean = \(\sqrt{MS_{Error}/n}\)
   \(n\) = number of replications

   For our greenhouse example we get:

   \(w=q_{.05(4,20)}\sqrt{(3.052⁄6)}=3.96(0.7132)=2.824\)

3.  

4. Step 2: Rank the means, calculate differences

   For the greenhouse example, we rank the group means as:

   29.20	28.6	25.87	21.00

   Start with the largest and second-largest means and calculate the difference:

   \(29.20 – 28.60 = 0.60\) which is less than our w of 2.824, so we indicate there is no significant difference between these two means by placing the letter “a” under each.

   29.20	28.6	25.87	21.00
   a	a

   Then calculate the difference between the largest and third-largest means:

   \(29.20 – 25.87 = 3.33\) which exceeds the critical w of 2.824, so we can label the third group mean with a “b” to show this difference is significant.

   29.20	28.6	25.87	21.00
   a	a	b

   Now we have to consider whether or not the second-largest and third-largest differ significantly. This is a step that sets up a back and forth process. Here

   \(28.6 – 25.87 = 2.73\), less than the critical w of 2.824, so these two means do not differ significantly. We need to add a “b” to the second group mean to show this:

   29.20	28.6	25.87	21.00
   a	ab	b

   Continuing down the line, we now calculate the next difference:

   \(28.60 – 21.00 = 7.60\), exceeding the critical w, so we now add a “c”:

   29.20	28.6	25.87	21.00
   a	ab	b	c

   Again, we need to go back and check to see if the third-largest also differs from the smallest:

   \(25.87 – 21.00 = 4.87\), which it does. So we are done.

These letters can be added to figures summarizing the results of the ANOVA.

The Tukey procedure explained above is valid only with equal sample sizes for each treatment level. In the presence of unequal sample sizes, more appropriate is the Tukey–Cramer Method, which calculates the standard deviation for each pairwise comparison separately. This method is available in SAS, R, and most other statistical software.

Although the Tukey procedure is the most widely used multiple comparison procedure, there are many other multiple comparison techniques.

An older approach, no longer offered in many statistical computing packages, is Fisher’s Protected Least Significant Difference (LSD). This is a method to compare all possible means, two at a time, as t-tests. Unlike an ordinary two-sample t-test, however, the method does rely on the experiment-wide error (the MSE). The LSD is calculated as:

\(LSD(\alpha)=t_{\alpha,df}s_{\bar{d}}\)

where \(t_\alpha\) is based on \(\alpha\) and df = error degrees of freedom from the ANOVA table. The standard error for the difference between two treatment means (\(s_{\bar{d}}\) or SE) is calculated as:

\(s_{\bar{d}}=\sqrt{\dfrac{2s^2}{n}}\)

Where n is the number of observations per treatment mean (replications) and \(s^2\) is the MSE from the ANOVA. As in the Tukey method, any pair of means that differ by more than the LSD value differ significantly. The major drawback of this method is that it does not control \(\alpha\) over an entire set of pairwise comparisons (the experiment-wise error rate) and hence is associated with Type 1 inflation.

The following multiple comparison procedures are much more assertive in dealing with Type 1 inflation. In theory, while we can set \(\alpha\) for a single test, the fact that we have T treatment levels means there are T(T - 1)/2 tests (the number of pairs of possible comparisons), and so we need to adjust \(\alpha\) to have the desired confidence level for the set of tests. The Tukey, Bonferroni, and Scheffé methods control the experiment-wise error, but in different ways. All three use a

“multiplier * SE”

but differ in the form of the multiplier.

Contrasts are comparisons involving two or more factor level means (discussed more in the following section). Pairwise mean comparisons can be thought of as a subset of possible contrasts among the means. If only pairwise comparisons are made, the Tukey method will produce the narrowest confidence intervals and is the recommended method. The Bonferroni and Scheffé methods are used for general tests of possible contrasts. The Bonferroni method is better when the number of contrasts being tested is about the same as the number of factor levels. The Scheffé method covers all possible contrasts, and as a result, is the most conservative of all the methods. The drawback for such a highly conservative test, however, is that it becomes more difficult to resolve differences among means, even though the ANOVA would indicate that they exist.

When treatment levels include a control and mean comparisons are restricted to only comparing treatment levels against a control level, Dunnett’s mean comparison method is appropriate. Because there are fewer comparisons made in this case, the test provides more power compared to a test (see Section 3.7) using the full set of all pairwise comparisons.

To illustrate these methods, the following output was obtained for the hypothetical greenhouse data of our example. We will be running these types of analyses later in the course.

Since the estimated means for F1 and F3 are covered by the same colored bar, they are not significantly different using the LSD approach.

Since the estimated means for F1 and F3 are covered by the same colored bar (red bar), they are not significantly different using Tukey's approach. Similarly, since F1 and F2 are covered by the same colored bar (blue bar) they are not significantly different using Tukey's approach.

Results from the Bonferroni approach are similar to the ones from Tukey's approach.

Results from the Scheffé approach are similar to the ones from Tukey's and Bonferroni's approaches.

We can see that the LSD method was the most liberal, that is, it indicated the largest number of significant differences between means. In this example, Tukey, Bonferroni, and Scheffé produced the same results. The Dunnett test was consistent with the other 4 methods, and this is not surprising given the small value of the control mean compared to the other treatment levels.

To get a closer look at the results of employing the different methods, we can focus on the differences between the means for each possible pair:

and compare the 95% confidence intervals produced:

You can see that the LSD produced the narrowest confidence intervals for the differences between means. Dunnett’s test had the next most narrow intervals (but only compares treatment levels to the control). The Tukey method produced intervals that were similar to those obtained for the LSD, and the Scheffé method produced the broadest confidence intervals.

What does this mean? When we need to be REALLY sure about our results, we should use conservative tests. If you are working in life-and-death situations such as in most clinical trials or bridge building you might want to be surer. If the consequences are less severe you can use a more liberal test, understanding there is more of a chance you might be incorrect (but still able to detect differences). In reality, you need to be consistent with the rigor used in your discipline. While we can't tell you which comparison to use, we can tell you the differences among the tests and the trade-offs for each one.

Unsurprisingly, paired comparison methods (presented in Sections 2.2 and 2.3) are limited to comparisons made only between treatment mean pairs. However, a contrast analysis procedure can be used to carry out comparisons of a much wider context, such as comparisons of treatment level groups or even testing of trends. In the context of a single-factor ANOVA model, a linear contrast can be defined as a linear combination of the treatment means such that their numerical coefficients add to zero. Mathematically, a contrast can be represented by

\(A=\sum_{i=1}^{T} a_{i} \bar{y}_{i}\)

where \(\bar{y}_{1}, \bar{y}_{2}, \ldots, \bar{y}_{T}\) represent the sample treatment means and \(\sum_{i=1}^{T} a_{i}=0\). The quantity A is a sample statistic and serves as an estimate for the parameter contrast \(\sum_{i=1}^{T} a_{i} \mu_{i}\). By choosing the numerical coefficients appropriately, linear contrasts can be used to make different comparisons among groups of treatment means, including but not limited to mean pairs. The table below gives 4 linear contrasts defined in terms of the 3 fertilizer levels, F1, F2, F3, and the Control in the greenhouse example.

Notice that values of each list of \(a_{i}\) (i = 1, 2, 3, 4) add to zero. The first contrast compares the first two fertilizer types in terms of their means (a pairwise comparison), the second compares the means of the 3 fertilizer types with the Control mean. The 3rd one is a comparison between the combined effect of fertilizer types 1 and 2 with fertilizer type 3, while the last contrast compares the 2nd and 3rd fertilizer types.

A pair of contrasts \(A=\sum_{i=1}^{T} a_{i} \bar{y}_{i}\), and \(B=\sum_{i=1}^{T} b_{i} \bar{y}_{i}\) is orthogonal if the products of their numerical coefficients add to zero. This can be expressed mathematically as

\(\sum_{i=1}^{T} a_{i}b_{i}=0\)

A set of contrasts is said to be orthogonal if every pair of contrasts in the set is orthogonal. Two orthogonal contrasts are not correlated. This means that if A and B are orthogonal, then the Covariance (A, B) = 0.  Furthermore, the sum of squares of the treatment displayed usually in the ANOVA table, can be partitioned into a set of (T-1) orthogonal contrasts each with 1 degree of freedom. Note that the maximal number of orthogonal contrasts associated with the treatment of T levels is (T-1) and each of them would be associated with one specific comparison independent of the other. In the table above, contrasts 1, 2, and 3 form an orthogonal set of (T-1) contrasts. 

The statistical significance of a linear contrast, which can be equated to testing for the zero contrast value, can be formulated using the null and alternative hypotheses:

\(H_0\colon \sum_{i=1}^{T} a_{i} \mu_{i}=0 \text { vs. } H_A\colon  \sum_{i=1}^{T} a_{i} \mu_{i} \neq 0 \text {, }\)

and can be tested using either,

\(t=\dfrac{\sum_{i=1}^{T} a_{i} \bar{y}_{i}}{\sqrt{\operatorname{MSE} \sum_{i=1}^{T} \frac{a_{i}^{2}}{n_i}}}\) with (N-T) degrees of freedom or \(F=\dfrac{\left(\sum_{i=1}^{T} a_{i} \bar{y}_{i}\right)^{2}}{\operatorname{MSE} \sum_{i=1}^{T} \frac{a_{i}^{2}}{n_i}}\)

with the numerator and denominator degrees of freedom equal to 1 and (N-T) respectively.

Note that MSE can be obtained from the ANOVA table. Applying the above formula, the t- statistic for testing contrast 2 above is

\(t=\dfrac{\sum_{i=1}^{T} a_{i} \bar{y}_{i}}{\sqrt{\operatorname{MSE} \sum_{i=1}^{T}\frac{a_{i}^{2}}{n_i}}}=\dfrac{28.6+25.867+29.2-(3 * 21)}{\sqrt{3.052 \times \frac{(1+1+1+9)}{6}}}=8.365\)

with df=20 and has a p-value of approximately 0. This indicates that the average plant height due to the combined treatment of the 3 fertilizer types differs significantly from the average plant height yielded by the control.

The above testing procedure is applicable to non-orthogonal contrasts as well. However, as non-orthogonal contrasts are not guaranteed to be uncorrelated, the conclusions arrived at may be "overlapping", leading to redundancies. In Lesson 3, examples are provided to illustrate how software can be used to conduct contrast testing. The hypothesis testing for trends using contrasts will be discussed in Lesson 10 ANCOVA II.

In this lesson, we became familiar with the ANOVA methodology to test for equality among treatment means. As follow-up procedures, we were also exposed to the Tukey method for paired mean comparisons which helped to identify significantly different treatment (factor) levels. The contrast analysis was also discussed as a means to compare differences among group means.