DA_ALDEx2 | R Documentation |
Fast run for the ALDEx2's differential abundance detection method. Support for Welch's t, Wilcoxon, Kruskal-Wallace, Kruskal-Wallace glm ANOVA-like, and glm tests.
DA_ALDEx2( object, assay_name = "counts", pseudo_count = FALSE, design = NULL, mc.samples = 128, test = c("t", "wilcox", "kw", "kw_glm", "glm"), paired.test = FALSE, denom = "all", contrast = NULL, verbose = TRUE )
object |
a phyloseq or TreeSummarizedExperiment object. |
assay_name |
the name of the assay to extract from the
TreeSummarizedExperiment object (default |
pseudo_count |
add 1 to all counts if TRUE (default
|
design |
a character with the name of a variable to group samples and
compare them or a formula to compute a model.matrix (when
|
mc.samples |
an integer. The number of Monte Carlo samples to use when estimating the underlying distributions. Since we are estimating central tendencies, 128 is usually sufficient. |
test |
a character string. Indicates which tests to perform. "t" runs Welch's t test while "wilcox" runs Wilcoxon test. "kw" runs Kruskal-Wallace test while "kw_glm" runs glm ANOVA-like test. "glm" runs a generalized linear model. |
paired.test |
A boolean. Toggles whether to do paired-sample tests.
Applies to |
denom |
An |
contrast |
character vector with exactly three elements: the name of a variable used in "design", the name of the level of interest, and the name of the reference level. If "kw" or "kw_glm" as test, contrast vector is not used. |
verbose |
an optional logical value. If |
A list object containing the matrix of p-values 'pValMat', the matrix of summary statistics for each tag 'statInfo', and a suggested 'name' of the final object considering the parameters passed to the function.
aldex
for the Dirichlet-Multinomial model
estimation. Several and more complex tests are present in the ALDEx2
framework.
set.seed(1) # Create a very simple phyloseq object counts <- matrix(rnbinom(n = 300, size = 3, prob = 0.5), nrow = 50, ncol = 6) metadata <- data.frame("Sample" = c("S1", "S2", "S3", "S4", "S5", "S6"), "group" = as.factor(c("A", "A", "A", "B", "B", "B"))) ps <- phyloseq::phyloseq(phyloseq::otu_table(counts, taxa_are_rows = TRUE), phyloseq::sample_data(metadata)) # Differential abundance with t test and denom defined by the user DA_t <- DA_ALDEx2(ps, design = "group", test = "t", denom = c(1,2,3), paired.test = FALSE, contrast = c("group", "B", "A")) # Differential abundance with wilcox test and denom = "iqlr" DA_w <- DA_ALDEx2(ps, design = "group", test = "wilcox", denom = "iqlr", paired.test = FALSE, contrast = c("group", "B", "A")) # Differential abundance with kw test and denom = "zero" # mc.samples = 2 to speed up (128 suggested) DA_kw <- DA_ALDEx2(ps, design = "group", test = "kw", denom = "zero", mc.samples = 2) # Differential abundance with kw_glm test and denom = "median" DA_kw_glm <- DA_ALDEx2(ps, design = "group", test = "kw", denom = "median", mc.samples = 2) # Differential abundance with glm test and denom = "all" DA_glm <- DA_ALDEx2(ps, design = ~ group, test = "glm", denom = "all", mc.samples = 2, contrast = c("group", "B", "A"))
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.