Re: [R] contrast matrix for aov

From: Prof Brian Ripley <ripley_at_stats.ox.ac.uk>
Date: Thu 10 Mar 2005 - 19:24:44 EST

On Wed, 9 Mar 2005, Darren Weber wrote:

> How do we specify a contrast interaction matrix for an ANOVA model?
>
> We have a two-factor, repeated measures design, with

Where does `repeated measures' come into this? You appear to have repeated a 2x2 experiment in each of 8 blocks (subjects). Such a design is usually analysed with fixed effects. (Perhaps you averaged over repeats in the first few lines of your code?)

> Cue Direction (2) x Brain Hemisphere(2)
>
> Each of these has 2 levels, 'left' and 'right', so it's a simple 2x2 design
> matrix. We have 8 subjects in each cell (a balanced design) and we want to
> specify the interaction contrast so that:
>
> CueLeft>CueRght for the Right Hemisphere
> CueRght>CueLeft for the Left Hemisphere.
>
> Here is a copy of the relevant commands for R:
>
> ########################################
> lh_cueL <- rowMeans( LHroi.cueL[,t1:t2] )
> lh_cueR <- rowMeans( LHroi.cueR[,t1:t2] )
> rh_cueL <- rowMeans( RHroi.cueL[,t1:t2] )
> rh_cueR <- rowMeans( RHroi.cueR[,t1:t2] )
> roiValues <- c( lh_cueL, lh_cueR, rh_cueL, rh_cueR )
>
> cuelabels <- c("CueLeft", "CueRight")
> hemlabels <- c("LH", "RH")
>
> roiDataframe <- data.frame( roi=roiValues, Subject=gl(8,1,32,subjectlabels),
> Hemisphere=gl(2,16,32,hemlabels), Cue=gl(2,8,32,cuelabels) )
>
> roi.aov <- aov(roi ~ (Cue*Hemisphere) + Error(Subject/(Cue*Hemisphere)),
> data=roiDataframe)

I think the error model should be Error(Subject). In what sense are `Cue' and `Cue:Hemisphere' random effects nested inside `Subject'?

Let me fake some `data':

set.seed(1); roiValues <- rnorm(32)
subjectlabels <- paste("V"1:8, sep = "")
options(contrasts = c("contr.helmert", "contr.poly")) roi.aov <- aov(roi ~ Cue*Hemisphere + Error(Subject), data=roiDataframe)

> roi.aov

Call:
aov(formula = roi ~ Cue * Hemisphere + Error(Subject), data = roiDataframe)

Grand Mean: 0.1165512

Stratum 1: Subject

Terms:

                 Residuals
Sum of Squares   4.200946
Deg. of Freedom         7

Residual standard error: 0.7746839

Stratum 2: Within

Terms:

                       Cue Hemisphere Cue:Hemisphere Residuals
Sum of Squares   0.216453   0.019712       0.057860 21.896872
Deg. of Freedom         1          1              1        21

Residual standard error: 1.021131
Estimated effects are balanced

Note that all the action is in one stratum, and the SSQs are the same as

aov(roi ~ Subject + Cue * Hemisphere, data = roiDataframe)

(and also the same as for your fit).

> print(summary(roi.aov))

It auto-prints, so you don't need print().

> ########################################
>
>
> I've tried to create a contrast matrix like this:
>
> cm <- contrasts(roiDataframe$Cue)
>
> which gives the main effect contrasts for the Cue factor. I really want to
> specify the interaction contrasts, so I tried this:
>
> ########################################
> # c( lh_cueL, lh_cueR, rh_cueL, rh_cueR )
> # CueRight>CueLeft for the Left Hemisphere.
> # CueLeft>CueRight for the Right Hemisphere
>
> cm <- c(-1, 1, 1, -1)
> dim(cm) <- c(2,2)

(That is up to sign what Helmert contrasts give you.)

> roi.aov <- aov( roi ~ (Cue*Hemisphere) + Error(Subject/(Cue*Hemisphere)),
> contrasts=cm, data=roiDataframe)
> print(summary(roi.aov))
> ########################################
>
> but the results of these two aov commands are identical. Is it the case that
> the 2x2 design matrix is always going to give the same F values for the
> interaction regardless of the contrast direction?

Yes, as however you code the design (via `contrasts') you are fitting the same subspaces. Not sure what you mean by `contrast direction', though.

However, you have not specified `contrasts' correctly:

contrasts: A list of contrasts to be used for some of the factors in

           the formula.

and cm is not a list, and an interaction is not a factor.

> OR, is there some way to get a summary output for the contrasts that is
> not available from the print method?

For more than two levels, yes: see `split' under ?summary.aov. Also, see se.contrasts which allows you to find the standard error for any contrast.

For the fixed-effects model you can use summary.lm:

> fit <- aov(roi ~ Subject + Cue * Hemisphere, data = roiDataframe)
> summary(fit)

                Df  Sum Sq Mean Sq F value Pr(>F)
Subject         7  4.2009  0.6001  0.5756 0.7677
Cue             1  0.2165  0.2165  0.2076 0.6533
Hemisphere      1  0.0197  0.0197  0.0189 0.8920
Cue:Hemisphere  1  0.0579  0.0579  0.0555 0.8161
Residuals      21 21.8969  1.0427

> summary.lm(fit)

Call:
aov(formula = roi ~ Subject + Cue * Hemisphere, data = roiDataframe)

Residuals:

     Min 1Q Median 3Q Max -1.7893 -0.4197 0.1723 0.5868 1.3033

Coefficients:

                  Estimate Std. Error t value Pr(>|t|)
[...]
Cue1             -0.08224    0.18051  -0.456    0.653
Hemisphere1       0.02482    0.18051   0.137    0.892
Cue1:Hemisphere1 -0.04252 0.18051 -0.236 0.816

where the F values are the squares of the t values.

-- 
Brian D. Ripley,                  ripley@stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford,             Tel:  +44 1865 272861 (self)
1 South Parks Road,                     +44 1865 272866 (PA)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595

______________________________________________
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Received on Thu Mar 10 19:53:16 2005

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:30:42 EST