Coventry University Logo Sigma Logo

Two-Way analysis of variance (ANOVA) in R (Using R Studio)

When to use a two-way ANOVA

A two-way analysis of variance (ANOVA) is a statistical test that can be used to assess whether the mean value of some outcome depends on two other factors.

Example: we might want to assess if weight depends on nationality and gender.

  • Here weight is our Dependent variable
  • Nationality and gender are our two Factors (Independent variables) that we think weight may depend on.

We might also be interested in whether the way that age affects maths scores is the same across all nationalities.

The observations/measurements on the Dependent variable, such as weight, must be Scale data.

For more help on “What test do I need” go to the sigma website statistical worksheets resources page.

Example

An educational researcher collected data on maths scores from adults of different ages living in different countries. They were grouped into nationality according to whether they were from the UK, Europe or the rest of the World. Their age was classified according to whether they fell into three pre-defined age groups: Younger, Middle-aged or Older. The data are shown below and can be downloaded in a CSV file called maths.csv, where Nationality was coded 1=UK, 2=Europe and 3=Rest of World, and Age was coded as 1=Younger, 2=Middle-aged and 3=Older.

A two-way Analysis of Variance (ANOVA) can be used here to test three different sets of hypotheses:

Test 1:

H0: There is no interaction between nationality and age

H1: There is an interaction between nationality and age

The presence of an interaction between nationality and age means that any age effects observed are not the same for all nationalities. Or equally that any differences due to nationality are not the same in every age group.

Test 2:

H0: There is no difference in mean maths score between the three nationality groups

H1: There is a difference in mean maths score between the three nationality groups

Test 3:

H0: There is no difference in mean maths score between the three age groups

H1: There is a difference in mean maths score between the three age groups

Note that in Tests 2 and 3 we assume that any age effects observed are the same for all nationalities. Or equally that any differences due to nationality are the same in every age group.

Hence, we should do Test 1 first since if conclude there is an interaction between nationality and age then we cannot do Tests 2 and 3. If there is an interaction we should skip tests 2 and 3 and look at the mean maths score for each nationality/age combination and report them (we will see how to plot these means in a profile plot).

If we conclude there is no interaction between nationality and age, then we can do Tests 2 and 3.

To get started with the analysis, first, bring the dataset into RStudio. To do this you can either run a read.csv() function if you know how to do this or alternatively you can follow these steps using the menus:

From the File menu select Import Dataset then From Text(base):

From the pop-up window navigate to the folder where you have saved the dataset, then once the file was selected click Open:

At the next dialogue box (see below), in the upper left corner in the “Name” field, amend the name of your dataset if you wish, in this example we named it as “maths”.

We should also make sure the Heading option below is set to Yes (otherwise the data will all be read in as text):

Finally, click on Import to complete the process. This imports the data set and is listed in the Environment in the top right of your RStudio screen as follows:

Using R

We will carry out the two way anova using the function aov(). First, we need to convert Nationality and Age.

  • Note: We indicate which variables we want to use by typing type the name of our dataset followed by a dollar sign “$” and the name of our column that contains the data.

Then, in the aov() function enter the dependent variable, followed by the independent.

maths$Nationality <- as.factor(maths$Nationality)
maths$Age <- as.factor(maths$Age)
anova <- aov(Maths_score ~ Nationality * Age, data=maths, type=III)
summary(anova)

Examining the Output

##                 Df Sum Sq Mean Sq F value   Pr(>F)    
## Nationality      2   2274  1137.0  11.255 0.000109 ***
## Age              2   1748   874.1   8.653 0.000661 ***
## Nationality:Age  4    320    79.9   0.791 0.537229    
## Residuals       45   4546   101.0                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

First, we undertake Test 1:

H0: There is no interaction between nationality and age

H1: There is an interaction between nationality and age

We examine the p-value in the row labelled Nationality*Age in the table, which assesses the interaction. In our case we have p=0.537. Since this is greater than 0.05, we conclude there is no evidence of an interaction between nationality and age. This means that any age effects observed can be assumed to be the same for all nationalities. Or equally that any differences due to nationality can be assumed to be the same in every age group.

Hence, we can proceed to undertake Tests 2 and 3:

Test 2 examines the nationality effect. To do this we examine the p-value in the row labelled Nationality in the table. Here we have p reported to be less than 0.001. Since this is less than 0.05, we conclude there is evidence of a nationality effect.

Test 3 examines the age effect. Examining the p-value in the row labelled Age we see this is also reported to be less than 0.001 and so we conclude there is evidence of an age effect.

Having established there are differences in maths marks due to both nationality and age, we need to undertake post hoc tests to understand where these differences are.

Post-hoc Tests

To obtain the post hoc tests, we apply the function TukeyHSD() to the ANOVA model previously fit:

TukeyHSD(anova)
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Maths_score ~ Nationality * Age, data = maths, type = III)
## 
## $Nationality
##          diff       lwr      upr     p adj
## 2-1  2.944444 -5.175603 11.06449 0.6562918
## 3-1 15.000000  6.879953 23.12005 0.0001494
## 3-2 12.055556  3.935509 20.17560 0.0022479
## 
## $Age
##          diff        lwr      upr     p adj
## 2-1  5.444444 -2.6756026 13.56449 0.2457384
## 3-1 13.833333  5.7132863 21.95338 0.0004503
## 3-2  8.388889  0.2688419 16.50894 0.0415005
## 
## $`Nationality:Age`
##               diff        lwr      upr     p adj
## 2:1-1:1  4.5000000 -14.401297 23.40130 0.9970285
## 3:1-1:1 17.3333333  -1.567964 36.23463 0.0954064
## 1:2-1:1 10.8333333  -8.067964 29.73463 0.6388973
## 2:2-1:1  9.6666667  -9.234630 28.56796 0.7632425
## 3:2-1:1 17.6666667  -1.234630 36.56796 0.0835706
## 1:3-1:1 12.3333333  -6.567964 31.23463 0.4698826
## 2:3-1:1 17.8333333  -1.067964 36.73463 0.0781349
## 3:3-1:1 33.1666667  14.265370 52.06796 0.0000275
## 3:1-2:1 12.8333333  -6.067964 31.73463 0.4160744
## 1:2-2:1  6.3333333 -12.567964 25.23463 0.9725853
## 2:2-2:1  5.1666667 -13.734630 24.06796 0.9924275
## 3:2-2:1 13.1666667  -5.734630 32.06796 0.3817025
## 1:3-2:1  7.8333333 -11.067964 26.73463 0.9101188
## 2:3-2:1 13.3333333  -5.567964 32.23463 0.3650439
## 3:3-2:1 28.6666667   9.765370 47.56796 0.0003571
## 1:2-3:1 -6.5000000 -25.401297 12.40130 0.9680020
## 2:2-3:1 -7.6666667 -26.567964 11.23463 0.9196002
## 3:2-3:1  0.3333333 -18.567964 19.23463 1.0000000
## 1:3-3:1 -5.0000000 -23.901297 13.90130 0.9939110
## 2:3-3:1  0.5000000 -18.401297 19.40130 1.0000000
## 3:3-3:1 15.8333333  -3.067964 34.73463 0.1669869
## 2:2-1:2 -1.1666667 -20.067964 17.73463 0.9999999
## 3:2-1:2  6.8333333 -12.067964 25.73463 0.9572087
## 1:3-1:2  1.5000000 -17.401297 20.40130 0.9999993
## 2:3-1:2  7.0000000 -11.901297 25.90130 0.9509443
## 3:3-1:2 22.3333333   3.432036 41.23463 0.0102053
## 3:2-2:2  8.0000000 -10.901297 26.90130 0.8999631
## 1:3-2:2  2.6666667 -16.234630 21.56796 0.9999354
## 2:3-2:2  8.1666667 -10.734630 27.06796 0.8891325
## 3:3-2:2 23.5000000   4.598703 42.40130 0.0056916
## 1:3-3:2 -5.3333333 -24.234630 13.56796 0.9906720
## 2:3-3:2  0.1666667 -18.734630 19.06796 1.0000000
## 3:3-3:2 15.5000000  -3.401297 34.40130 0.1874808
## 2:3-1:3  5.5000000 -13.401297 24.40130 0.9886122
## 3:3-1:3 20.8333333   1.932036 39.73463 0.0209771
## 3:3-2:3 15.3333333  -3.567964 34.23463 0.1984041

Looking down the p adj column in the Nationality table, we can see there are some statistically significant differences, where p is less than 0.05. There is a statistically significant difference between the UK (coded 1) and the World (coded 3), p reported as less than 0.001, and between Europe (coded 2) and the World (coded 3), p=0.002. But there is no evidence of a statistically significant difference between the UK (coded 1) and Europe (coded 2), p=0.656. Hence the UK and Europe are similar, but both are different to the rest of the world.

This time looking down the p adj column in the Age table, we can see there are again some statistically significant differences. There is a statistically significant difference between younger (coded 1) and older people (coded 3), p reported as less than 0.001, and between middle-aged (coded 2) and older people (coded 3), p=0.042. But there is no evidence of a statistically significant difference between the younger (coded 1) and middle-aged people (coded 2), p=0.246. Hence it seems younger and middle-aged people have similar ability in maths, but these are both different to older people.

We also want marginal mean tables for nationality, age, and their interation. This can be done using the aggreate() function.

aggregate(Maths_score ~ Nationality, data=maths, mean)
aggregate(Maths_score ~ Age, data=maths, mean)
aggregate(Maths_score ~ Nationality + Age, data=maths, mean)

In the output we have the means for each nationality and age group, as well as their combinations. For example, the mean mark for those in the UK (coded 1) that are in the younger age group (coded 1) was 50.167.

We can also create a profile plot using the results previously obtained.

means <- aggregate(Maths_score ~ Nationality + Age, data=maths, mean)

interaction.plot(means$Age, means$Nationality, means$Maths_score,
                 xlab="Age", ylab="Maths Score")

Note: To generate plots in R, we need to indicate the data or variables we want to plot. We can do that by indicating which columns in our dataset contain that data, for the case of Rstudio, we can type the name of our dataset followed by a dollar sign “$” and the name of our column that contains the data as in the snippet of code above. Additionally, to add labels, we can add them by using the xlab and ylab options inside the function as seen in the code provided.

Reporting Results

We could report the results as:

“A two-way ANOVA was used to assess the effect of nationality and age on maths scores. There was no evidence of an interaction effect between nationality and age (p=0.537), but there was evidence of differences due to both nationality (p<0.001) and age (p<0.001). There was a statistically significant difference between the UK (M=57.89), p<0.001, and the World (M=72.89) and between Europe (M=60.83) and the World, p=0.002, but no evidence of difference between the UK and Europe (p=0.656). There was also a statistically significant difference between younger people (M=57.44) and older people (M=71.28), p<0.001, and between middle-aged people (M=62.89) and older people, p=0.042. But there is no evidence of a statistically significant difference between the younger and middle-aged people (p=0.246).”

Further Work

Checking assumptions: Normality

For the test to be valid it should be reasonable to assume that the measurements in each age group and each nationality are approximately normally distributed. If we have 30 or more measurements in each nationality or each age group, then we can safely make that assumption, and we need not check this any further. Since we only had 18 measurements in each nationality group or 18 in each age group, we need to do some further assessments.

The assumption of normality is based on what are called the residuals. The residuals are the differences between each observed value and the mean value for that nationality/age sub-group. We therefore just need to check if this single set of residuals can be assumed to be normally distributed.

Normality could be judged by examining the shape of the histogram to see if it makes a roughly symmetric bell-shaped curve.

To obtain the residuals from our two-way ANOVA model, we can use the residuals() function. We can then plot these using the hist() function.

residuals <- residuals(anova)
hist(residuals)

We can also assess normality using the Shapiro-Wilk test, which can be done using shapiro.test() function. We need a non-significant result, i.e. the sig value (p-value) needs to be greater than 0.05 to be able to assume normality. In our example, p = 0.861 and so we can assume normality.

shapiro.test(residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals
## W = 0.98795, p-value = 0.8614

Checking assumptions: Equal variances

The ANOVA test we undertook also assumes homogeneity (equality) of variance. This essentially means, can we assume the amount by which maths marks vary within each age group is similar, and also the variation in maths marks within each nationality group is similar.

Levene’s test is used in R to evaluate the homogeneity of variance assumption.

First, install the package car. Then, apply the library() function to this. This will give us access to the function LeveneTest() which we will use next.

Note: to add additional functions to Rstudio we can install external packages. These can add more capabilities not present in the default Rstudio or when there is no other way to carry on an analysis. After installing the package we need to call it by using library(). This will indicate Rstudio to use the functions inside this package to run calculations or procedures requested by the user.

install.packages("car", repos="https://cloud.r-project.org")  
library(car)

We can then use LeveneTest().

leveneTest(Maths_score ~ Nationality * Age, data = maths)

We examine the p-value. You need the p-value to be greater than 0.05 to be able to assume homogeneity of variances. In our case p=0.8105 so this assumption is fine in our case.

For more resources, see sigma.coventry.ac.uk Adapted from material developed by Coventry University Creative Commons License