[BioC] Limma :: decideTest question

[Lorenzo Bomba]

Hi all,

I would like to ask a question about the function decideTest:

I'm analysing microarray data and I'm using a the following commands to get the differentially express genes:

#RDE is a  MAlist object 

target<-read.table(file="TargetALL.txt", header=TRUE)
TS <- factor(target$Diet, levels=c("CTR0","CTR","TRT"))
design <- model.matrix(~0+TS)
colnames(design) <- c("CTR0","CTR","TRT")
fit3_1 <- lmFit(RDE,design,weights=RDE$weights)
cont.matrix <- makeContrasts(CTR0vsCTR = CTR0-CTR, CTR0vsTRT= CTR0-TRT, CTRvsTRT= CTR-TRT,levels=design)
fit3_1c <- contrasts.fit(fit3_1, cont.matrix)
fit3 <- eBayes(fit3_1c)
TableDEGfit3T <- topTable(fit3, adjust="BH", number= 29767)

# follow the "TargetALL.txt" file that I use to create the design matrix
FileName	Diet
S1M	CTR0
S2	CTR0
S3	CTR0
S4	CTR0
S6	CTR0
S10M	CTR
S11	CTR
S12	CTR
S13	CTR
S14	CTR
S15	CTR
S17rt	TRT
S18	TRT
S19	TRT
S20	TRT
S21	TRT
S22	TRT

#after I make a SEPARATE decide test I got the number of gene in common between CTR0vsCTR and CTR0vsTRT contrast like in the venn #diagramm following

results <- decideTests(fit3,p.value=0.01)
#Venn Diagrams
jpeg("Venndigramsnocov.jpg")
vennDiagram(results)
dev.off()

-------------- next part --------------

after that I tried to do separate analysis of the two control using this commands

#CONTRAST:  CTR0vsCTR
target<-read.table(file="TargetALL.txt", header=TRUE)
TS <-factor(target$Diet,levels=c("CTR0","CTR"))
design <- model.matrix(~-1+TS)
colnames(design) <- c("CTR0","CTR")
fit4_1 <- lmFit(RDEctr0ctr,design,weights=RDEctr0ctr$weights)
cont.matrix <- makeContrasts(CTR0vsCTR= CTR0-CTR,levels=design)
fit4_1c <- contrasts.fit(fit4_1, cont.matrix)
fit4 <- eBayes(fit4_1c)
TableDEGfit4T <- topTable(fit4, adjust="BH", number= 29767)


#CONTRAST:  CTR0vsTRT
target<-read.table(file="TargetALL.txt", header=TRUE)
TS <-factor(target$Diet,levels=c("CTR0","TRT"))
design <- model.matrix(~-1+TS)
colnames(design) <- c("CTR0","TRT")
fit5_1 <- lmFit(RDEctr0ctr,design,weights=RDEctr0ctr$weights)
cont.matrix <- makeContrasts(CTR0vsTRT= CTR0-TRT,levels=design)
fit5_1c <- contrasts.fit(fit5_1, cont.matrix)
fit5 <- eBayes(fit5_1c)
TableDEGfit5T <- topTable(fit5, adjust="BH", number= 29767)

and than I tried to figured out what are the gene in common that I've suppose to be the same number of the vennDiagram: 6472 genes in common

p001_4<-which(TableDEGfit4T$adj.P.Val <= 0.01)
p001_5<-which(TableDEGfit5T$adj.P.Val <= 0.01)
Descr4_001<-TableDEGfit4T$Description[p001_4]
Descr5_001<-TableDEGfit5T$Description[p001_5]
com <- match(Descr4_001, Descr5_001)
length(which(com != "NA"))

after This I got this number 4176 instead 6472
I cannot understand because in the decideTest the contrast are take separately .........
Thanks in advance for the help!

Lorenzo Bomba