[R] limma contrast matrix

[James W. MacDonald]

Hi John,

The limma package is part of Bioconductor, so you should be asking this 
question on the BioC listserv <bioconductor at r-project.org> rather than 
R-help.

In addition, please see the posting guide. The example you give is not 
self-contained, and the term 'does not work' is so unspecific as to be 
useless.

On 8/4/2011 8:53 AM, Belmont, John W wrote:
> I am trying to correct for the effect of 2 covariates in a gene expression data set.
>
>
> design<-model.matrix(~0 + Factor + cov1 + cov2)
>
>
> QUESTION:
> How to set up the contrast matrix?
>
> The usual commands
>
> fit<- lmFit(selDataMatrix, design)
> cont.matrix<- makeContrasts(FacCont=Group1-Group2, levels=design)

Why would your design matrix have columns labeled 'Group1' and 'Group2'? 
Your code doesn't indicate that you did any such naming.

> fit2<- contrasts.fit(fit, cont.matrix)
>
> does not work because the original design matrix includes the covariates.
>
> I think I don't really understand how the contrast matrix works.

I agree. So let's make some wild assumptions, and maybe I can help.

Factor <- factor(rep(1:2, each=6))
cov1 <- factor(rep(1:2, times=2, each=3))
## let's say you had two batches.
cov2 <- rnorm(12)
## and the other covariate is continuous

mat <- model.matrix(~ 0 + Factor + cov1 + cov2)

 > mat
    Factor1 Factor2 cov12        cov2
1        1       0     0  0.75675940
2        1       0     0 -0.80761696
3        1       0     0  0.61228480
4        1       0     1 -1.13920820
5        1       0     1 -0.24367358
6        1       0     1  0.32244694
7        0       1     0  0.51438468
8        0       1     0 -2.23587057
9        0       1     0 -0.06560733
10       0       1     1 -0.11273432
11       0       1     1  0.33398074
12       0       1     1  0.70900581

Now we have a model that controls for the fact that you did things in 
two batches (naughty boy) and that you have some random continuous 
covariate you want to control for. So Factor1 is the mean of the first 
group after controlling for the other two covariates, and Factor2 is the 
same, for the other group. Can you see that a contrast comparing these 
two factors is what you are after?

So something like

makeContrasts(whatIwant=Factor1-Factor2, groups=design)

should get you where you want to go.

Best,

Jim


>
>
> John W. Belmont, M.D.,Ph.D.
> Professor
> Department of Molecular and Human Genetics
> Baylor College of Medicine
> 1100 Bates, Room 8070
> Houston, TX 77030
> 713-798-4634
> fax: 713-798-7187
>
> ______________________________________________
> 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.

-- 
James W. MacDonald, M.S.
Biostatistician
Douglas Lab
University of Michigan
Department of Human Genetics
5912 Buhl
1241 E. Catherine St.
Ann Arbor MI 48109-5618
734-615-7826
**********************************************************
Electronic Mail is not secure, may not be read every day, and should not be used for urgent or sensitive issues