[BioC] limma: NA coefficients due to missing values

[Gordon K Smyth]

Dear Tao,

The coefficients for days 2-4 have already been estimated without NAs, see 
fit1$coef[,11:13].

Using contrasts will introduce NAs unless you remove the missing patient 
effect.  This can be done by,

  > fit.day <- fit1[,11:13]
  > colnames(fit.day)
  [1] "day2" "day3" "day4"
  > head(fit.day$coef)
            day2        day3        day4
  [1,]  0.1018951  0.14728872  0.01348279
  [2,]  0.8222967  0.95183289  0.02500825
  [3,]  0.6601298 -0.27436185  0.16624044
  [4,] -0.7082047  0.04495111  0.25636112
  [5,] -0.0628469 -0.04033322 -0.59622820
  [6,] -0.0296513 -0.31283864  0.14225534

Then you can happily take as many contrasts as you like of the day 
effects, using fit3 instead of fit1.  For example:

   cont.mat <- makeContrasts(day3-day2, levels=colnames(fit.day))

Best wishes
Gordon

---------------------------------------------
Professor Gordon K Smyth,
Bioinformatics Division,
Walter and Eliza Hall Institute of Medical Research,
1G Royal Parade, Parkville, Vic 3052, Australia.
http://www.statsci.org/smyth

On Tue, 14 May 2013, Shi, Tao wrote:

> Dear Gordon and list,

Please see the example below and R output.  The coefficients for Probe #2 
can't be estimated b/c patient 10 data were missing.  But in fact, they 
can be estimated just using data from other 9 patients.  Is there a way to 
tell limma to do that?

Thanks!

Tao



library(limma)
dat <- matrix(rnorm(400), ncol=40)
ptID <- factor(rep(1:10, each=4))
day <- factor(rep(1:4,10))
dat[2, 37:40] <- NA
design <- model.matrix(~0+ptID+day)
contrast.matrix <- makeContrasts(day2=day2, day3=day3,day4=day4, levels=design)
fit1 <- lmFit(dat, design=design)
fit2 <- eBayes(contrasts.fit(fit1, contrast.matrix))


Warning message:
Partial NA coefficients for 1 probe(s) 
> fit2$p.value
       Contrasts
              day2      day3      day4
   [1,] 0.40027892 0.4917185 0.0387787
   [2,]         NA        NA        NA
   [3,] 0.05281279 0.4914354 0.6088744
   [4,] 0.72648798 0.6816108 0.7551084
   [5,] 0.35528206 0.7962423 0.9827632
   [6,] 0.09552328 0.3458950 0.1327465
   [7,] 0.76918217 0.5945065 0.4511735
   [8,] 0.40125232 0.1066155 0.1269416
   [9,] 0.56923379 0.2981261 0.9285374
  [10,] 0.04276991 0.4693320 0.3139434
> fit2$coef
       Contrasts
              day2       day3         day4
   [1,]  0.3737307  0.3054398  0.921284513
   [2,]         NA         NA           NA
   [3,] -0.8628711 -0.3056396 -0.227256710
   [4,] -0.1553395  0.1821986  0.138509886
   [5,] -0.4107858  0.1146613 -0.009593154
   [6,]  0.7421089 -0.4188821  0.668951727
   [7,] -0.1303084 -0.2364268 -0.334736056
   [8,] -0.3729581 -0.7182350 -0.679192601
   [9,]  0.2528092  0.4624636  0.039823091
  [10,] -0.9030535 -0.3214420 -0.447554386


> sessionInfo()
R version 3.0.0 (2013-04-03)
Platform: x86_64-w64-mingw32/x64 (64-bit)

locale:
[1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United States.1252    LC_MONETARY=English_United States.1252 LC_NUMERIC=C                          
[5] LC_TIME=English_United States.1252   

attached base packages:
[1] grDevices datasets  splines   graphics  stats     tcltk     utils     methods   base    

other attached packages:
[1] limma_3.16.3    svSocket_0.9-55 TinnR_1.0-5     R2HTML_2.2.1    Hmisc_3.10-1.1  survival_2.37-4

loaded via a namespace (and not attached):
[1] cluster_1.14.4  grid_3.0.0      lattice_0.20-15 svMisc_0.9-69   tools_3.0.0   
>

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