Entering edit mode
2.8 years ago
vinishavvenugopal
▴
30
Hi,
The plotMA in my analysis shows only one significant gene. What can I infer from this and can I even go ahead with plotting volcano plot and functional analysis using genes with pvalue<0.05 rather than adjpvalue<0.05, also considering the fact there are no genes having high normalized counts with pvalue<0.05.
How do I proceed with this analysis?
dds <- DESeqDataSetFromTximport(txi, met, ~Type)
dds <- dds[rowMeans(counts(dds)) > 20,]
dds <- estimateSizeFactors(dds)
dds <- DESeq(dds)
plotMA(dds, alpha=0.05,ylim=c(-2,2), lfcThreshold = 1.5, colSig = "blue")
res <- results(dds, contrast=c("Type","WT","P154S"))
res.out <- as.data.frame(res)
res <- res[!is.na(res$pvalue),]
resOrdered <- res[(res$pvalue < 0.05),]
ressig <- resOrdered[(resOrdered$log2FoldChange < -1.5) | (resOrdered$log2FoldChange >1.5) ,]
Any help is appreciated. Thanks, Vinisha
Take the following with a huge grain of salt. Your experiment is showing you that whatever condition you are measuring, it doesn't have a large impact on gene expression space (alternatively, the experiment could be noisy enough that you're not able to easily see the effect given the statistics). One message here is that you need to do a better or different experiment to understand the relationship between your condition and gene expression. From an experimentalist's perspective, I would say the experiment is not completely dead yet. You could ignore the magnitude of the statistics (i.e. a p-value cutoff) and simply take the top x number of genes up or down as ranked by p-value, and look for enrichment or clues in those sets. For instance, of all genes with positive logFC, rank by p-value and take the top 300. The idea here is to see if your experiment has pushed any biology towards the edges. If you see something that makes any sense given your system, this might be insightful and help you design the next experiment. This experiment itself will not be publishable, but may help you take a next step, or rule something out, like many experiments with negative results. The risk is that there really is nothing useful here, and you could waste of a lot of time sifting through a negative result, time that would be better spent elsewhere. The danger with genomic data is that if you look hard enough you can always find something, even when those somethings are red herrings.
That makes a lot of sense. Thanks!
When I use edgeR, I could come up with more significant genes relative to deseq. Is it fine if I proceed with edgeR?
How do I determine which type of differential expression analysis to use? (deseq,edgeR)
One of my friends likes pizza more than the other. Which one do I believe? You can use either one. EdgeR is known for calling more genes DE than DESeq2, but this is simply the difference in where and how to draw the line between false positives and false negatives under the hood. If you search the literature you can find some comparisons, but I'd say it doesn't really matter. They are both fantastic, and your biology will not, and should not, depend on one versus the other. At the end of the day, true positives and true negatives should not depend on software.
Thank you for such an amazing answer, and I am forced to appreciate the beauty/simplicity of this answer. Finally, after so long, someone has come up with an excellent explaination of lack of DEGs from RNA seq data.
you have limited the y-axis to
ylim=c(-2,2). Probably you see some more genes by relaxing that limitAnyway, you can try to select
lfcThreshold=1.0(two-fold up/down-regulated) if you don't have enough significant genes.No. adjPval is there to do
multiple hypothesis testing corrections(search here or google), and you can't avoid thatThanks Santosh. But I think the ylim doesn't matter here. If you can see the blue triangle ones, that's basically when it's beyond the limit you have used.
I really appreciate your response. I'll check adjPvalue.