31.3.1 Splitting variability

ANOVA does this by splitting up the total variability (how much the data differ on the whole, called SSTotal) into two components:

• The variability between groups (how much the group means differ, called SSG)
• The variability within groups (how much the values within each group differ, called SSE)

Here, the SS stands for sum of squares, and gives us a hint about how it is calculated. Those calculations are long, confusing, and tedious, so (in a moment) we will let R handle them. The G rather clearly stands for groups, but the E stands, perhaps confusingly at first, for “error” instead of something else. This is for historical reasons, and is related to ANOVA’s emergence from linear modelling, which we will cover in the next tutorial.

31.3.2 Scaling variability

However, as we learned with the t-statistic, when we estimate variability from a smaller group, we are likely to underestimate it. In the plots above, there are 30 data points in each group and only 3 groups. We need to find some way to standardize those scores before we compare them.

ANOVA does this by calulating something called the “Mean Square”, which is a measure of the “mean” variability. To scale these though, the sums of squares (SS) are divided by the degrees of freedom, rather than the total samples. This accounts for the risk of underestimation, and builds in some protections against over-interpretting poorly designed tests.

So, we now need to calculate degrees of freedom both for the groups and for the errors. The groups should be simple: they follow the same rule as chi-square, the number of groups minus one (often denoted as k - 1). We might think that we can now do the same for the errors (using n in place of k), but here we run into a problem.

The total degrees of freedom for the test should only be n - 1, so, we need to account for the fact that we are already using some of those for the groups. In ANOVA, this is accomplished by setting the degrees of freedom for the errors (within groups) to n - k, yielding the correct total degress of freedom. Don’t worry if that was a bit confusing, R will calculate these for you, you just need to know where they come from.

So, we can now scale each of our SS terms to calculate the mean sum of squares for groups (MSG) and the mean sum of squares for errors (MSE) using these formulas:

\[ MSG = \frac{SSG}{k-1} \]

\[ MSE = \frac{SSE}{n-k} \]

31.3.3 Our new test statistic

Now, as we calculated a z-statistic for proportions, a t-statistic for comparing one or two means, and a chi-square-statistic for comparing distributions, we now need a new (complicated) statistic to test differences in means between multiple groups. I mentioned that ANOVA is comparing the source of variability, so it should be clear that we are some how going to compare our (now appropriately scaled) MSE and MSG values. But … how?

Well, like the chi-square, the distribution we will use (the f-distribution) treats big numbers as less likely under the null (that the groups are all equal). So, we need to make sure that if MSG is bigger than MSE we get a bigger value. However, as we saw in panels A and C above, we need this test to work similarly no matter how big the values are. So, what we need is a ratio of the two sources of variability. This ratio is called the F-statistic, and is calculated as:

\[ F = \frac{MSG}{MSE} \]

So, the bigger the group differences are compared to the variability within groups, the bigger the F-statistic.

31.3.4 A distribution to compare to

Finally, we need to know whether or not our F-statistic is unusual or not. As always, we could simulate, but here (because the calculations are beyond the scope of this tutorial) we will jump straight to the approximation. The F-distribution has the same set of functions as all of the others we have encountered. The only main difference is that it requires two specificied degrees of freedom (one for the groups and one for the errors).

Just to see what it looks like (here, using 3, as if we had 4 groups (4 - 1 = 3); and 20, as if we had 24 samples (24 - 4 = 20) ).

# Example f- distribution
# If there were 4 groups and 24 samples
curve(df(x,
         df1 = 4 - 1 ,
         df2 = 24 - 4),
      from = 0, to = 5,
      add = FALSE, xpd = TRUE,
      col = "green3", lwd = 4,
      ylab = "Density")

You won’t need to worry much about using these functions directly, as we will leave the rest up to R.

31.3.5 Assumptions of ANOVA

As before, when we use these shortcuts, they come at the cost of making (and meeting) some assumptions. For ANOVA, the first (normality) is the same one we have encountered before. As always, these tests rely on the Central Limit Theorem, and so require that either the underlying data is roughly normal or that the sample size is large. So, we always need to plot our data before running ANOVA. Second, ANOVA assumes that the variances are equal in each group.

The assumption of equal variances comes from the history of the ANOVA as it was designed to test outcomes of experiments, where if the null (no effect of treatment) is true, then all of the groups will have the same variability. In practice, as long as the standard deviation of one group isn’t double that of another, we are generally ok. Plotting the data is usually enough to test this too.

31.4.1 Plot the result

Finally, we always want to be able to visualize our results. For this, just like when using the t-test, we want a comparison of means with some indication of their variability, such as a confidence interval. So, we also want ot use the plotmeans() function here. Just like when using t.test(), the argument to plotmeans() is the same as the argument to lm(), so we can just copy and paste it in, though (as before), I also like to turn off some of the plot details.

# Plot the differences
plotmeans(rLegInjury ~ size,
          data = crash,
          n.label = FALSE,
          connect = FALSE)

As expected, we see that the confidence intervals are substantially overlapping: just like we expected given our large p-value.

31.4.2 Try it out

Is the mean chest deceleration (column chestDecel, a measure of how quickly the body slows) the same in all sizes of cars?

State your null and alternative hypotheses, run and interpret a test, then visualize the result.

You should get a p value of 0.0000000074, and a plot like this:

31.5.1 Next steps

There are two general approaches to determining which group(s) actually differ. The first, which many introductory books recommend, is to manually make direct pairwise comparisons and see which ones are significant. Only doing this following a significant ANOVA result reduces the likelihood of Type I error (we already believe there is at least one real difference), but it doesn’t elminate it. Those introductory texts leave dealing with this (very major) concern, beyond the scope of the text, but give you just enough information to be dangerous.

Instead, I am going to give you just a little more information, but ask you to do less interpretation from it. There are a number of systems that handle the issue of multiple comparisions after an ANOVA, but one of the most common is called Tukey’s Honest Significant Differences. In R, this is implemented by passing a (slightly different) ANOVA object to the function TukeyHSD(). For this, we use that aov() function, instead of anova().

The TukeyHSD() function reports a confidence interval (95% by default, but can be changed with the conf.level = argument) and a p-value for each pairwise comparison of groups. These p-values have been adjusted for the multiple comparisons, so we can (more) safely interpret them. If the p adj is below your alpha threshold, you can conclude that the two groups are different. To see it in action:

# Run a Tukey "posthoc" test
TukeyHSD( aov(hotdogCalLM) )

##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = hotdogCalLM)
## 
## $Type
##                    diff       lwr       upr     p adj
## Meat-Beef      1.855882 -16.82550  20.53726 0.9688129
## Poultry-Beef -38.085294 -56.76667 -19.40391 0.0000277
## Poultry-Meat -39.941176 -59.36515 -20.51720 0.0000239

Here, we see that (as expected) the Poultry group differs from both Beef and Meat, but Beef and Meat do not differ from each other.

Finally, you will almost always want to plot something to show the results of a comparison like this.

# Plot the means and CI's
plotmeans(Calories ~ Type,
          data = hotdogs,
          n.label = FALSE,
          connect = FALSE )

We can also get the means from the anova object if we want to report them directly. Here, rep gives us the sample size in each group, and the values are the mean of each group (Grand mean is the overall mean).

# Get the means
model.tables( aov(hotdogCalLM), type = "means" )

## Tables of means
## Grand mean
##          
## 145.4444 
## 
##  Type 
##      Beef  Meat Poultry
##     156.8 158.7   118.8
## rep  20.0  17.0    17.0

31.5.2 Try it out

Is the mean measure of right leg injury (column rLegInjury) the same in various types of cars? Set up, run, interpret, and plot an appropriate test. Which pair-wise differences are significant?

You should get a p value of 0.0214, and a plot like this:

The comparison between pickups and compact cars is the only significant difference identified. This suggests that we may not have sufficient power to accurately identify any of the smaller differences between groups.