[R] puzzle about contrasts
Prof Brian Ripley
ripley at stats.ox.ac.uk
Tue Sep 9 18:10:57 CEST 2008
-0.5*(A+B) is not a contrast, which is the seat of your puzzlement.
All you can get from y ~ x is an intercept (a column of ones) and a single
'contrast' column for 'x'.
If you use y ~ 0+x you can get two columns for 'x', but R does not give
you an option of what columns in the case: see the source of contrasts().
So you would need to replace contrasts(), which I think will be hard as
model.matrix.default will look in the 'stats' namespace. It would
probably be easier to create the model matrix yourself.
On Tue, 9 Sep 2008, Kenneth Knoblauch wrote:
> 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 at 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.
--
Brian D. Ripley, ripley at 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
More information about the R-help
mailing list