Hi, I got a quick question about my results for differential expression analysis based on miRNA dataset. I'm using limma package to do the DE analysis and the script I used was attached below. However, I ended up identifying a few miRNA that has very high log fold change (from 2000-80000). Since I've only worked with mRNA dataset before, and never saw a fold change this large, I'm wondering if I made any mistakes here or if there is any additional QC steps I need to do. I took a more detailed look at the specific miRNA with the extremely high fold change value across my treatment and control group and its value looks highly variable. I'm wondering if there's any QC step I should do to get rid of miRNA with high variation before doing DE analysis?
Thanks!
design <- model.matrix(~ -1 + factor(RNA.CLEAN.na$gose_max))
colnames(design) = c('Good','Bad')
fit <- lmFit(t(RNA.CLEAN.na[1:800]), design)
TBI.Baseline.contrast <- makeContrasts(gose = Good - Bad,
levels = design)
fit2 <- contrasts.fit(fit, TBI.Baseline.contrast)
fit2 <- eBayes(fit2)
TbI.baseline.results <- topTable(fit2, coef = "gose", adjust = "BH", num = Inf)
head(TbI.baseline.results)
TbI.baseline.results$ID = row.names(TbI.baseline.results)
head(TbI.baseline.results)
write.table(TbI.baseline.results, "./TBI_BL_DE_results.txt", row.names = F,
col.names = T, quote = F, sep = "\t")
you might want to look at apeglm to shrink those fold changes - https://github.com/search?q=lmfit+apeglm&type=code