Asking, are simple effects different from 0
On 3/5/2008 10:09 AM, jebyrnes wrote:
Huh. Very interesting. I haven't really worked with manipulating contrast matrices before, save to do a prior contrasts. Could you explain the matrix you laid out just a bit more so that I can generalize it to my case?
Each column corresponds to one of the coefficients in the model, and
each row specifies a particular contrast. The numbers in the matrix
indicate how the model coefficients are combined to indicate a
particular difference in means.
For example, the first row indicates that the third coefficient
(woolB) is multiplied by -1. The baseline categories are A and L for
the wool and tension factors, so the woolB effect in fm is the simple
effect of B vs. A in the baseline category of the tension factor.
Multiplying this coefficient by -1 produces an A vs. B comparison in the
baseline category of the tension factor.
When I want to check that contrasts are as intended, I use contrast()
in the contrast package (by Steve Weston, Jed Wing, & Max Kuhn). That
function allows you to specify factor levels by name to construct the
contrast. For example:
library(contrast)
# M vs. H at B
contrast(fm, a=list(tension = "M", wool = "B"),
b=list(tension = "H", wool = "B"))
lm model parameter contrast
Contrast S.E. Lower Upper t df Pr(>|t|)
10 5.157299 -0.3694453 20.36945 1.94 48 0.0584
It also allows you to print the design matrix for a contrast:
contrast(fm, a=list(tension = "M", wool = "B"),
b=list(tension = "H", wool = "B"))$X
(Intercept) tensionM tensionH woolB tensionM:woolB tensionH:woolB
1 0 1 -1 0 1 -1
Chuck Cleland wrote:
One approach would be to use glht() in the multcomp package. You
need to work out how to formulate the matrix of coefficients that give
the desired contrasts. Here is an example using the warpbreaks data
frame:
fm <- lm(breaks ~ tension*wool, data=warpbreaks)
# names(coef(fm))
# (Intercept) tensionM tensionH woolB tensionM:woolB tensionH:woolB
cm <- rbind(
"A vs. B at L" = c(0, 0, 0,-1, 0, 0),
"A vs. B at M" = c(0, 0, 0,-1,-1, 0),
"A vs. B at H" = c(0, 0, 0,-1, 0,-1),
"M vs. L at A" = c(0, 1, 0, 0, 0, 0),
"M vs. H at A" = c(0, 1,-1, 0, 0, 0),
"L vs. H at A" = c(0, 0,-1, 0, 0, 0),
"M vs. L at B" = c(0, 1, 0, 0, 1, 0),
"M vs. H at B" = c(0, 1,-1, 0, 1,-1),
"L vs. H at B" = c(0, 0,-1, 0, 0,-1))
library(multcomp)
summary(glht(fm, linfct = cm), test = adjusted(type="none"))
Simultaneous Tests for General Linear Hypotheses
Fit: lm(formula = breaks ~ tension * wool, data = warpbreaks)
Linear Hypotheses:
Estimate Std. Error t value p value
A vs. B at L == 0 16.3333 5.1573 3.167 0.002677 **
A vs. B at M == 0 -4.7778 5.1573 -0.926 0.358867
A vs. B at H == 0 5.7778 5.1573 1.120 0.268156
M vs. L at A == 0 -20.5556 5.1573 -3.986 0.000228 ***
M vs. H at A == 0 -0.5556 5.1573 -0.108 0.914665
L vs. H at A == 0 20.0000 5.1573 3.878 0.000320 ***
M vs. L at B == 0 0.5556 5.1573 0.108 0.914665
M vs. H at B == 0 10.0000 5.1573 1.939 0.058392 .
L vs. H at B == 0 9.4444 5.1573 1.831 0.073270 .
---
Signif. codes: 0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1
(Adjusted p values reported -- none method)
Chuck Cleland, Ph.D. NDRI, Inc. 71 West 23rd Street, 8th floor New York, NY 10010 tel: (212) 845-4495 (Tu, Th) tel: (732) 512-0171 (M, W, F) fax: (917) 438-0894