estimateDominance {mia}R Documentation

Estimate dominance measures

Description

This function calculates community dominance indices. This includes the ‘Absolute’, ‘Berger-Parker’, ‘Core abundance’, ‘Gini’, ‘McNaughton’s’, ‘Relative’, and ‘Simpson's’ indices.

Usage

estimateDominance(
  x,
  abund_values = "counts",
  index = c("absolute", "dbp", "core_abundance", "gini", "dmn", "relative",
    "simpson_lambda"),
  ntaxa = 1,
  aggregate = TRUE,
  name = index,
  ...,
  BPPARAM = SerialParam()
)

## S4 method for signature 'SummarizedExperiment'
estimateDominance(
  x,
  abund_values = "counts",
  index = c("absolute", "dbp", "core_abundance", "gini", "dmn", "relative",
    "simpson_lambda"),
  ntaxa = 1,
  aggregate = TRUE,
  name = index,
  ...,
  BPPARAM = SerialParam()
)

Arguments

x

a SummarizedExperiment object

abund_values

A single character value for selecting the assay to calculate the sample-wise estimates.

index

a character vector, specifying the indices to be calculated.

ntaxa

Optional and only used for the Absolute and Relative dominance indices: The n-th position of the dominant taxa to consider (default: ntaxa = 1). Disregarded for the indices “dbp”, “core_abundance”, “Gini”, “dmn”, and “Simpson”.

aggregate

Optional and only used for the Absolute, dbp, Relative, and dmn dominance indices: Aggregate the values for top members selected by ntaxa or not. If TRUE, then the sum of relative abundances is returned. Otherwise the relative abundance is returned for the single taxa with the indicated rank (default: aggregate = TRUE). Disregarded for the indices “core_abundance”, “gini”, “dmn”, and “simpson”.

name

A name for the column(s) of the colData where the calculated Dominance indices should be stored in.

...

additional arguments currently not used.

BPPARAM

A BiocParallelParam object specifying whether calculation of estimates should be parallelized. (Currently not used)

Details

A dominance index quantifies the dominance of one or few species in a community. Greater values indicate higher dominance.

Dominance indices are in general negatively correlated with alpha diversity indices (species richness, evenness, diversity, rarity). More dominant communities are less diverse.

estimateDominance calculates the following community dominance indices:

Value

x with additional colData named *name*

Author(s)

Leo Lahti and Tuomas Borman. Contact: microbiome.github.io

References

Berger WH & Parker FL (1970) Diversity of Planktonic Foraminifera in Deep-Sea Sediments. Science 168(3937):1345-1347. doi: 10.1126/science.168.3937.1345

Gini C (1921) Measurement of Inequality of Incomes. The Economic Journal 31(121): 124-126. doi: 10.2307/2223319

McNaughton, SJ and Wolf LL. (1970). Dominance and the niche in ecological systems. Science 167:13, 1–139

Simpson EH (1949) Measurement of Diversity. Nature 163(688). doi: 10.1038/163688a0

See Also

Examples

data(esophagus)

# Calculates Simpson's lambda (can be used as a dominance index)
esophagus <- estimateDominance(esophagus, index="simpson_lambda")

# Shows all indices
colData(esophagus)

# Indices must be written correctly (e.g. dbp, not dbp), otherwise an error
# gets thrown
## Not run: esophagus <- estimateDominance(esophagus, index="dbp")
# Calculates dbp and Core Abundance indices
esophagus <- estimateDominance(esophagus, index=c("dbp", "core_abundance"))
# Shows all indices
colData(esophagus)
# Shows dbp index
colData(esophagus)$dbp
# Deletes dbp index
colData(esophagus)$dbp <- NULL
# Shows all indices, dbp is deleted
colData(esophagus)
# Deletes all indices
colData(esophagus) <- NULL

# Calculates all indices
esophagus <- estimateDominance(esophagus)
# Shows all indices
colData(esophagus)
# Deletes all indices
colData(esophagus) <- NULL

# Calculates all indices with explicitly specified names
esophagus <- estimateDominance(esophagus,
    index = c("dbp", "dmn", "absolute", "relative",
              "simpson_lambda", "core_abundance", "gini"),
    name  = c("BergerParker", "McNaughton", "Absolute", "Relative",
              "SimpsonLambda", "CoreAbundance", "Gini")
)
# Shows all indices
colData(esophagus)


[Package mia version 1.2.2 Index]