Specifying contrasts for lme in 3x3 repeated-measures design

Question

I would like to test the effects of two categorical variables, each with 3 levels, on some continuous data. Three groups of participants completed the same task for 3 types of stimuli, so it's a 3x3 repeated measures design. I'm using lme. The fixed-effects part of the model is the interaction between participant group ("play") and stimuli ("ins"), and the random-effects part is subject grouped within stimulus category.

contrasts(iqr$ins) <- contr.sum(3)
contrasts(iqr$play) <- contr.sum(3)
iqr.lme <- lme(iqr ~ play*ins, random = ~1|subject/ins, data = iqr)

I'm using sum contrasts and would like to see comparisons between each level of the "play" and "ins" variables. In other words, I would like to use a contrast matrix for each variable that looks like the one below. If I understand correctly, -1 is the reference level and 1 is the level compared to it. Specifying contr.sum only results in contrasts for the first two columns being given, however.

          [,1]  [,2]  [,3]
clarinet    1     0     1
piano       0     1    -1
violin     -1    -1     0

Here is the fixed effects part of the summary output:

Fixed effects: iqr ~ play + ins + play * ins 
                 Value  Std.Error DF   t-value p-value
(Intercept)  2.9768519 0.11599299 66 25.664066  0.0000
play1       -0.0185185 0.16403886 33 -0.112891  0.9108
play2        0.2037037 0.16403886 33  1.241801  0.2231
ins1         0.1342593 0.09005357 66  1.490882  0.1408
ins2        -0.3240741 0.09005357 66 -3.598681  0.0006
play1:ins1  -0.3009259 0.12735498 66 -2.362891  0.0211
play2:ins1   0.0601852 0.12735498 66  0.472578  0.6381
play1:ins2   0.1574074 0.12735498 66  1.235974  0.2208
play2:ins2  -0.1481481 0.12735498 66 -1.163269  0.2489

If I understand correctly, ins1 is the contrast for clarinet vs. violin (using the contrast matrix above), and ins2 is the contrast for piano vs. violin. The problem is that I would also like to see a contrast for clarinet vs. piano.

My solution was to run a second model using a contrast matrix that I specified myself. I thought this would give me the missing contrast and replicate the clarinet vs. violin values from the original output.

cst1 <- cbind(c(-1,1,0), c(-1,0,1))
contrasts(iqr$ins) <- cst1
contrasts(iqr$play) <- cst1
iqr.lme <- lme(iqr ~ play*ins, random = ~1|subject/ins, data = iqr)
contrasts(iqr$ins)

          [,1]  [,2]
clarinet   -1    -1
piano       1     0
violin      0     1

Here, ins1 should correspond to clarinet vs. piano, right? But the values are the same as the ins2 contrast in my original output, indicating that the contrast matrix I specified was blatantly ignored and ins1 instead corresponds to piano vs. violin. The values for play2 and ins2 are new, and I could just assume that they indicate the contrast that was missing from the original output, but I would really rather understand what is going on than base my results on dubious assumptions!

Fixed effects: iqr ~ play * ins 
                 Value  Std.Error DF   t-value p-value
(Intercept)  2.9768519 0.11599299 66 25.664066  0.0000
play1        0.2037037 0.16403886 33  1.241801  0.2231
play2       -0.1851852 0.16403886 33 -1.128910  0.2671
ins1        -0.3240741 0.09005357 66 -3.598681  0.0006
ins2         0.1898148 0.09005357 66  2.107799  0.0389
play1:ins1  -0.1481481 0.12735498 66 -1.163269  0.2489
play2:ins1  -0.0092593 0.12735498 66 -0.072704  0.9423
play1:ins2   0.0879630 0.12735498 66  0.690691  0.4922
play2:ins2  -0.2314815 0.12735498 66 -1.817608  0.0737

Am I misunderstanding how contrast matrices are read or specified? Or is there actually a way to get all the contrasts I want listed at once? I am fairly new to R and very new to nlme, so I could very well be missing something obvious. I hope this is enough information for someone to be able to help. Thanks a lot!

Answer 1

If you have a factor with $n$ levels, you can test $n - 1$ effects (at maximum) associated with this factor in one model.

Have a look at the sum contrasts. The function contr.sum(3) returns the following matrix:

  [,1] [,2]
1    1    0
2    0    1
3   -1   -1

The columns denote the contrasts. For a factor with $3$ levels, you obtain $2$ contrasts. Your interpretation of the contrasts is not correct. The first one tests the first level of the factor against levels 2 and 3. The second contrasts tests the second levels against levels 1 and 3.

You manually specified a contrast matrix with cbind(c(-1,1,0), c(-1,0,1)):

     [,1] [,2]
[1,]   -1   -1
[2,]    1    0
[3,]    0    1

This matrix also allows two tests. You test level 2 against levels 1 and 3. This is identical to the second test in the first matrix. Hence, you obtain the identical coefficient in both models. The second test of this contrast matrix is level 3 against levels 1 and 2.

If you want to compare levels against each other, there are two possibilities:

• Treatment contrasts:

  contr.treatment(3)
  
    2 3
  1 0 0
  2 1 0
  3 0 1
  

  Here, you test the first level against (a) the second and (b) the third one.

• Repeated contrasts (sliding differences):

  library(MASS)
  contr.sdif(3)
  
           2-1        3-2
  1 -0.6666667 -0.3333333
  2  0.3333333 -0.3333333
  3  0.3333333  0.6666667
  

  Here, you test (a) the first level against the second one and (b) the second level against the third one. This matches the column names.

Please note that no contrast matrix will allow you to have more than $2$ contrasts for a factor with $3$ levels.

If you want to manually create an appropriate contrast matrix based on specific tests, you can use the ginv function from the MASS package:

cst1 <- cbind(c(-1,1,0), c(-1,0,1)) # your contrast matrix

library(MASS)
t(ginv(cst1))

       [,1]       [,2]
[1,] -0.3333333 -0.3333333
[2,]  0.6666667 -0.3333333
[3,] -0.3333333  0.6666667

The resulting matrix allows the tests you are looking for, i.e., level 2 against level 1 and level 3 against level 1.

For more details on how you can create a contrast matrix based on manually defined comparisons, have a look at this answer.