divergence {microbiome}R Documentation

Diversity within a Sample Group

Description

Quantify microbiota divergence (heterogeneity) within a given sample set.

Usage

divergence(x, method = "bray", coreset = NULL)

Arguments

x

phyloseq object

method

dissimilarity method: any method available via stats::cor or phyloseq::distance function. Note that some methods ("jsd" for instance) do not work with the group divergence.

coreset

phyloseq object; the samples to be used to define the centroid

Details

Microbiota divergence (heterogeneity / spread) within a given sample set can be quantified by the average sample dissimilarity or beta diversity. Taking average over all pairwise dissimilarities is sensitive to sample size and heavily biased as the similarity values are not independent. To reduce this bias, the dissimilarity of each sample against the group mean is calculated. This generates one value per sample. These can be compared between groups in order to compare differences in group homogeneity.

Note that this measure is still affected by sample size. Subsampling or bootstrapping can be applied to equalize sample sizes between comparisons.

The spearman mode is a simple indicator that returns average spearman correlation between samples of the input data and the overall group-wise average. The inverse of this measure (ie rho instead of 1-rho as in here) was used in Salonen et al. (2014) to quantify group homogeneity.

Value

Vector with dissimilarities; one for each sample, quantifying the dissimilarity of the sample from the group-level mean.

Author(s)

Contact: Leo Lahti microbiome-admin@googlegroups.com

References

The inter- and intra-individual homogeneity measures used in Salonen et al. ISME J. 8:2218-30, 2014 were obtained as 1 - beta where beta is the group diversity as quantified by the spearman method.

To cite this R package, see citation('microbiome')

See Also

the vegdist function from the vegan package provides many standard beta diversity measures

Examples

# Assess beta diversity among the African samples
# in a diet swap study (see \code{help(dietswap)} for references)
data(dietswap)
b <- divergence(subset_samples(dietswap, nationality == 'AFR'))

[Package microbiome version 1.2.1 Index]