MicrobiomeStat is an R package for statistical analysis and visualization of microbiome and multi-omics data. It supports cross-sectional, paired-sample, and longitudinal study designs through a unified interface, with built-in automated report generation.

• Unified interface -- test_alpha(), test_beta(), test_taxa() and their plot_*() counterparts automatically dispatch to the correct method based on study design.
• LinDA -- Linear models for differential abundance analysis of compositional data (Zhou et al. 2022), with an enhanced linda2() variant supporting detection-depth weighting and tree-guided smoothing.
• Longitudinal analysis -- Trend tests, volatility assessment, spaghetti plots, and per-time-point comparisons for time-series microbiome data.
• Automated reports -- mStat_generate_report_single(), _pair(), and _long() produce publication-ready PDF/HTML reports from a single function call.
• Broad data import -- Native support for phyloseq, QIIME2, DADA2, mothur, BIOM, DESeq2, edgeR, SummarizedExperiment, and MultiAssayExperiment objects.
• 120+ functions covering alpha diversity, beta diversity, taxonomic composition, and their change/trend/volatility variants across all three study designs.

Install from CRAN:

install.packages("MicrobiomeStat")

Or install the development version from GitHub for the full feature set:

# install.packages("devtools")
devtools::install_github("cafferychen777/MicrobiomeStat")

library(MicrobiomeStat)

# Load example data
data(peerj32.obj)

# Alpha diversity
alpha.obj <- mStat_calculate_alpha_diversity(
  peerj32.obj$feature.tab,
  alpha.name = c("shannon", "observed_species")
)

generate_alpha_boxplot_single(
  data.obj = peerj32.obj,
  alpha.obj = alpha.obj,
  alpha.name = "shannon",
  subject.var = "subject",
  time.var = "time",
  t.level = "2",
  group.var = "group",
  strata.var = "sex"
)

# Beta diversity ordination
generate_beta_ordination_single(
  data.obj = peerj32.obj,
  dist.obj = NULL,
  pc.obj = NULL,
  subject.var = "subject",
  time.var = "time",
  t.level = "2",
  group.var = "group",
  strata.var = "sex",
  dist.name = "BC"
)

# Differential abundance with LinDA
generate_taxa_test_single(
  data.obj = peerj32.obj,
  subject.var = "subject",
  time.var = "time",
  t.level = "2",
  group.var = "group",
  adj.vars = "sex",
  feature.level = "Genus",
  feature.dat.type = "count"
)

# From phyloseq
data.obj <- mStat_convert_phyloseq_to_data_obj(phyloseq_obj)

# From QIIME2
data.obj <- mStat_import_qiime2_as_data_obj(
  otu.qza = "table.qza",
  taxa.qza = "taxonomy.qza",
  sam.tab = "metadata.txt"
)

# From DADA2
data.obj <- mStat_import_dada2_as_data_obj(
  seq.tab = seqtab, tax.tab = taxtab, sam.tab = metadata
)

See ?mStat_import_biom_as_data_obj, ?mStat_import_mothur_as_data_obj, and the Bioconductor converters for other formats.

• Tutorials: www.microbiomestat.wiki
• Function reference: cafferychen777.github.io/MicrobiomeStat
• Vignettes: vignette("MicrobiomeStat") and vignette("longitudinal-analysis")

If you use MicrobiomeStat in your research, please cite:

Zhou, H., He, K., Chen, J., & Zhang, X. (2022). LinDA: linear models for differential abundance analysis of microbiome compositional data. Genome Biology, 23(1), 95. https://doi.org/10.1186/s13059-022-02655-5

citation("MicrobiomeStat")

• GitHub Issues -- bug reports and feature requests
• Discord -- questions and discussion