4

Within the context of a 2-way ANOVA, I would like to enter contrast levels for one factor (10, 20, 30 degrees) within levels of another factor (SpeciesA, SpeciesB). This assumes the interaction effect is significant so the contrasts cannot simply be done on the temperature main effect.

I have attempted to do this within EZanova, multcomp and phia, but have been unable to figure out a syntax that works. Have I missed a package that allows this?

Contrast and coefficients

Contrast1: Does 10 degrees differ from 20 and 30 for Species A? Species:Temp (2,-1,-1,0,0,0) Contrast2: Does 10 degrees differ from 20 and 30 for Species B? Species:Temp (0,0,0,2,-1,-1) Contrast3: Does 10 degrees for Species A differ from 20 and 30 degrees for Species B? Species:Temp (2,0,0,0,-1,-1)

species<-rep(c("speciesA","speciesB"),each=12)
temp<-rep(rep(c("10","20","30"),each=4),2)
y<-rnorm(24,5,3)
(result<-anova(lm(y~species*temp)))
B. Bingham
  • 41
  • 1
  • 2
  • Check out `contrasts` in the `stats` package. There are different functions that you can use for creating a contrast matrix(e.g. `contrasts`, `contr.treatment`, `contr.poly`, `contr.helmert`, etc.) – small_data88 Sep 12 '15 at 00:00
  • I tried: `anova(lm(y ~ species*temp + C(temp, contr = matrix(c(2,-1,-1,0,0,0,0,0,0,2,-1,-1,2,0,0,0,-1,-1), nrow = 3, ncol = 6,byrow = T), how.many = 6)+ C(species, contr = matrix(c(2,-1,-1,0,0,0,0,0,0,2,-1,-1,2,0,0,0,-1,-1), nrow = 2, ncol = 6,byrow = T), how.many = 6)))` , but the output was still the same as without adding the contrast matrices with `C()`. – small_data88 Sep 12 '15 at 15:30

1 Answers1

2

It's part of base-R. See the ?contrasts page. I will illustrate for the first requested contrast; the other two should be very straightforward. I'm assuming you want these contrast p-values calculated separately, so I use the how.many argument to limit the contrast construction to a single comparison. Otherwise you would have gotten 5 contrasts, since by default the contrasts function tries to build a spanning set of orthogonal contrasts:

set.seed(123)
species <- rep(c("speciesA","speciesB"),each=12)
temp <- rep(rep(c("10","20","30"),each=4),2)
y<-rnorm(24,5,3)

intvar <- interaction(species, temp)  # create an interaction variable
?contrasts

Need to know which levels are which to properly order the contrast values:

> levels(intvar)
[1] "speciesA.10" "speciesB.10" "speciesA.20" "speciesB.20" "speciesA.30"
[6] "speciesB.30"

So this should build a contrast for 10 versus 20 or 30 angles among the A species:

 contrasts(intvar, how.many=1) <- c(2,0,-1,0,-1,0)
 anova(lm(y~intvar) )
#------------
Analysis of Variance Table

Response: y
          Df  Sum Sq Mean Sq F value Pr(>F)
intvar     1   0.013  0.0129  0.0015 0.9695
Residuals 22 190.306  8.6503

You seem to be at a fairly early stage in learning R so I would recommend that you learn to build dataframes to pass as data=-arguments to regression functions, rather than working on "loose" objects in your work-space and please do not pick up the nasty habit of using attach.

IRTFM
  • 258,963
  • 21
  • 364
  • 487
  • Thanks for the response/helpful information. The code you provided gives a clear path to contrasts on the full set of interaction means. Ideally, I would like to do that contrast within the 2-way ANOVA, meaning that the denominator for the contrast F ratio is the MSResidual from the full ANOVA (i.e., with the residuals from the full model and 18 df). I can do that by hand using the 2-way ANOVA output combined with your code. Apparently there isn't yet an easy way to do that in R without some manual calculations. – B. Bingham Sep 14 '15 at 19:11
  • There were a total of 21 df, so a "full model" with 5 contrasts would have 21-5=16 df in the residuals. If you were only going to have 3 contrasts in the "full model" (which doesn't seem full in any sense that I understand), then the number 18 would make sense in reference to the df in residuals. No manual calculations need.... just a clear specification of what is desired. This sounds like a homework problem which was why I left it with parts for the "reader" to fill in. I have never enjoyed completing someone else's homework for them (and not getting any credit.) – IRTFM Sep 15 '15 at 02:26