I hope this isn't a silly question, I'd like some advice on following up a threeway interaction in a mixed effects model. I've been building my models incrementally, like this:
lmer0<- lmer(predval ~ 1 + (1 + Quantile|Subj), data = d)
lmer1<- update(lmer0, .~. + epoch)
lmer2<- update(lmer1, .~. + Quantile)
lmer3<- update(lmer2., .~. + Congruence)
lmer4<- update(lmer3, .~. + FixDur)
lmer5<- update(lmer4, .~. + Quantile:Congruence)
lmer6<- update(lmer5, .~. + Quantile:FixDur)
lmer7<- update(lmer6, .~. + Congruence:FixDur)
lmer8<- update(lmer7, .~. + Quantile:Congruence:FixDur)
Here's what I get when I test these models against each other. All interaction terms improve the model (indicated by decreasing AIC/BIC values). This printout indicates the fully interactive model (lmer8) is the best model so far.
Df AIC BIC logLik Chisq Chi Df Pr(>Chisq)
lmer0 5 2062522 2062573 -1031256
lmer1 6 2053156 2053217 -1026572 9368.2240 1 < 2.2e-16 ***
lmer2 7 2053144 2053215 -1026565 14.2730 1 0.0001581 ***
lmer3 8 2051857 2051938 -1025921 1288.7862 1 < 2.2e-16 ***
lmer4 11 2051604 2051716 -1025791 259.3176 3 < 2.2e-16 ***
lmer5 12 2051606 2051728 -1025791 0.0284 1 0.8662848
lmer6 15 2051168 2051320 -1025569 443.8425 3 < 2.2e-16 ***
lmer7 18 2051152 2051335 -1025558 21.9526 3 6.673e-05 ***
lmer8 21 2050697 2050911 -1025327 460.9543 3 < 2.2e-16 ***
It's the last term in the model that really interests me - the three way interaction. Unfortunately I can't embed the image here (newbie), but you can find a plot of the data shows Congruence (shape) by Quantile (Colour gradient/x-axis) by Fixation Duration (facet_grid) here:
What I think I should do now is as follows:
- Test the Congruence x Quantile 2-way at each level of FixDur (i.e. Begin new model comparison at each level of FixDur, that's each column in the plot above).
- For those two-way interactions that are significant (they all are), follow with t-tests at each Quantile (i.e. 20 paired t-tests per FixDur). I'm not too sure how to conduct these using lmer, I'm basically just using the t.test function.
The obvious issue with this will be multiple comparisons. Essentially I want to do 20 t-tests at each level of Fixation Duration (so I can say something about the timecourse of my effects). That's a lot of t-tests, so using a Bonferroni or Scheffe correction is going to kill me!
Does anyone have a suggestion on what method or approach would be suitable for dealing with this problem? It's a hunch, but I feel as though the co-occurence of significant p-vals should count for something. That is, if the quantiles that yield a t-vals significant at .05 occur all in a row (i.e. show the effect emerging, as in the first column of the graph), that should give us confidence that they're not type 1 errors. Am I going down a dangerous road?
This is a terribly long post i'm sorry. I'd welcome all thoughts though on any part of it! My thanks to you in advance.
p.adjust
with (saytype="BH"
to control false discovery rate -- very popular nowadays because it's not as conservative as methods that control experimentwise error ...) $\endgroup$