ancombc | R Documentation |
Determine taxa whose absolute abundances, per unit volume, of
the ecosystem (e.g., gut) are significantly different with changes in the
covariate of interest (e.g., group). The current version of
ancombc
function implements Analysis of Compositions of Microbiomes
with Bias Correction (ANCOM-BC) in cross-sectional data while allowing
for covariate adjustment.
ancombc(
data = NULL,
taxa_are_rows = TRUE,
assay.type = NULL,
assay_name = "counts",
rank = NULL,
tax_level = NULL,
aggregate_data = NULL,
meta_data = NULL,
formula,
p_adj_method = "holm",
prv_cut = 0.1,
lib_cut = 0,
group = NULL,
struc_zero = FALSE,
neg_lb = FALSE,
tol = 1e-05,
max_iter = 100,
conserve = FALSE,
alpha = 0.05,
global = FALSE,
n_cl = 1,
verbose = TRUE
)
data |
the input data. The |
taxa_are_rows |
logical. Whether taxa are positioned in the rows of the feature table. Default is TRUE. |
assay.type |
alias for |
assay_name |
character. Name of the count table in the data object
(only applicable if data object is a |
rank |
alias for |
tax_level |
character. The taxonomic level of interest. The input data
can be agglomerated at different taxonomic levels based on your research
interest. Default is NULL, i.e., do not perform agglomeration, and the
ANCOM-BC anlysis will be performed at the lowest taxonomic level of the
input |
aggregate_data |
The abundance data that has been aggregated to the desired
taxonomic level. This parameter is required only when the input data is in
|
meta_data |
a |
formula |
the character string expresses how microbial absolute
abundances for each taxon depend on the variables in metadata. When
specifying the |
p_adj_method |
character. method to adjust p-values. Default is "holm".
Options include "holm", "hochberg", "hommel", "bonferroni", "BH", "BY",
"fdr", "none". See |
prv_cut |
a numerical fraction between 0 and 1. Taxa with prevalences
(the proportion of samples in which the taxon is present)
less than |
lib_cut |
a numerical threshold for filtering samples based on library
sizes. Samples with library sizes less than |
group |
character. the name of the group variable in metadata.
The |
struc_zero |
logical. whether to detect structural zeros based on
|
neg_lb |
logical. whether to classify a taxon as a structural zero using its asymptotic lower bound. Default is FALSE. |
tol |
numeric. the iteration convergence tolerance for the E-M algorithm. Default is 1e-05. |
max_iter |
numeric. the maximum number of iterations for the E-M algorithm. Default is 100. |
conserve |
logical. whether to use a conservative variance estimator for the test statistic. It is recommended if the sample size is small and/or the number of differentially abundant taxa is believed to be large. Default is FALSE. |
alpha |
numeric. level of significance. Default is 0.05. |
global |
logical. whether to perform the global test. Default is FALSE. |
n_cl |
numeric. The number of nodes to be forked. For details, see
|
verbose |
logical. Whether to generate verbose output during the ANCOM-BC fitting process. Default is FALSE. |
A taxon is considered to have structural zeros in some (>=1)
groups if it is completely (or nearly completely) missing in these groups.
For instance, suppose there are three groups: g1, g2, and g3.
If the counts of taxon A in g1 are 0 but nonzero in g2 and g3,
then taxon A will be considered to contain structural zeros in g1.
In this example, taxon A is declared to be differentially abundant between
g1 and g2, g1 and g3, and consequently, it is globally differentially
abundant with respect to this group variable.
Such taxa are not further analyzed using ANCOM-BC, but the results are
summarized in the overall summary. For more details about the structural
zeros, please go to the
ANCOM-II paper.
Setting neg_lb = TRUE
indicates that you are using both criteria
stated in section 3.2 of
ANCOM-II
to detect structural zeros; otherwise, the algorithm will only use the
equation 1 in section 3.2 for declaring structural zeros. Generally, it is
recommended to set neg_lb = TRUE
when the sample size per group is
relatively large (e.g. > 30).
a list
with components:
feature_table
, a data.frame
of pre-processed
(based on prv_cut
and lib_cut
) microbial count table.
zero_ind
, a logical data.frame
with TRUE
indicating the taxon is detected to contain structural zeros in
some specific groups.
samp_frac
, a numeric vector of estimated sampling
fractions in log scale (natural log).
delta_em
, estimated sample-specific biases
through E-M algorithm.
delta_wls
, estimated sample-specific biases through
weighted least squares (WLS) algorithm.
res
, a list
containing ANCOM-BC primary result,
which consists of:
lfc
, a data.frame
of log fold changes
obtained from the ANCOM-BC log-linear (natural log) model.
se
, a data.frame
of standard errors (SEs) of
lfc
.
W
, a data.frame
of test statistics.
W = lfc/se
.
p_val
, a data.frame
of p-values. P-values are
obtained from two-sided Z-test using the test statistic W
.
q_val
, a data.frame
of adjusted p-values.
Adjusted p-values are obtained by applying p_adj_method
to p_val
.
diff_abn
, a logical data.frame
. TRUE if the
taxon has q_val
less than alpha
.
res_global
, a data.frame
containing ANCOM-BC
global test result for the variable specified in group
,
each column is:
W
, test statistics.
p_val
, p-values, which are obtained from two-sided
Chi-square test using W
.
q_val
, adjusted p-values. Adjusted p-values are
obtained by applying p_adj_method
to p_val
.
diff_abn
, A logical vector. TRUE if the taxon has
q_val
less than alpha
.
Huang Lin
kaul2017analysisANCOMBC
\insertReflin2020analysisANCOMBC
ancom
ancombc2
library(ANCOMBC)
if (requireNamespace("microbiome", quietly = TRUE)) {
data(atlas1006, package = "microbiome")
# subset to baseline
pseq = phyloseq::subset_samples(atlas1006, time == 0)
# run ancombc function
set.seed(123)
out = ancombc(data = pseq, tax_level = "Family",
formula = "age + nationality + bmi_group",
p_adj_method = "holm", prv_cut = 0.10, lib_cut = 1000,
group = "bmi_group", struc_zero = TRUE, neg_lb = FALSE,
tol = 1e-5, max_iter = 100, conserve = TRUE,
alpha = 0.05, global = TRUE, n_cl = 1, verbose = TRUE)
} else {
message("The 'microbiome' package is not installed. Please install it to use this example.")
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.