[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