[BioC] Design and not estimable coefficients

[Gordon K Smyth]

Dear Dave,

This is a common type of design, with two groups and a batch effect.  To 
find genes DE in mu vs wt adjusting for exp differences:

   mutation <- factor(target$mutation)
   mutation <- relevel(mutation, ref="wt")
   exp <- factor(target$exp)
   design <- model.matrix(~exp+mutation)
   fit <- lmFit(x,design)
   fit <- eBayes(fit)
   topTable(fit, coef=3)

Note that the design matrix has only 3 columns, not four as yours do. 
The batch effect requires only one extra column.

Your design matrix has a redundant column because you've defined the 
intercept term twice.  The first and second columns (wt and mu) add up to 
the intercept column, and so do the 3rd and 4th columns (exp1 and exp2).

You say that limma failed on your design matrix, but really it didn't. 
It correctly removed the redundant column from your design matrix, and 
gave you a warning to let you know.

Best wishes
Gordon

---------------- original message ----------------
[BioC] Design and not estimable coefficients
Dave Canvhet dcanvhet at gmail.com
Wed Apr 25 18:08:59 CEST 2012

Dear all,

I have the following target

target

filename mutation exp
sample1.cel   wt     1
sample2.cel   wt     1
sample3.cel  mu     1
sample4.cel  mu     1
sample5.cel   wt     2
sample6.cel   mu     2
sample7.cel  mu     2
sample8.cel  mu     2


I'm intersting in mu vs wt, but taking into account the exp factor (which 
seems having most impact on signal as reveal by a PCA : PC1 seems to be 
assoiated to exp)

So I've set the following design so as to model the resulting signal as a
combination of these two factor:

design = matrix(0, ncol=4, nrow = 8)
design[which(target$mutation == "wt"),1] = 1
design[which(target$mutation == "mu"),2] = 1
design[which(target$exp == 1),3] = 1
design[which(target$exp == 2),4] = 1
colnames(design) = c("wt","mu","exp1","exp2")


design
      wt mu exp1 exp2
[1,]  1  0    1    0
[2,]  1  0    1    0
[3,]  0  1    1    0
[4,]  0  1    1    0
[5,]  1  0    0    1
[6,]  0  1    0    1
[7,]  0  1    0    1
[8,]  0  1    0    1

fit = lmFit(x, design)
but it failed (even before doinf any contrast matrix), raising the
following error :
Coefficients not estimable: exp2

Can you please tell me what is wrong is such design ?


SO I've used the approach described in Limma User guide :
f <- paste(target$mutation,target$exp,sep="")
f <- factor(f)
[1] wt1 wt1 mu1 mu1 wt2 mu2 mu2 mu2
Levels: mu1 mu2 wt1 wt2


design <- model.matrix(~0+f)
colnames(design) <- levels(f)
design

   mu1 mu2 wt1 wt2
1   0   0   1   0
2   0   0   1   0
3   1   0   0   0
4   1   0   0   0
5   0   0   0   1
6   0   1   0   0
7   0   1   0   0
8   0   1   0   0
attr(,"assign")
[1] 1 1 1 1
attr(,"contrasts")
attr(,"contrasts")$f
[1] "contr.treatment"

using such a design matrix, can I use correctly the following contrast
matrix to get the gene differentially expressed between mutant and wild
type ?
cont.matrix <- makeContrasts(mu_vs_wt = (mu2+mu1) - 
(wt2+wt1),levels=design)


many thanks for you answer.

==
Dave Canvhet

______________________________________________________________________
The information in this email is confidential and intend...{{dropped:4}}