# How to perform a planned comparison test (in R, using gmodels)

## Task

Suppose that ANOVA reveals a significant difference between treatment levels, and you wish to explore further through post hoc analysis by comparing two specific treatment levels. How can we perform perform planned comparisons, also called a contrast test?

## Solution

Usually, you have data you wish to compare, but we will use example data here.
We load the “oats” dataset from R’s `MASS`

package, about the yield of oats
from a split-plot field trial using three varieties (V) and four levels of
manurial treatment (N). The experiment was laid out in 6 blocks (B) of 3 main
plots, each split into 4 sub-plots. The varieties were applied to the main plots
and the manurial treatments to the sub-plots.

1
2
3

# install.package('MASS') # if you have not already done so, and want this data
library(MASS)
df <- oats

Before we perform the contrast test, let’s verify that the yield of oats `Y`

depends on the nitrogen manurial treatment given to it `N`

.

1
2

aov1 <- aov(Y ~ N, data = df)
summary(aov1)

1
2
3
4
5

Df Sum Sq Mean Sq F value Pr(>F)
N 3 20020 6673 14.2 2.78e-07 ***
Residuals 68 31965 470
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

The $p$-value in the `Pr(>F)`

column is below $\alpha=0.05$.
So at the 5% significance level, the yield differs according to the
nitrogen manurial treatment. We assume that the model assumptions are met but
do not verify them here.

We now want to perform a planned comparison test (or contrast test) on the data
to see whether there is a difference between the $N<0.5$ levels and the
$N>0.5$ levels. We will use the `fit.contrast`

function in the `gmodels`

package.
Since the order of the levels is 0, 0.2, 0.4 and 0.6, the
contrast coefficients will be $-0.5$, $-0.5$, $0.5$, $0.5$, respectively.

1
2
3

# install.package('gmodels') # if you have not already done so
library(gmodels)
fit.contrast(aov1, "N", coeff=c(-1/2,-1/2,1/2,1/2))

1
2
3
4

Estimate Std. Error t value Pr(>|t|)
N c=( -0.5 -0.5 0.5 0.5 ) 29.66667 5.110338 5.805225 1.855598e-07
attr(,"class")
[1] "fit_contrast"

The $p$-value in the `Pr(>|t|)`

column is below $\alpha=0.05$.
This tells us that there is a significant difference between the average yields
of the $N<0.5$ and $N>0.5$ levels.

Content last modified on 24 July 2023.

See a problem? Tell us or edit the source.

Contributed by:

- Krtin Juneja (KJUNEJA@falcon.bentley.edu)
- Nathan Carter (ncarter@bentley.edu)