Hi,
I'm trying to redefine the contrasts for a linear model.
With a 2 level factor, x, with levels A and B, a two level
factor outputs A and B - A from an lm fit, say
lm(y ~ x). I would like to set the contrasts so that
the coefficients output are -0.5 (A + B) and B - A,
but I can't get the sign correct for the first coefficient
(Intercept).
Here is a toy example,
set.seed(12161952)
y <- rnorm(10)
x <- factor(rep(letters[1:2], each = 5))
## so A and B =
tapply(y, x, mean)
a b
-0.7198888 0.8323837
## and with treatment contrasts
coef(lm(y ~ x)) ## A and B - A
(Intercept) xb
-0.7198888 1.5522724
Then, I try to redefine the contrasts
### would like contrasts: -0.5 (A + B) and B - A
D1 <- matrix( c(-0.5, -0.5,
-1, 1),
2, 2, byrow = TRUE)
C1 <- solve(D1)
Cnt <- C1[, -1]
contrasts(x) <- Cnt
coef(lm(y ~ x))
(Intercept) x1
0.05624745 1.55227241
but note that the desired value is
-0.5 * sum(tapply(y, x, mean))
[1] -0.05624745
I note that the first column of C1 is -1's not +1's
and that working by hand, if I tamper with the model matrix
mm <- model.matrix(y ~ x)
mm[, 1] <- -1
mm
(Intercept) x1
1 -1 -0.5
2 -1 -0.5
3 -1 -0.5
4 -1 -0.5
5 -1 -0.5
6 -1 0.5
7 -1 0.5
8 -1 0.5
9 -1 0.5
10 -1 0.5
attr(,"assign")
[1] 0 1
attr(,"contrasts")
attr(,"contrasts")$x
[,1]
a -0.5
b 0.5
solve(t(mm) %*% mm) %*% t(mm) %*% y ##Yes, I know. Use QR
[,1]
(Intercept) -0.05624745
x1 1.55227241
gives the correct sign.
So, I guess my question reduces to how one would set the
contrasts for the model.matrix to be correct
for this to work out correctly?
Thank you.
Ken
--
Ken Knoblauch
Inserm U846
Institut Cellule Souche et Cerveau
Département Neurosciences Intégratives
18 avenue du Doyen Lépine
69500 Bron
France
tel: +33 (0)4 72 91 34 77
fax: +33 (0)4 72 91 34 61
portable: +33 (0)6 84 10 64 10
http://www.sbri.fr
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
