9
$\begingroup$

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?

$\endgroup$
  • $\begingroup$ This sounds like a Helmert contrast (c is the first level, then b, then a). $\endgroup$ – Michael M Aug 25 '15 at 8:03
  • $\begingroup$ I tried Helmert too, but the numbers are not the means I'm looking for. I've edited the question to include Helmert contrasts, thanks. $\endgroup$ – M4RT1NK4 Aug 25 '15 at 8:11
12
$\begingroup$

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.

$\endgroup$
  • $\begingroup$ Great, thanks! Just to make sure I understand this now - if I wanted to compare the first level to the rest of the levels in a 4 level variable, mat would be c(1, -1/3, -1/3, -1/3)? So I always set the numbers as they would be in the formula (a + (b+c+d)/3) and then ginv scales it appropriately so that the coefficients directly reflect the difference. And when you changed the order in the Helmert example, that was just to match the question? Otherwise, the results should be the same, regardless of the order of contrasts, right? $\endgroup$ – M4RT1NK4 Aug 27 '15 at 8:59
  • $\begingroup$ @M4RT1NK4 Your formula and the corresponding contrast ist correct. The order of columns was just changed to match the order of columns in the question. The order of rows, however, is important, because the first level is the reference level. In your example, the reference level is the third level. $\endgroup$ – Sven Hohenstein Aug 27 '15 at 10:33
  • $\begingroup$ @SvenHohenstein I had a related question based on this answer, mind having a look? stats.stackexchange.com/questions/357781/… $\endgroup$ – mat Jul 20 '18 at 9:57

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.