Russel88
commented
6 years ago Hi jjscarpa,
There is no google group, or similar, but it might be a good idea to start one.
The intention with DAtest was to choose one method to use. I don't think that using an ensemble approach with allDA is wrong, it's just harder to interpret. Also, running testDA gives you an idea of the specificity and sensitivity of single methods, not a combination. allDA was meant as a sanity check; do the top methods find nearly the same significant features. My advice would therefore be to simply use ds2.
Regarding DESeq2, the default with multi-class predictors is for ds2(x) to run a likelihood ratio test. While you only get one p-value, you can get log fold-changes for all comparisons:
Assuming your species are called "mouse" and "cat":
res <- DA.ds2(data.new, predictor = "Sp", allResults = TRUE)
res_mouse.vs.cat <- results(res, contrast = c("predictor","mouse","cat")@listData
Alternatively, if you want a p-value for each comparison you can set out.all = FALSE in both testDA, DA.ds2 and so on. ds2 will then run Wald tests and you will get a p-value for each comparison between species, but I think this is intended for two-class predictors or multi-class predictors with a common reference you wish to compare to. But you could in principle run it for each comparison:
res <- DA.ds2(data.new, predictor = "Sp", allResults = TRUE, out.all = FALSE)
res_mouse.vs.cat <- results(res, contrast = c("predictor","mouse","cat")@listData
But note that p-values are then adjusted for each comparison independently (unless you correct them yourself), which might inflate the false positive rate.
Hi! Not an issue at all. I was just hoping to start a dialog on best practices. If there is a more appropriate venue for this discussion (e.g., a google group) please let me know.
In brief I am working with 16S rRNA data--5 host species, 10 individuals per species. I am interested in identifying differentially abundant OTUs across host species. I previously used LEfSe but am interested in exploring DAtest. I processed the data using dada2 and now have it in a phyloseq object.
I would really appreciate feedback on whether I am on the right track and how best to identify OTUs that are differentially abundant by host species.
Here is my DAtest workflow.
1) first I pre-process the data to reduce the number of features tested
data.new <- preDA(gpt_gg_no_cyano, min.samples = 5, min.reads = 100, min.abundance = 0)2) Then I use host species as the predictor to run the test.
mytest <- testDA(data.new, predictor = "Sp", effectSize = 20, out.all = TRUE)3) Results:
summary(mytest)-Method AUC FPR Spike.detect.rate Rank -DESeq2 man. geoMeans (ds2) 0.768 0.005 0.017 1.25 -Log Linear reg. 2 (llm2) 0.690 0.019 0.025 1.50 -DESeq2 (ds2x) 0.736 0.010 0.017 1.75 -ANOVA (aov) 0.732 0.023 0.017 2.25 -Log LIMMA 2 (lli2) 0.674 0.024 0.017 3.25 -Log ANOVA 2 (lao2) 0.663 0.027 0.008 5.00 -Log LIMMA (lli) 0.574 0.033 0.000 6.75 -Log ANOVA (lao) 0.569 0.034 0.000 7.50 -Log Linear reg. (llm) 0.569 0.034 0.000 7.50 -Kruskal-Wallis (kru) 0.500 0.021 0.000 8.25
The rest were NA
4) I focus here on the top 5 methods. I then run posthoc analysis on
llm2andaovand compare all 5 methods:res.allData <- allDA(data.new, predictor = "Sp", paired = NULL, covars = NULL, tests = c("ds2", "llm2", "aov", "lli2", "ds2x"), relative = TRUE, cores = (detectCores() - 1), rng.seed = 123, p.adj = "fdr", args = list(), out.all = NULL, alpha = 0.1, core.check = TRUE)5) I end up with a table of p-values by feature. My thinking here is to select OTUs that are significant for several (>3?) of the 5 methods.
6) the other issue I am having is which host species the OTUs are enriched in. My understanding is that
ds2andds2xonly do one vs. all comparisons? Is this correct?Thanks!