Differential analysis using limma-voom
run_limma_voom(
ps,
group,
confounders = character(0),
contrast = NULL,
taxa_rank = "all",
transform = c("identity", "log10", "log10p"),
norm = "none",
norm_para = list(),
voom_span = 0.5,
p_adjust = c("none", "fdr", "bonferroni", "holm", "hochberg", "hommel", "BH", "BY"),
pvalue_cutoff = 0.05,
...
)
ps a phyloseq::phyloseq
object.
character, the variable to set the group, must be one of the var of the sample metadata.
character vector, the confounding variables to be adjusted.
default character(0)
, indicating no confounding variable.
this parameter only used for two groups comparison while there are multiple groups. For more please see the following details.
character to specify taxonomic rank to perform
differential analysis on. Should be one of
phyloseq::rank_names(phyloseq)
, or "all" means to summarize the taxa by
the top taxa ranks (summarize_taxa(ps, level = rank_names(ps)[1])
), or
"none" means perform differential analysis on the original taxa
(taxa_names(phyloseq)
, e.g., OTU or ASV).
character, the methods used to transform the microbial
abundance. See transform_abundances()
for more details. The
options include:
"identity", return the original data without any transformation (default).
"log10", the transformation is log10(object)
, and if the data contains
zeros the transformation is log10(1 + object)
.
"log10p", the transformation is log10(1 + object)
.
the methods used to normalize the microbial abundance data. See
normalize()
for more details.
Options include:
"none": do not normalize.
"rarefy": random subsampling counts to the smallest library size in the data set.
"TMM": trimmed mean of m-values. First, a sample is chosen as reference. The scaling factor is then derived using a weighted trimmed mean over the differences of the log-transformed gene-count fold-change between the sample and the reference.
"RLE", relative log expression, RLE uses a pseudo-reference calculated using the geometric mean of the gene-specific abundances over all samples. The scaling factors are then calculated as the median of the gene counts ratios between the samples and the reference.
"CSS": cumulative sum scaling, calculates scaling factors as the cumulative sum of gene abundances up to a data-derived threshold.
arguments passed to specific normalization methods. Most users will not need to pass any additional arguments here.
width of the smoothing window used for the lowess
mean-variance trend for limma::voom()
. Expressed as a proportion
between 0 and 1.
method for multiple test correction, default none
,
for more details see stats::p.adjust.
cutoff of p value, default 0.05.
extra arguments passed to limma::eBayes()
.
a microbiomeMarker
object.
contrast
must be a two length character or NULL
(default). It is only
required to set manually for two groups comparison when there are multiple
groups. The order determines the direction of comparison, the first element
is used to specify the reference group (control). This means that, the first
element is the denominator for the fold change, and the second element is
used as baseline (numerator for fold change). Otherwise, users do required
to concern this parameter (set as default NULL
), and if there are
two groups, the first level of groups will set as the reference group; if
there are multiple groups, it will perform an ANOVA-like testing to find
markers which difference in any of the groups.
Law, C. W., Chen, Y., Shi, W., & Smyth, G. K. (2014). voom: Precision weights unlock linear model analysis tools for RNA-seq read counts. Genome biology, 15(2), 1-17.
data(enterotypes_arumugam)
mm <- run_limma_voom(
enterotypes_arumugam,
"Enterotype",
contrast = c("Enterotype 3", "Enterotype 2"),
pvalue_cutoff = 0.01,
p_adjust = "none"
)
mm
#> microbiomeMarker-class inherited from phyloseq-class
#> normalization method: [ none ]
#> microbiome marker identity method: [ limma_voom ]
#> marker_table() Marker Table: [ 14 microbiome markers with 5 variables ]
#> otu_table() OTU Table: [ 244 taxa and 39 samples ]
#> sample_data() Sample Data: [ 39 samples by 9 sample variables ]
#> tax_table() Taxonomy Table: [ 244 taxa by 1 taxonomic ranks ]