Hello! All.
I'm using DESeq to check differential gene expression , but I got in doubt recent days and meet strange result which different from DESeq munual's demo , below I post my code and wish your kindly help:
library(DESeq2)
workDir <- "C:/Users/Administrator/Desktop/rawcounts"
setwd(workDir)
directory<-workDir
sampleFiles <- grep(".reads",list.files(directory),value=TRUE)
stage <- factor(c("B","CD4","CD8","CLP","CMP","EryA","EryB","GMP","Granulocyte","HSC","LT_HSC",
"MEP","MF","Mono","MPP","NK"))
sampleTable<-data.frame(sampleName=sampleFiles, fileName=sampleFiles,stage = stage)
dds<-DESeqDataSetFromHTSeqCount(sampleTable=sampleTable, directory=directory,design =~ stage)
dds <- dds[rowSums(counts(dds)) > 1,]
dds<-DESeq(dds)
res<-results(dds)
resOrdered <- res[order(res$padj),]
summary(resOrdered)
And I got the result:
Note:My datas have no replicates, so I wonder if this is the problem.
Thank you again for your attention!
1 answer
IMHO, you need some replicates to have a meaningful comparison. And yes, the FDR/p-values depend very much on replicates, so this might explain the things. You may see this post from Simon Anders for further details regarding running DESeq2 w/o replicates
http://seqanswers.com/forums/showpost.php?p=107433&postcount=2
Apart from that, you also need to choose which comparisons are you trying to make. In the default setting, the results() returns the comparison of the last level of the last variable in the design formula over the first level of this variable. For example, for a simple two-group comparison, this would return the log2 fold changes of the second group over the first group (the reference level).
see ?results for details
Also, see this for multifactor design https://support.bioconductor.org/p/67600/#67612
Newer and faster alternative pipeline is to use transcript abundance quantification methods such as Salmon (Patro et al. 2017), Sailfish (Patro, Mount, and Kingsford 2014), kallisto (Bray et al. 2016), or RSEM (B. Li and Dewey 2011), to estimate abundances without aligning reads.
Log in to answer this question.
