[BioC] limma - different parametrization and weights

Hans-Ulrich Klein h.klein at uni-muenster.de
Tue Nov 4 11:56:17 CET 2008


Dear All,

I used limma with two different parametrizations. Both approaches should 
be equivalent in my opinion. However, if I use weights, the results 
differ. (Results of both approaches are equal without weights.) I 
attached an example below. Does someone know the reason for this?

Regards,
Hans-Ulrich



Here is the example:

library("limma")

n=30
m=100
data = matrix(rnorm(n*m, mean=8, sd=1), ncol=n, nrow=m)
W = matrix(rbinom(n*m, 1, p=0.8), ncol=n, nrow=m)
W[W==0] = 1/2
disease = factor(c(rep("D1", n/3), rep("D2", n/3), rep("D3", n/3)))
batch = factor(sample(c("B1", "B2"), n, replace=TRUE))

D1 = model.matrix(~ disease + batch)
D2 = model.matrix(~ 0 + disease + batch)

fit1 = lmFit(data, D1, weights=W)
fit2 = lmFit(data, D2, weights=W)

contrs = makeContrasts(contrasts=
    c("diseaseD2-diseaseD1", "diseaseD3-diseaseD1"), levels=D2)
fit2c = contrasts.fit(fit2, contrs)

eFit1 = eBayes(fit1)
eFit2 = eBayes(fit2c)

topTable(eFit1, coef="diseaseD2")
        logFC  AveExpr         t    P.Value adj.P.Val         B
71  1.0992458 8.414207  2.547634 0.01148965 0.6240313 -4.461615
79  1.0701491 8.306736  2.412824 0.01660306 0.6240313 -4.482316
65  0.9004096 7.956182  2.061818 0.04033472 0.6240313 -4.515367
73  0.8769716 7.822919  2.014672 0.04508839 0.6240313 -4.519390
97 -0.8838333 8.059430 -1.969606 0.05006854 0.6240313 -4.526330
28 -0.9103940 8.163924 -2.000665 0.04658917 0.6240313 -4.527588
92 -0.8426017 7.948664 -1.885870 0.06055682 0.6240313 -4.530451
78  0.8518031 7.921508  1.928056 0.05506380 0.6240313 -4.531817
96 -0.8062308 8.137124 -1.833017 0.06807629 0.6240313 -4.542318
76  0.7613192 8.176955  1.771061 0.07785820 0.6240313 -4.543365

topTable(eFit2, coef="diseaseD2-diseaseD1")
        logFC  AveExpr         t    P.Value adj.P.Val         B
71  1.0992458 8.414207  2.525839 0.01220672 0.6062981 -4.466368
79  1.0701491 8.306736  2.344510 0.01989314 0.6062981 -4.495402
28 -0.9103940 8.163924 -2.104442 0.03641183 0.6062981 -4.510339
65  0.9004096 7.956182  2.083941 0.03825580 0.6062981 -4.511464
73  0.8769716 7.822919  2.003999 0.04622820 0.6062981 -4.521191
78  0.8518031 7.921508  1.993785 0.04734180 0.6062981 -4.521311
92 -0.8426017 7.948664 -1.905554 0.05793940 0.6062981 -4.527255
97 -0.8838333 8.059430 -1.942207 0.05331756 0.6062981 -4.530611
45  0.7744936 7.745213  1.769790 0.07807041 0.6062981 -4.544968
8  -0.7040262 7.806587 -1.693483 0.09170026 0.6062981 -4.547372


sessionInfo()
R version 2.8.0 (2008-10-20)
x86_64-pc-linux-gnu

locale:
C

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

other attached packages:
[1] limma_2.16.2


-- 
Hans-Ulrich Klein
Department of Medical Informatics and Biomathematics
University of Münster, Germany
Tel.: +49 (0)251 83-58405



More information about the Bioconductor mailing list