Introduction

We recommend switching from phyloseq to the TreeSummarizedExperiment based methods described in the Orchestrating Microbiome Analysis (OMA) online book, which is where method development now takes place. The mia package provides the analysis methods for that container.

To ease the transition, most functions in this package accept a TreeSummarizedExperiment (or any other SummarizedExperiment-derived object) wherever they accept a phyloseq object. This vignette walks through the standard operations on a TreeSummarizedExperiment. The companion vignette, vignette("vignette", package="microbiome"), gives the general introduction to the package.

Almost everything is unchanged: the same function called on a TreeSummarizedExperiment returns the same value it would return for the equivalent phyloseq object. There are two things to know before starting, both covered in their own sections below:

transform() stores its result as a new named assay rather than overwriting the counts, so the result is read back with abundances(x, assay.type=).
The functions that manipulate the taxonomy table remain phyloseq only.

library(microbiome)

A note on masking

Attaching mia also attaches its dependencies, two of which export functions that share a name with a microbiome function. Since they are attached later, they take precedence:

IRanges::transform() masks microbiome::transform()
Biostrings::coverage() masks microbiome::coverage()

Neither raises an error when called on a microbiome object, so the symptom is a silently wrong result rather than a failure. Call these two with their package prefix whenever mia is attached, as this vignette does throughout:

microbiome::transform(tse, "compositional")
microbiome::coverage(tse)

The rest of the package is unaffected.

Getting a TreeSummarizedExperiment

mia converts an existing phyloseq object. Here we use the dietswap data set (O’Keefe et al. 2015) that ships with this package.

library(mia)

data(dietswap)
tse <- convertFromPhyloseq(dietswap)
tse

## class: TreeSummarizedExperiment 
## dim: 130 222 
## metadata(0):
## assays(1): counts
## rownames(130): Actinomycetaceae Aerococcus ... Xanthomonadaceae
##   Yersinia et rel.
## rowData names(3): Phylum Family Genus
## colnames(222): Sample-1 Sample-2 ... Sample-221 Sample-222
## colData names(8): subject sex ... timepoint.within.group bmi_group
## reducedDimNames(0):
## mainExpName: NULL
## altExpNames(0):
## rowLinks: NULL
## rowTree: NULL
## colLinks: NULL
## colTree: NULL

mia also reads the common file formats directly, with importBiom(), importMothur(), importQIIME2() and others, so a conversion step is not needed for new analyses. See the OMA book for those.

Accessors

abundances(), meta() and taxa() are the three accessors the rest of the package is built on. They work the same way on both containers, returning the abundance matrix, the sample metadata and the feature names.

a <- abundances(tse)
dim(a)

## [1] 130 222

a[seq_len(4), seq_len(3)]

##                  Sample-1 Sample-2 Sample-3
## Actinomycetaceae        0        1        0
## Aerococcus              0        0        0
## Aeromonas               0        0        0
## Akkermansia            18       97       67

head(meta(tse), 3)

##          subject  sex nationality group   sample timepoint
## Sample-1     byn male         AAM    DI Sample-1         4
## Sample-2     nms male         AFR    HE Sample-2         2
## Sample-3     olt male         AFR    HE Sample-3         2
##          timepoint.within.group  bmi_group
## Sample-1                      1      obese
## Sample-2                      1       lean
## Sample-3                      1 overweight

head(taxa(tse), 5)

## [1] "Actinomycetaceae"             "Aerococcus"                  
## [3] "Aeromonas"                    "Akkermansia"                 
## [5] "Alcaligenes faecalis et rel."

For a TreeSummarizedExperiment these read the assay, colData and rownames respectively. The number of features and samples comes from the object itself.

c(features=nrow(tse), samples=ncol(tse))

## features  samples 
##      130      222

Choosing the assay

A SummarizedExperiment can hold several assays at once. abundances() takes an assay.type argument to pick between them. It defaults to the counts assay when one is present, and otherwise to the first assay.

assayNames(tse)

## [1] "counts"

identical(abundances(tse), abundances(tse, assay.type="counts"))

## [1] TRUE

Asking for an assay that does not exist is an error rather than a silent fallback:

abundances(tse, assay.type="nonexistent")

## Error in `.se_assay()`:
## ! assay.type 'nonexistent' not found. Available assays: counts

Transformations

This is the one place where the result differs from the phyloseq version, and it is deliberate. For a phyloseq object, transform() overwrites the abundance table. For a SummarizedExperiment, transforming is additive: the result is stored as a new assay named after the transformation, and the original counts are left untouched. This follows the mia::transformAssay() convention.

tse <- microbiome::transform(tse, "compositional")
assayNames(tse)

## [1] "counts"        "compositional"

The transformed values are read back by naming the assay:

abundances(tse, assay.type="compositional")[seq_len(4), seq_len(3)]

##                     Sample-1     Sample-2    Sample-3
## Actinomycetaceae 0.000000000 4.742483e-05 0.000000000
## Aerococcus       0.000000000 0.000000e+00 0.000000000
## Aeromonas        0.000000000 0.000000e+00 0.000000000
## Akkermansia      0.002127911 4.600209e-03 0.002319704

The counts are still there, unchanged:

abundances(tse)[seq_len(4), seq_len(3)]

##                  Sample-1 Sample-2 Sample-3
## Actinomycetaceae        0        1        0
## Aerococcus              0        0        0
## Aeromonas               0        0        0
## Akkermansia            18       97       67

Use the name argument to control the assay name, for instance to follow the mia naming convention:

tse <- microbiome::transform(tse, "compositional", name="relabundance")
assayNames(tse)

## [1] "counts"        "compositional" "relabundance"

Transformations chain, since each one adds an assay:

tse <- microbiome::transform(tse, "clr")
tse <- microbiome::transform(tse, "Z")
assayNames(tse)

## [1] "counts"        "compositional" "relabundance"  "clr"          
## [5] "Z"

When you only need the transformed matrix and not a modified object, pass transform to abundances() directly. This never round-trips through the object and works identically for both containers:

clr <- abundances(tse, transform="clr")
clr[seq_len(4), seq_len(3)]

##                    Sample-1  Sample-2  Sample-3
## Actinomycetaceae -2.6385281 -2.016809 -3.114983
## Aerococcus       -2.6385281 -3.115421 -3.114983
## Aeromonas        -2.6385281 -3.115421 -3.114983
## Akkermansia       0.9723898  2.157578  1.790291

Alpha diversity and related indices

These take the object and return a per-sample data frame or vector, as they do for phyloseq.

head(alpha(tse, index=c("shannon", "gini_simpson")), 3)

##          diversity_shannon diversity_gini_simpson
## Sample-1          2.940698              0.8677550
## Sample-2          2.822472              0.8765881
## Sample-3          2.407962              0.7669587

head(richness(tse), 3)

##          observed    chao1
## Sample-1      104 113.2308
## Sample-2      110 120.0000
## Sample-3      103 106.6000

head(microbiome::diversity(tse, index="shannon"), 3)

##           shannon
## Sample-1 2.940698
## Sample-2 2.822472
## Sample-3 2.407962

head(evenness(tse, index="pielou"), 3)

##             pielou
## Sample-1 0.6331719
## Sample-2 0.6004646
## Sample-3 0.5195476

head(dominance(tse, index="absolute"), 3)

##          absolute
## Sample-1     2774
## Sample-2     5121
## Sample-3    13271

head(rarity(tse, index="log_modulo_skewness"), 3)

##          log_modulo_skewness
## Sample-1            2.057691
## Sample-2            2.057552
## Sample-3            2.052324

The most abundant taxon per sample, and the overall top taxa:

head(dominant(tse), 3)

## [1] "Bacteroides vulgatus et rel."      "Prevotella melaninogenica et rel."
## [3] "Prevotella melaninogenica et rel."

top_taxa(tse, n=5)

## [1] "Prevotella melaninogenica et rel."  "Oscillospira guillermondii et rel."
## [3] "Bacteroides vulgatus et rel."       "Clostridium cellulosi et rel."     
## [5] "Prevotella oralis et rel."

Divergence from a reference sample, here the median across samples:

ref <- apply(abundances(tse), 1, median)
head(divergence(tse, ref, method="bray"), 3)

##  Sample-1  Sample-2  Sample-3 
## 0.5245336 0.5426288 0.6134272

Sample-wise summaries:

head(readcount(tse), 3)

## Sample-1 Sample-2 Sample-3 
##     8459    21086    28883

head(microbiome::coverage(tse), 3)

## Sample-1 Sample-2 Sample-3 
##        4        3        2

head(low_abundance(tse), 3)

##   Sample-1   Sample-2   Sample-3 
## 0.02896323 0.02997249 0.03396462

Core microbiota and prevalence

head(prevalence(tse, detection=1, sort=TRUE), 5)

##                      Vibrio       Uncultured Mollicutes 
##                           1                           1 
## Uncultured Clostridiales II  Uncultured Clostridiales I 
##                           1                           1 
##          Tannerella et rel. 
##                           1

core_members() lists the taxa passing a detection and prevalence threshold, and core() returns the object filtered down to them.

core_members(tse, detection=1, prevalence=50/100)

##  [1] "Akkermansia"                          
##  [2] "Alcaligenes faecalis et rel."         
##  [3] "Allistipes et rel."                   
##  [4] "Anaerostipes caccae et rel."          
##  [5] "Anaerotruncus colihominis et rel."    
##  [6] "Anaerovorax odorimutans et rel."      
##  [7] "Bacteroides fragilis et rel."         
##  [8] "Bacteroides intestinalis et rel."     
##  [9] "Bacteroides ovatus et rel."           
## [10] "Bacteroides plebeius et rel."         
## [11] "Bacteroides splachnicus et rel."      
## [12] "Bacteroides stercoris et rel."        
## [13] "Bacteroides uniformis et rel."        
## [14] "Bacteroides vulgatus et rel."         
## [15] "Bifidobacterium"                      
## [16] "Bryantella formatexigens et rel."     
## [17] "Bulleidia moorei et rel."             
## [18] "Butyrivibrio crossotus et rel."       
## [19] "Campylobacter"                        
## [20] "Clostridium (sensu stricto)"          
## [21] "Clostridium cellulosi et rel."        
## [22] "Clostridium colinum et rel."          
## [23] "Clostridium difficile et rel."        
## [24] "Clostridium leptum et rel."           
## [25] "Clostridium nexile et rel."           
## [26] "Clostridium orbiscindens et rel."     
## [27] "Clostridium ramosum et rel."          
## [28] "Clostridium sphenoides et rel."       
## [29] "Clostridium stercorarium et rel."     
## [30] "Clostridium symbiosum et rel."        
## [31] "Collinsella"                          
## [32] "Coprobacillus catenaformis et rel."   
## [33] "Coprococcus eutactus et rel."         
## [34] "Desulfovibrio et rel."                
## [35] "Dialister"                            
## [36] "Dorea formicigenerans et rel."        
## [37] "Eggerthella lenta et rel."            
## [38] "Enterobacter aerogenes et rel."       
## [39] "Enterococcus"                         
## [40] "Escherichia coli et rel."             
## [41] "Eubacterium biforme et rel."          
## [42] "Eubacterium cylindroides et rel."     
## [43] "Eubacterium hallii et rel."           
## [44] "Eubacterium rectale et rel."          
## [45] "Eubacterium siraeum et rel."          
## [46] "Eubacterium ventriosum et rel."       
## [47] "Faecalibacterium prausnitzii et rel." 
## [48] "Fusobacteria"                         
## [49] "Helicobacter"                         
## [50] "Klebisiella pneumoniae et rel."       
## [51] "Lachnobacillus bovis et rel."         
## [52] "Lachnospira pectinoschiza et rel."    
## [53] "Lactobacillus gasseri et rel."        
## [54] "Lactobacillus plantarum et rel."      
## [55] "Megasphaera elsdenii et rel."         
## [56] "Mitsuokella multiacida et rel."       
## [57] "Oceanospirillum"                      
## [58] "Oscillospira guillermondii et rel."   
## [59] "Outgrouping clostridium cluster XIVa" 
## [60] "Oxalobacter formigenes et rel."       
## [61] "Papillibacter cinnamivorans et rel."  
## [62] "Parabacteroides distasonis et rel."   
## [63] "Peptococcus niger et rel."            
## [64] "Peptostreptococcus micros et rel."    
## [65] "Phascolarctobacterium faecium et rel."
## [66] "Prevotella melaninogenica et rel."    
## [67] "Prevotella oralis et rel."            
## [68] "Prevotella ruminicola et rel."        
## [69] "Prevotella tannerae et rel."          
## [70] "Proteus et rel."                      
## [71] "Roseburia intestinalis et rel."       
## [72] "Ruminococcus bromii et rel."          
## [73] "Ruminococcus callidus et rel."        
## [74] "Ruminococcus gnavus et rel."          
## [75] "Ruminococcus lactaris et rel."        
## [76] "Ruminococcus obeum et rel."           
## [77] "Sporobacter termitidis et rel."       
## [78] "Streptococcus bovis et rel."          
## [79] "Streptococcus intermedius et rel."    
## [80] "Streptococcus mitis et rel."          
## [81] "Subdoligranulum variable at rel."     
## [82] "Sutterella wadsworthia et rel."       
## [83] "Tannerella et rel."                   
## [84] "Uncultured Clostridiales I"           
## [85] "Uncultured Clostridiales II"          
## [86] "Uncultured Mollicutes"                
## [87] "Veillonella"                          
## [88] "Vibrio"                               
## [89] "Yersinia et rel."

tse.core <- core(tse, detection=1, prevalence=50/100)
c(features=nrow(tse.core), samples=ncol(tse.core))

## features  samples 
##       89      222

Subsetting preserves the class, so the result is a TreeSummarizedExperiment that can be fed back into any of these functions:

class(tse.core)

## [1] "TreeSummarizedExperiment"
## attr(,"package")
## [1] "TreeSummarizedExperiment"

rare() and rare_members() are the complements, and core_abundance() gives the fraction of the community the core accounts for:

head(core_abundance(tse, detection=1, prevalence=50/100), 3)

## [1] NA

length(rare_members(tse, detection=1, prevalence=50/100))

## [1] 41

Subsetting

remove_taxa() and remove_samples() drop features and samples by name, again returning the same class they were given.

drop <- taxa(tse)[seq_len(5)]
tse.sub <- remove_taxa(drop, tse)
c(before=nrow(tse), after=nrow(tse.sub))

## before  after 
##    130    125

drop.s <- colnames(tse)[seq_len(5)]
tse.sub <- remove_samples(drop.s, tse)
c(before=ncol(tse), after=ncol(tse.sub))

## before  after 
##    222    217

Visualization

The plotting functions that go through the accessors take a TreeSummarizedExperiment and return a ggplot object as usual.

plot_core(tse, prevalences=seq(.1, 1, .2), detections=3)

Core microbiota heatmap.

boxplot_abundance(tse, x="nationality", y=taxa(tse)[[1]])

Abundance of a single taxon across groups.

boxplot_alpha(tse, x_var="nationality", index="shannon")

Alpha diversity across groups.

plot_density(tse, variable=taxa(tse)[[1]])

Abundance density for a single taxon.

plot_atlas() and plot_tipping() accept a TreeSummarizedExperiment in the same way.

Time series

Given a time field in the sample metadata, the time series helpers work as they do for phyloseq. time_normalize() writes back into the object and so returns a TreeSummarizedExperiment.

tse.t <- tse
tse.t$time <- tse.t$timepoint

length(timesplit(tse.t))

## [1] 4

head(time_sort(tse.t)$time, 3)

## NULL

class(time_normalize(tse.t))

## [1] "TreeSummarizedExperiment"
## attr(,"package")
## [1] "TreeSummarizedExperiment"

baseline() and collapse_replicates() likewise accept and return the same container.

What is not supported

The functions that build or manipulate the taxonomy table are phyloseq only, and error on a TreeSummarizedExperiment:

aggregate_taxa(), aggregate_rare()
plot_composition()
map_levels()
psmelt2(), otu_tibble(), tax_tibble(), sample_tibble()
add_besthit(), add_refseq()
the read_*() family

mia covers this ground for TreeSummarizedExperiment objects, and does it better: use mia::agglomerateByRank() in place of aggregate_taxa(), mia::meltSE() in place of psmelt2(), and the miaViz package for composition plots. The OMA book documents these.

Session info

sessionInfo()

## R version 4.6.1 (2026-06-24)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 26.04.1 LTS
## 
## Matrix products: default
## BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
## LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.32.so;  LAPACK version 3.12.0
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
##  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
##  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
##  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
## 
## time zone: Etc/UTC
## tzcode source: system (glibc)
## 
## attached base packages:
## [1] stats4    stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] mia_1.21.11                     TreeSummarizedExperiment_2.21.0
##  [3] Biostrings_2.81.9               XVector_0.53.0                 
##  [5] SingleCellExperiment_1.35.2     MultiAssayExperiment_1.39.1    
##  [7] SummarizedExperiment_1.43.0     Biobase_2.73.2                 
##  [9] GenomicRanges_1.65.4            Seqinfo_1.3.2                  
## [11] IRanges_2.47.5                  S4Vectors_0.51.10              
## [13] BiocGenerics_0.59.12            generics_0.1.4                 
## [15] MatrixGenerics_1.25.0           matrixStats_1.5.0              
## [17] microbiome_1.35.1               ggplot2_4.0.3                  
## [19] phyloseq_1.57.0                 BiocStyle_2.41.0               
## 
## loaded via a namespace (and not attached):
##   [1] DBI_1.3.0                   gridExtra_2.3.1            
##   [3] permute_0.9-10              rlang_1.3.0                
##   [5] magrittr_2.0.5              ade4_1.7-24                
##   [7] scater_1.41.2               otel_0.2.0                 
##   [9] compiler_4.6.1              mgcv_1.9-4                 
##  [11] DelayedMatrixStats_1.35.0   vctrs_0.7.3                
##  [13] reshape2_1.4.5              stringr_1.6.0              
##  [15] pkgconfig_2.0.3             crayon_1.5.3               
##  [17] fastmap_1.2.0               labeling_0.4.3             
##  [19] scuttle_1.23.2              rmarkdown_2.32             
##  [21] ggbeeswarm_0.7.3            DirichletMultinomial_1.55.1
##  [23] purrr_1.2.2                 xfun_0.61                  
##  [25] bluster_1.23.1              cachem_1.1.0               
##  [27] beachmat_2.29.3             jsonlite_2.0.0             
##  [29] biomformat_1.41.0           DelayedArray_0.39.8        
##  [31] BiocParallel_1.47.0         irlba_2.4.1                
##  [33] parallel_4.6.1              cluster_2.1.8.3            
##  [35] R6_2.6.1                    bslib_0.12.0               
##  [37] stringi_1.8.9               RColorBrewer_1.1-3         
##  [39] jquerylib_0.1.4             iterators_1.0.14           
##  [41] Rcpp_1.1.2                  knitr_1.52                 
##  [43] DECIPHER_3.9.4              Matrix_1.7-6               
##  [45] splines_4.6.1               igraph_2.3.4               
##  [47] tidyselect_1.2.1            abind_1.4-8                
##  [49] yaml_2.3.12                 viridis_0.6.5              
##  [51] vegan_2.7-6                 codetools_0.2-20           
##  [53] lattice_0.23-1              tibble_3.3.1               
##  [55] plyr_1.8.9                  withr_3.0.3                
##  [57] treeio_1.37.1               S7_0.2.2                   
##  [59] Rtsne_0.17                  evaluate_1.0.5             
##  [61] survival_3.8-12             pillar_1.11.1              
##  [63] BiocManager_1.30.27         foreach_1.5.2              
##  [65] sparseMatrixStats_1.25.0    scales_1.4.0               
##  [67] tidytree_0.4.8              glue_1.8.1                 
##  [69] lazyeval_0.2.3              maketools_1.3.2            
##  [71] tools_4.6.1                 data.table_1.18.6.1        
##  [73] BiocNeighbors_2.7.3         sys_3.4.3                  
##  [75] ScaledMatrix_1.21.0         buildtools_1.0.0           
##  [77] fs_2.1.0                    grid_4.6.1                 
##  [79] tidyr_1.3.2                 ape_5.8-1                  
##  [81] nlme_3.1-171                beeswarm_0.4.0             
##  [83] BiocSingular_1.29.2         vipor_0.4.7                
##  [85] cli_3.6.6                   rsvd_1.0.5                 
##  [87] rappdirs_0.3.4              S4Arrays_1.13.2            
##  [89] viridisLite_0.4.3           dplyr_1.2.1                
##  [91] gtable_0.3.6                yulab.utils_0.2.5          
##  [93] sass_0.4.10                 digest_0.6.39              
##  [95] SparseArray_1.13.4          ggrepel_0.9.8              
##  [97] decontam_1.33.0             farver_2.1.2               
##  [99] multtest_2.69.0             htmltools_0.5.9            
## [101] lifecycle_1.0.5             MASS_7.3-66

References

O’Keefe, SJD, JV Li, L Lahti, et al. 2015. “Fat, Fiber and Cancer Risk in African, Americans and Rural Africans.” Nature Communications 6 (April): 6342. https://doi.org/10.1038/ncomms7342.

Using microbiome with TreeSummarizedExperiment