[R-sig-ME] testing simple main effects by contrasts

Reinhold Kliegl reinhold.kliegl at gmail.com
Sat Jun 13 20:19:49 CEST 2009

>From what I understand you want to test the effect of Drug within each
of the two Groups. Here is a simple way to do this:

# Convert 2 Drug x 2 Group design to 1 x 4 Cond design
df$Cond <- factor(paste(df$Group, df$Drug, sep="_"))

# ANOVA: one main, two nested fixed effects
cmat <- matrix(c( -1, -1, +1, +1,             # Main effect of Group
	                     -1, +1,  0,  0,             # Effect of Drug | Group==A
	                       0,  0, -1, +1),  4,  3)    # Effect of Drug | Group==B
colnames(cmat) <- c(".group", ".drug_A", ".drug_B")

contrasts(df$Cond) <- cmat/2
print(mod <- lmer(value ~ Cond + (1 |Subj), data=df), cor=FALSE)
print(mod <- lmer(value ~ Cond + Block + (1 |Subj), data=df),
cor=FALSE)  # Block as fixed effect
print(mod <- lmer(value ~ Cond + (1 |Subj) + (1|Block), data=df),
cor=FALSE) # Block as random effect

The fixed-effect coefficients reported correspond to
(Intercept)  = mean of four conditions
Group  = difference between the two groups
.drug_A = difference between Drug levels in Group A
.drug_B = difference between Drug levels in Group B

and associated standard errors and "t-values".  Estimates line up with
table means for balanced designs.

You mention 2 Days of testing but this factor did not appear in the example.

The models you tested with lme and lmer where structurally different
with respect to Block. Given that you have only 3  Blocks, I would
include Block as a fixed effect, not as a random effect, but opinions
differ on this. Anyway, this has not much of a bearing on the
estimates of the Group difference and the estimates of Drug effects
within Groups.

Reinhold Kliegl

On Sat, Jun 13, 2009 at 4:57 PM, Erich
Studerus<erich.studerus at bli.uzh.ch> wrote:
> Hi,
> I have the following study design: Each of two groups of subjects (group A
> and group B) was tested repeatedly on two experimental days. At each
> experimental day, placebo or drug were administered in a randomized order
> and subjects were tested under the the influence of the drug in 3
> consecutive blocks. I have therefore one between subjects factor (Group A vs
> Group B) and two within subjects factors (Drug and Block). Here's a
> reproducible example:
> x <- expand.grid(Block= paste('B', 1:3, sep=''),
>  Drug = c('Placebo', 'Drug'), Subj=1:24)
> df <- data.frame(x, Group = rep(c('A','B'), each=72),
>  value=rnorm(144, sd=10))
> m <- model.matrix(lm(value~Group*Drug+Block,df))
> df$value <- df$value+m[,2]*7+ m[,3]*2 + m[,4]*5 + m[,5]*8 + m[,6]-2
> #delete some observations to simulate an unbalanced design
> df <- df[-as.numeric(sample(rownames(df),20)),]
> mod <- lme(value ~ Group*Drug+Block, random = ~1|Subj, data = df)
> anova(mod)
> Now, since the interaction between Drug and Group is significant, the main
> effects of Drug and Group can not be interpreted. I therefore define
> contrasts to test the simple main effects, that is, the effect of the Drug
> in each Group, averaged over all levels of the Block factor.
> I use the contrast package to set up the contrasts:
> library(contrast)
> con1 <- contrast(mod, list(Block=c('B1','B2','B3'), Group = 'A', Drug =
> 'Pla'),
>     list(Block=c('B1','B2','B3'), Group = 'A', Drug = 'Drug'),
> type='average')
> con2 <- contrast(mod, list(Block=c('B1','B2','B3'), Group = 'B', Drug =
> 'Pla'),
>     list(Block=c('B1','B2','B3'), Group = 'B', Drug = 'Drug'),
> type='average')
> Now, I use the multcomp-package to correct the p-values for multiple
> comparisons:
> library(multcomp)
> summary(glht(mod, linfct = rbind(con1$X,con2$X)))
> This approach works fine with lme, but not with lmer. In the model above I
> modeled the block factor as a fixed effect. However, according to my
> understanding, it would be more parsimonious to model the block factor as a
> random effect. I'm not directly interested in its effect and it's just an
> additional source of variance due to habituation effects . Unfortunately,
> the Block and Subj factors are crossed and therefore can not easily be
> modeled by lme. After comparing different models, the best fitting model
> with lmer looks like this:
> lmer(value ~ Group*Drug+(Drug|Subj)+(1|Block))
> How can I set up contrasts to test simple main effects for this lmer model?
> I know, that lmer has a contrasts argument. However, as far as I can see,
> one can only test contrasts for the levels of one factor at a time. Any
> comments are highly appreciated.
> Erich
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

More information about the R-sig-mixed-models mailing list