Using R to calculate between- by within-subj ANOVA interaction contrasts using car or lme

Question

Psychology student here. As part of my thesis I have some data from an experiment with two between-subjects IVs, classification and condition, and a within-subjects IV, trial_type. The data are of the form:

test = data.frame(
  ID=rep(c(1,2,3,4,5,6), each=3), 
  condition=rep(c('comp', 'seq', 'comp_text'), each=3, times=2),
  classification=rep(c('rule', 'exemplar'), each=3, length.out=18),
  trial_type=rep(c('ambig', 'unambig', 'trained'), length.out=18),
  value = c(0.25, 0.75, 1.00, 1.00, 1.00, 1.00, 0.00, 1.00, 0.75, 1.00, 1.00, 1.00, 0.25, 0.75, 0.75, 0.25, 0.75, 0.50))

I've analysed these as an anova using nlme, and generating type 3 SS, as this is what my course requires

model <- lme(value~condition*classification*trial_type, 
   random=~1|ID, 
   correlation = corCompSymm(form = ~1|ID)
)
anova(model, type='marginal')

I'd like to investigate an interaction present in my data. Specifically I want to see if the values for 'trained' trials are higher in the 'seq' condition than the other two. In the past I've used contrast() from the 'contrast' package to generate a contrast matrix, which I would then send to glht() from the multcomp package. In this instance I've tried:

cntr1 <- contrast(model, 
   a=list(condition=c('seq'), classification=c('rule','exemplar'), trial_type=c('trained')), 
   b=list(condition=c('comp','comp_text'), classification=c('rule', 'exemplar'), trial_type=c('trained'))
)

However, this returns

Error in testStatistic(fit, X, modelCoef, covMat, conf.int = conf.int) : 
Non-positive definite approximate variance-covariance

This appears to be something to do with the variance-covariance matrix, but I'm not exactly sure what the issue is.

My question is, is there a way to resolve this issue? If not, is there a way to run a similar contrast in another package, for instance, 'car'? I can generate the same omnibus F-test using Anova from the car package by converting the data so that each trial type has its own column, but I'm not sure how to test this sort of interaction using Anova().

Answer 1

Does this help you (let me know if I misunderstood you somewhere):

1. There were too few combinations in your provided test to perform the mixed ANOVA. There is only exactly one instance of each of all possible combinations. For illustration purposes I adjusted this:

   test = data.frame( ID=rep(c(1:12), each=3), condition=rep(c('comp', 'seq', 'comp_text'), each=6, times=2), classification=rep(c('rule', 'exemplar'), each=6, length.out=18), trial_type=rep(c('ambig', 'unambig', 'trained'), length.out=18), value = round(rnorm(36),2))

2. I would use ezANOVA for this, which provides all the interactions you want:

First re-arrange the data more suitable for ezANOVA:

ambig <- test[test$trial_type == "ambig", 'value']
unambig <- test[test$trial_type == "unambig", 'value']
trained <- test[test$trial_type == "trained", 'value']
df1 <- data.frame(ambig, unambig, trained)
df1 <- stack(df1)
subj <-  as.numeric(rep(unique(test$ID), 3))
condition <- rep(unique(test$condition), 3)
classification <- rep(unique(test$classification), 3)
df1[3] <-  subj
df1[4] <-  condition
df1[5] <- classification
colnames(df1) <- c("value", "trial_type", "id", "condition", "classification")

Then perform your the 3 (Trial_type; within-subjects) by 2 (Condition; between-subjects) by 2 (Classification; between-subjects) mixed ANOVA on value.

library(ez)
mixed_aov <- ezANOVA(
  data = df1
  , dv = value
  , wid = id
  , within = trial_type
  , within_covariates = NULL
  , between = .(classification, condition)
  , observed = NULL
  , diff = NULL
  , reverse_diff = FALSE
  , type = 3
  , white.adjust = FALSE
  , detailed = F
  , return_aov = T
)

See for details on the arguments of ezANOVA() see the documentation here.

The output gives you all your interactions.

For contrasts you would then use follow-up t-tests.

P.S. Also watch for your typo in the lme specification (trial_type) .