DESeq2 with different replicates
0
0
Entering edit mode
5.8 years ago
robertorun • 0

There are 104 samples, 38 samples are treated, 66 samples are control. The following script was used to do DGE analysis. But there is no any different expression genes, even though using the argument of dds$padj<0.05, alpha = 0.05, lfcThreshold=1. But there are more than 6000 DGEs, if you use edgeR, DEGSeq or GFOLD et al.

Would you like to tell me where is wrong? Or which soft package is more fit to do this analysis (with different replicates)? Thanks very much!

different-replicate DESeq RNA-Seq • 2.3k views
ADD COMMENT
0
Entering edit mode

code is the following:

library('DESeq2')
directory <-"./At_Count/"
sampleFiles <- grep("*.txt",list.files(directory),value=TRUE)
sampleFiles

sampleCondition <- c("A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.I","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W","A.W")
sampleTable<-data.frame(sampleName=sampleFiles, fileName=sampleFiles, condition=sampleCondition)
sampleTable

ddsHTSeq<-DESeqDataSetFromHTSeqCount(sampleTable=sampleTable, directory=directory, design=~condition)
ddsHTSeq

colData(ddsHTSeq)$condition<-factor(colData(ddsHTSeq)$condition, levels=c("A.I","A.W"))
#ddsHTSeq <- ddsHTSeq[rowSums(counts(ddsHTSeq)) >520,]
ddsHTSeq

dds<-DESeq(ddsHTSeq)
summary(dds)
table(dds$padj<0.001)

res<-results(dds,alpha = 0.05, lfcThreshold=2)

summary(res)
table(res$padj<0.05)

write.csv(as.data.frame(res),file="./At.IW.deseq2.MAplot.LFC2.Alp05.csv")

png(file="./At.IW.deseq2.MAplot.LFC2.Alp05.png",width=10,height=7.5,units="in",res=600)

layout(matrix(c(1,2,3,4,5,6), 1, 1, byrow=TRUE))
plotMA(res,ylim=c(-10,10),main="Improved VS Wild of At")

dev.off()
ADD REPLY
1
Entering edit mode

Maybe drop the lfcThreshold argument from results. When this is supplied, DESeq2 performs a test to determine if the absolute value of the fold change is greater than the value of lfcThreshold. This is totally different than asking if a gene is simply differentially expressed.

When lfcThreshold is supplied, the null hypothesis is "between these two conditions the absolute value of the fold change of a gene is not greater than lfcThreshold". When you're just looking looking at p-values, your null hypothesis is "gene x is not differentially expressed."

Just drop it out and filter the results on log2FC after.

ADD REPLY

Login before adding your answer.

Traffic: 2923 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6