How to set custom contrasts with lmer in R

Question

I am using lmer in R to check the effect of condition (cond) on some result. Here are some made up data, where s is the subject identifier and a, b and c are conditions.

library("tidyr")
library("dplyr")
set.seed(123)
temp <- data.frame(s = paste0("S", 1:30), 
                   a = rnorm(30, -2, 1), 
                   b = rnorm(30, -3, 1), 
                   c = rnorm(30, -4, 1)) 

I would like to compare

1. level a to the mean of levels b and c and
2. level b to level c.

My question is, how do I set the contrasts to do this in such a way that the intercept reflects the mean of the three conditions and the two computed estimates directly reflect differences as defined in 1. and 2.?

I tried with

c1 <- cbind(c(-0.5, 0.25, 0.25), c(0, -0.5, 0.5))
gather(temp, cond, result, a, b, c) %>%
  lmer(result ~ cond + (1|s), data = ., contrasts = list(cond = c1))

where cond2 seems to be OK, but cond1 is not.

Following How to interpret these custom contrasts?, I tried to use the generalized inverse instead, but these estimates don't make sense either.

c2 <- t(ginv(c1))
gather(temp, cond, result, a, b, c) %>%
  lmer(result ~ cond + (1|s), data = ., contrasts = list(cond = c2))

I tried Helmert contrasts too, but the means still don't match up.

gather(temp, cond, result, a, b, c) %>%
  mutate(cond = factor(cond, levels = c("c", "b", "a"))) %>%
  lmer(result ~ cond + (1|s), data = ., contrasts = list(cond = contr.helmert))

What is the correct way to do this?

Answer 1

For the following steps, we need the data frame in the long format. The data frame dat contains the dependent variable result, the categorical predictor cond (levels: a, b, and c), and the random factor s.

library(tidyr)
dat <- gather(temp, cond, result, a, b, c)

In the following, I will illustrate two approaches to create a contrast matrix corresponding to the conditions you want to compare:

1. $a - \frac{b+c}{2}$
2. $b - c$

Custom contrasts

The matrix mat corresponds to the level differences.

mat <- rbind(c(1, -0.5, -0.5),     # a vs. (b + c) / 2
             c(0, 1, -1))          # b vs. c

To create the actual contrast matrix, we compute the generalized inverse with ginv (from MASS).

library(MASS)
cMat <- ginv(mat)
#            [,1]          [,2]
# [1,]  0.6666667 -7.130169e-17
# [2,] -0.3333333  5.000000e-01
# [3,] -0.3333333 -5.000000e-01

This contrast matrix cMat can be used in lmer.

library(lme4)
res <- lmer(result ~ cond + (1|s), data = dat, 
            contrasts = list(cond = cMat))
coef(summary(res))    
#              Estimate Std. Error    t value
# (Intercept) -2.948115  0.0946025 -31.163182
# cond1        1.351517  0.2006822   6.734612
# cond2        1.153918  0.2317279   4.979625

As you can see, the fixed-effect estimates correspond to the differences specified above. Furthermore, the intercept represents the overall mean.

Helmert contrast with contr.helmert

You can also use the built-in contr.helmert function to create the contrast matrix.

cHelmert <- contr.helmert(3)
#   [,1] [,2]
# 1   -1   -1
# 2    1   -1
# 3    0    2

However, the order does not correspond to the one you specified in the question. Hence, we have to reverse the order of columns and rows. The first column corresponds to b vs. a and the second one corresponds to c vs. the mean of b and a.

cHelmert2 <- cHelmert[c(3:1), 2:1]
#   [,1] [,2]
# 3    2    0
# 2   -1    1
# 1   -1   -1

Compare the contrast matrix cHelmert2 to cMat. You will notice that the columns are scaled versions of the other matrix.

The result of lmer is:

library(lme4)
res2 <- lmer(result ~ cond + (1|s), data = dat, 
             contrasts = list(cond = cHelmert2))
coef(summary(res2))    
#               Estimate Std. Error    t value
# (Intercept) -2.9481150 0.09460250 -31.163182
# cond1        0.4505056 0.06689407   6.734612
# cond2        0.5769590 0.11586393   4.979625

This contrast matrix allows fo the same comparisons as the custom contrast matrix. However, since the values in the matrix are different, the fixed-effects coefficients are different too. Not suprisingly, the $t$-values are the same.