getPrevalence | R Documentation |
These functions calculate the population prevalence for taxonomic ranks in a
SummarizedExperiment-class
object.
getPrevalence(x, ...)
## S4 method for signature 'ANY'
getPrevalence(
x,
detection = 0,
include_lowest = FALSE,
sort = FALSE,
na.rm = TRUE,
...
)
## S4 method for signature 'SummarizedExperiment'
getPrevalence(
x,
assay.type = assay_name,
assay_name = "counts",
as_relative = FALSE,
rank = NULL,
...
)
getPrevalentFeatures(x, ...)
## S4 method for signature 'ANY'
getPrevalentFeatures(x, prevalence = 50/100, include_lowest = FALSE, ...)
## S4 method for signature 'SummarizedExperiment'
getPrevalentFeatures(
x,
rank = NULL,
prevalence = 50/100,
include_lowest = FALSE,
...
)
getPrevalentTaxa(x, ...)
## S4 method for signature 'ANY'
getPrevalentTaxa(x, ...)
getRareFeatures(x, ...)
## S4 method for signature 'ANY'
getRareFeatures(x, prevalence = 50/100, include_lowest = FALSE, ...)
## S4 method for signature 'SummarizedExperiment'
getRareFeatures(
x,
rank = NULL,
prevalence = 50/100,
include_lowest = FALSE,
...
)
getRareTaxa(x, ...)
## S4 method for signature 'ANY'
getRareTaxa(x, ...)
subsetByPrevalentFeatures(x, ...)
## S4 method for signature 'SummarizedExperiment'
subsetByPrevalentFeatures(x, rank = NULL, ...)
subsetByPrevalentTaxa(x, ...)
## S4 method for signature 'ANY'
subsetByPrevalentTaxa(x, ...)
subsetByRareFeatures(x, ...)
## S4 method for signature 'SummarizedExperiment'
subsetByRareFeatures(x, rank = NULL, ...)
subsetByRareTaxa(x, ...)
## S4 method for signature 'ANY'
subsetByRareTaxa(x, ...)
getPrevalentAbundance(
x,
assay.type = assay_name,
assay_name = "relabundance",
...
)
## S4 method for signature 'ANY'
getPrevalentAbundance(
x,
assay.type = assay_name,
assay_name = "relabundance",
...
)
## S4 method for signature 'SummarizedExperiment'
getPrevalentAbundance(x, assay.type = assay_name, assay_name = "counts", ...)
x |
a
|
... |
additional arguments
|
detection |
Detection threshold for absence/presence. Either an
absolute value compared directly to the values of |
include_lowest |
logical scalar: Should the lower boundary of the
detection and prevalence cutoffs be included? (default: |
sort |
logical scalar: Should the result be sorted by prevalence?
(default: |
na.rm |
logical scalar: Should NA values be omitted when calculating
prevalence? (default: |
assay.type |
A single character value for selecting the
|
assay_name |
a single |
as_relative |
logical scalar: Should the detection threshold be applied
on compositional (relative) abundances? (default: |
rank |
a single character defining a taxonomic rank. Must be a value of
|
prevalence |
Prevalence threshold (in 0 to 1). The
required prevalence is strictly greater by default. To include the
limit, set |
getPrevalence
calculates the relative frequency of samples that exceed
the detection threshold. For SummarizedExperiment
objects, the
prevalence is calculated for the selected taxonomic rank, otherwise for the
rows. The absolute population prevalence can be obtained by multiplying the
prevalence by the number of samples (ncol(x)
). If as_relative =
FALSE
the relative frequency (between 0 and 1) is used to check against the
detection
threshold.
The core abundance index from getPrevalentAbundance
gives the relative
proportion of the core species (in between 0 and 1). The core taxa are
defined as those that exceed the given population prevalence threshold at the
given detection level as set for getPrevalentFeatures
.
subsetPrevalentFeatures
and subsetRareFeatures
return a subset of x
.
The subset includes the most prevalent or rare taxa that are calculated with
getPrevalentFeatures
or getRareFeatures
respectively.
getPrevalentFeatures
returns taxa that are more prevalent with the
given detection threshold for the selected taxonomic rank.
getRareFeatures
returns complement of getPrevalentTaxa
.
subsetPrevalentFeatures
and subsetRareFeatures
return subset of x
.
All other functions return a named vectors:
getPrevalence
returns a numeric
vector with the
names being set to either the row names of x
or the names after
agglomeration.
getPrevalentAbundance
returns a numeric
vector with
the names corresponding to the column name of x
and include the
joint abundance of prevalent taxa.
getPrevalentTaxa
and getRareFeatures
return a
character
vector with only the names exceeding the threshold set
by prevalence
, if the rownames
of x
is set.
Otherwise an integer
vector is returned matching the rows in
x
.
Leo Lahti
For getPrevalentAbundance
: Leo Lahti and Tuomas Borman.
Contact: microbiome.github.io
A Salonen et al. The adult intestinal core microbiota is determined by analysis depth and health status. Clinical Microbiology and Infection 18(S4):16 20, 2012. To cite the R package, see citation('mia')
agglomerateByRank
,
getTopTaxa
data(GlobalPatterns)
tse <- GlobalPatterns
# Get prevalence estimates for individual ASV/OTU
prevalence.frequency <- getPrevalence(tse,
detection = 0,
sort = TRUE,
as_relative = TRUE)
head(prevalence.frequency)
# Get prevalence estimates for phylums
# - the getPrevalence function itself always returns population frequencies
prevalence.frequency <- getPrevalence(tse,
rank = "Phylum",
detection = 0,
sort = TRUE,
as_relative = TRUE)
head(prevalence.frequency)
# - to obtain population counts, multiply frequencies with the sample size,
# which answers the question "In how many samples is this phylum detectable"
prevalence.count <- prevalence.frequency * ncol(tse)
head(prevalence.count)
# Detection threshold 1 (strictly greater by default);
# Note that the data (GlobalPatterns) is here in absolute counts
# (and not compositional, relative abundances)
# Prevalence threshold 50 percent (strictly greater by default)
prevalent <- getPrevalentFeatures(tse,
rank = "Phylum",
detection = 10,
prevalence = 50/100,
as_relative = FALSE)
head(prevalent)
# Gets a subset of object that includes prevalent taxa
altExp(tse, "prevalent") <- subsetByPrevalentFeatures(tse,
rank = "Family",
detection = 0.001,
prevalence = 0.55,
as_relative = TRUE)
altExp(tse, "prevalent")
# getRareFeatures returns the inverse
rare <- getRareFeatures(tse,
rank = "Phylum",
detection = 1/100,
prevalence = 50/100,
as_relative = TRUE)
head(rare)
# Gets a subset of object that includes rare taxa
altExp(tse, "rare") <- subsetByRareFeatures(tse,
rank = "Class",
detection = 0.001,
prevalence = 0.001,
as_relative = TRUE)
altExp(tse, "rare")
# Names of both experiments, prevalent and rare, can be found from slot altExpNames
tse
data(esophagus)
getPrevalentAbundance(esophagus, assay.type = "counts")
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.