[BioC] contrasts in limma and missing values

[Gordon Smyth]

Dear Hans-Ulrich,

This behaviour follows directly from arithmetic in R, for example

    c(1,0) %*% c(1,NA)

produces NA rather than 0. The reason for this is that the unknown 
value corresponding to the NA could in principle have been Inf, 
meaning that 0*NA can't be inferred to be 0.

You could argue that limma should assume that infinite expression 
values are impossible, therefore it should treat 0*NA as 0, but I 
haven't so far been convinced that this finesse is worth the effort.

In your example, you could get the contrast you want without NAs 
simply by telling limma explicitly that the 3rd coefficient doesn't 
appear in your contrast, which you can do by

   fit2 <- contrasts.fit(fit[,-3],c(1,-1))

Best wishes
Gordon

>Date: Mon, 29 Oct 2007 19:42:27 +0100
>From: Hans-Ulrich Klein <h.klein at uni-muenster.de>
>Subject: [BioC] contrasts in limma and missing values
>To: BioC Mailing List <bioconductor at stat.math.ethz.ch>
>Message-ID: <47262993.8030201 at uni-muenster.de>
>Content-Type: text/plain; charset=ISO-8859-1; format=flowed
>
>Dear List,
>
>I have a question concerning missing values in limma. Consider the
>following short example:
>
>M = rbind(c(10,9,5,4,5,12,13,7,8),c(10,9,5,4,5,12,13,7,8))
>w = rbind(c(1,1,1,1,1,1,1,1,1),c(1,1,1,1,1,0,0,0,0))
>class = factor(c(1,1,2,2,2,1,1,2,2))
>lab = factor(c(1,1,1,1,1,2,2,2,2))
>
>design =  model.matrix(~ 0 + class + lab)
>contrasts = makeContrasts(contrasts="class1-class2", levels=design)
>
>fit = lmFit(M, design=design, weights=w)
>fit2 = contrasts.fit(fit, contrasts)
>
>Obviously, the effect of "lab" can not be estimated for gene 2. All
>measurements of gene 2 in lab2 have weights = 0. However, the contrast
>of interest (here class1 - class2) can be estimated by using the
>observations from lab1, but in the limma-object (fit2 in the example
>above) the contrast is NA. Can someone explain me whether/why this is
>sensible?
>
>Thanks for your help,
>Hans-Ulrich