Hello all,
I have been working on my 16S amplcon data for a while now and I have gotten to the last of the downstream analysis where I am stuck and I dont know hwo to move forward. I have data set that I woud say loks like a full factorial; Genotype (4 levels; G1, G2, G3, & G4), Day (3 levels; D1, D2 & D3), Treatment (6 levels; Control, Atrazin, PFOS, Diclo, Arsenic, wastewater) and Replicates (3 biolgical replicates of the genotypes across the time points and treatment).
I have run a differential abundance analysis using the function "ancombc2" that uses the lmerTest in its model. This i think suites my kind of data because it will allow me look for interaction among the variabels and I will also have a nested model with replicates as random effect. Please see below my code:
set.seed(123)
output2 = ancombc2(data = ps, assay_name = "counts", tax_level = "Genus",
fix_formula = "Treatment * Genotype * Day", rand_formula = "(1|Replicates)",p_adj_method = "holm", pseudo = 0, pseudo_sens = FALSE,prv_cut = 0.10, lib_cut = 0, s0_perc = 0.05,group = "Treatment", struc_zero = FALSE, neg_lb = FALSE,alpha = 0.05, n_cl = 2, verbose = FALSE,global = TRUE, pairwise = TRUE, dunnet = TRUE, trend = FALSE,iter_control = list(tol = 1e-2, max_iter = 20, verbose = FALSE),
em_control = list(tol = 1e-5, max_iter = 100),lme_control = lme4::lmerControl(),
mdfdr_control = list(fwer_ctrl_method = "holm", B = 100),
trend_control = list(contrast = NULL, node = NULL, solver = "ECOS", B = 100))
# ps = phyloseq object
I assume that the pairwise comparison will be agaisnt the base "Treatment", am not too famiiar with the meaning of the ancombc output.
The "output" has several files: global, prim pairs, and Dunn test. I can see in the 'prim' output interactions but most are false in terms of p-val but the 'global' has a different table structure with diff_abun column, W, adj_pval and the taxon. I other to move forward with this analysis, my aim is to identify ASVs,/ kegg genes that are enriched and then visualise this. but at this point I dont know how to selct the diff_adun ASVs to create a list that will be use for enrichement analysis. To clarify, I am using the amcombe package to run differential abundance analysis on both picurst2 kegg output and phyloseq object for ASVs
I would be grateful if anyone could share their thoughts on this. Thank you
res_global2.csv (5.4 KB)