Two-way ANOVA

In this exercise, we are going to fit a two-way ANOVA using the data set ‘soil’ in the provided soil.csv file. This data set contains soil phosphorus concentrations (\(mg\ kg^{-1}\)) derived from two different soil types: an alluvial soil (soil in proximity to a river experiencing regular flooding) and soil from an upland terrace (no soil nutrient replenishment through flooding). Soil samples were taken from the A-, B- and C-horizon. In statistical terms we have a continuous response variable and two categorical predictors (soil type ‘loc’ and soil horizon ‘horizon’). Because we have two predictor variables, the resulting model is called a two-way ANOVA.

  1. Fit a two-way ANOVA model using the explanatory variables ‘loc’ and ‘horizon’ including their interaction.
  2. Inspect the model diagnostic plots and comment on whether the variance homogeneity and the normality assumption are met.
  3. If there is a significant interaction term, follow up with a post-ho comparison, i.e. slice the interaction using the functions in the emmeans package.
  4. Report the overall significance of the interaction term as you would in a paper and produce a publication-quality figure with lower case letters indicating the results of the post-hoc test.
## First: Set your working directory (in the menu bar click on 'Session' => 'Set Working Directory' => 'Choose Directory...'. Copy working directory from the console and paste it into the RMarkdown script.)

## Read in the data set
soilp <- read.csv("soil.csv")

## Check structure and get data summary
head(soilp)
##        loc horizon      p depth.cm
## 1 alluvial       A 345.45     0-10
## 2 alluvial       A 220.64     0-10
## 3 alluvial       A 224.98     0-10
## 4 alluvial       A 246.72     0-10
## 5 alluvial       A 286.60     0-10
## 6 alluvial       A 283.31     0-10
str(soilp)
## 'data.frame':    54 obs. of  4 variables:
##  $ loc     : Factor w/ 2 levels "alluvial","terrace": 1 1 1 1 1 1 1 1 1 1 ...
##  $ horizon : Factor w/ 3 levels "A","B","C": 1 1 1 1 1 1 1 1 1 2 ...
##  $ p       : num  345 221 225 247 287 ...
##  $ depth.cm: Factor w/ 3 levels "0-10","20-30",..: 1 1 1 1 1 1 1 1 1 3 ...
summary(soilp)
##        loc     horizon       p            depth.cm 
##  alluvial:27   A:18    Min.   : 36.72   0-10  :18  
##  terrace :27   B:18    1st Qu.:117.07   20-30 :18  
##                C:18    Median :156.34   Oct-20:18  
##                        Mean   :164.05              
##                        3rd Qu.:201.25              
##                        Max.   :345.45