[BioC] contrast.fit and different p-values for same comparison

Lakshmanan Iyer laxvid at gmail.com
Fri Feb 12 20:17:54 CET 2010


Hi
I used limma... topTable on the "same contrast" for a complete dataset
(36 arrays, 12 conditions, 3 replicates per condition) and a subset of
the same (18 arrays).
I used lmFit, followed by contrasts.fit and eBayes in both cases.

However, I get different Ave Expression and P-values in the two cases.

Why is this happening? I cannot readily explain this.

Thanks for any help/clues.

The first row is full dataset
second row is the subset

> rbind (x[x[,1]=="ILMN_2772632",],x1[x1[,1]=="ILMN_2772632",])
                 ID    logFC  AveExpr        t      P.Value    adj.P.Val
37561  ILMN_2772632 4.198562 10.40451 22.86259 2.007820e-21 9.091611e-17
375611 ILMN_2772632 4.198562  9.85615 24.89322 1.224698e-17 1.631045e-14
              B
37561  37.44148
375611 29.95781
----------------------------------------
-best
-Lax
Here is the snippet of code with some output.

library (limma);
eset <- read.delim("results/expressions.tsv",sep="\t",stringsAsFactors=F)
rownames(eset) <- eset[,1]
eset <- eset[,-1]
names(eset) <- gsub("^X","",names(eset))
samples <- read.delim("mySamples.tsv",sep="\t",stringsAsFactors=F)
TS <- factor(samples$Condition, levels=unique(samples$Condition))
design <- model.matrix(~0 + TS)
colnames(design) <- levels(TS)
cont.matrix <- makeContrasts(Day1AgedPnxVsCtrl = Day1AgedPnx -
Day1AgedControl, Day3AgedPnxVsCtrl = Day3AgedPnx - Day3AgedControl,
Day7AgedVsCtrl = Day7AgedPnx - Day7AgedControl, Day1YoungPnxVsCtrl =
Day1YoungPnx - Day1YoungControl, Day3YoungPnxVsCtrl = Day3YoungPnx -
Day3YoungControl, Day7YoungPnxVsCtrl =Day7YoungPnx - Day7YoungControl,
levels=design)

fit <- lmFit(eset, design)
fit2 <- contrasts.fit(fit,cont.matrix)
fit2 <- eBayes(fit2)
x <- topTable(fit2, coef="Day1AgedPnxVsCtrl",num=dim(eset)[1])
#
TS1 <- factor(samples$Condition[1:18], levels=unique(samples$Condition[1:18]))
design1 <- model.matrix(~0 + TS1)
colnames(design1) <- levels(TS1)
cont.matrix1 <- makeContrasts(Day1AgedPnxVsCtrl = Day1AgedPnx -
Day1AgedControl, Day1YoungPnxVsCtrl = Day1YoungPnx - Day1YoungControl,
levels=design1)
fit1 <- lmFit(eset[,1:18], design1)
fit21 <- contrasts.fit(fit1,cont.matrix1)
fit21 <- eBayes(fit21)
topTable(fit21, coef="Day1AgedPnxVsCtrl")
x1 <- topTable(fit21, coef="Day1AgedPnxVsCtrl",num=dim(eset)[1])

> dim (eset)
[1] 45281    36
> dim (design)
[1] 36 12
> dim (cont.matrix)
[1] 12  6
> dim (design1)
[1] 18  6
> dim (cont.matrix1)
[1] 6 2
> rbind (x[x[,1]=="ILMN_2772632",],x1[x1[,1]=="ILMN_2772632",])
                 ID    logFC  AveExpr        t      P.Value    adj.P.Val
37561  ILMN_2772632 4.198562 10.40451 22.86259 2.007820e-21 9.091611e-17
375611 ILMN_2772632 4.198562  9.85615 24.89322 1.224698e-17 1.631045e-14
              B
37561  37.44148
375611 29.95781


-------------
> sessionInfo()
R version 2.10.1 (2009-12-14)
i486-pc-linux-gnu

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8
 [5] LC_MONETARY=C              LC_MESSAGES=en_US.UTF-8
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C
 [9] LC_ADDRESS=C               LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

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

other attached packages:
[1] limma_3.2.1

loaded via a namespace (and not attached):
[1] tools_2.10.1
>



More information about the Bioconductor mailing list