Merge pull request #65 from yangchen2/master

get_pvalue_tool_columns using p-value columns instead of corrected P-value columns

qadabra/workflow/rules/common.smk

@@ -50,14 +50,24 @@ def get_pvalue_tool_columns(wildcards):
     target = d["target_level"]
     reference = d["reference_level"]
 
+    # columns = {
+    #     "edger": "PValue",
+    #     "deseq2": "pvalue",
+    #     "ancombc": "pvals",
+    #     "aldex2": f"model.{covariate}{target} Pr(>|t|)",
+    #     "maaslin2": "pval",
+    #     "metagenomeseq": "pvalues",
+    #     "corncob": "fit.p",
+    # }
+
     columns = {
-        "edger": "PValue",
-        "deseq2": "pvalue",
-        "ancombc": "pvals",
-        "aldex2": f"model.{covariate}{target} Pr(>|t|)",
-        "maaslin2": "pval",
-        "metagenomeseq": "pvalues",
-        "corncob": "fit.p",
+        "edger": "PValue_BH_adj",
+        "deseq2": "padj",
+        "ancombc": "qvals",
+        "aldex2": f"model.{covariate}{target} Pr(>|t|).BH",
+        "maaslin2": "pval_BH_adj",
+        "metagenomeseq": "adjPvalues",
+        "corncob": "adjusted_p_values",
     }
     return columns[wildcards.tool]
 

qadabra/workflow/scripts/R/ancombc.R

@@ -57,15 +57,15 @@ coef_col <- paste("coefs", covariate, target, sep = ".")
 pval_col <- paste("pvals", covariate, target, sep = ".")
 coefs <- ancombc.results$res$beta
 pvals <- ancombc.results$res$p_val
-# qvals <- ancombc.results$res$q_val
+qvals <- ancombc.results$res$q_val
 
 # Modify result names
 row.names(coefs) <- gsub("^F_", "", row.names(coefs))
 row.names(pvals) <- gsub("^F_", "", row.names(pvals))
-# row.names(qvals) <- gsub("^F_", "", row.names(qvals))
+row.names(qvals) <- gsub("^F_", "", row.names(qvals))
 
-results_all <- data.frame(coefs=coefs, pvals=pvals)
-colnames(results_all) <- c("coefs", "pvals")
+results_all <- data.frame(coefs=coefs, pvals=pvals, qvals=qvals)
+colnames(results_all) <- c("coefs", "pvals", "qvals")
 
 # Save results to output file
 write.table(results_all, file=snakemake@output[[1]], sep="\t")