In cmmr/rbiom: Read/Write, Analyze, and Visualize 'BIOM' Data

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  R.options = list(
    pillar.print_min = 5,
    pillar.print_max = 5 ))

Create

The general purpose as_rbiom() function can convert most data types into an rbiom object - see importing for details. Here we'll import a dataset from a BIOM file.

library(rbiom)

file <- system.file(package = "rbiom", "extdata", "hmp50.bz2")
biom <- as_rbiom(file)

biom

Inspect

The rbiom object has many helpful accessors.

| Accessor | Content | |--------------------------|--------------------------------------------------------| | $counts | Abundance of each OTU in each sample. | | $metadata | Sample mappings to metadata (treatment, patient, etc). | | $taxonomy | OTU mappings to taxonomic ranks (genus, phylum, etc). | | $otus, $n_otus | OTU names. | | $samples, $n_samples | Sample names. | | $fields, $n_fields | Metadata field names. | | $ranks, $n_ranks | Taxonomic rank names. | | $tree, $sequences | Phylogenetic tree / sequences for the OTUs, or NULL. | | $id, $comment | Arbitrary strings for describing the dataset, or "". | | $depth | Rarefaction depth, or NULL if unrarefied. |

biom$counts[1:4,1:8] %>% as.matrix()

biom$fields

# Use pull() to automatically setNames().
pull(biom, 'Age') %>% head()

pull(), sample_sums(), taxa_matrix(), taxa_means(), taxa_sums()

Clone

Rbiom objects are passed by reference. The common <- assignment operator creates a second reference to the same object - it does not create a second object. To create a copy of an rbiom object, use the object's $clone() method.

wzxhzdk:3

wzxhzdk:4

Modify

There are seven components of an rbiom object which you can modify directly. Assigning new values to these components will trigger validation checks and inter-component OTU/sample synchronization. See Working with rbiom Objects for additional details.

| Component | What can be assigned. | |--------------|--------------------------------------------------------| | $counts | matrix of abundances; OTUs (rows) by samples (columns) | | $metadata | data.frame with '.sample' as the first column | | $taxonomy | data.frame with '.otu' as the first column | | $tree | phylo object with the phylogenetic tree for OTUs | | $sequences | character vector of reference sequences for OTUs | | $id | string with a title for the dataset | | $comment | string with additional information about the dataset |

Rarefy Counts

A common way to normalize microbiome count data is to rarefy it. This process drops samples with too few observations, and randomly removes observations from the remaining samples, so that all samples have the same "rarefaction depth".

sample_sums(biom) %>% head()

biom <- rarefy(biom)

sample_sums(biom) %>% head()

Add Metadata

Additional sample metadata columns can be added to biom$metadata (a tibble data.frame). The first column, '.sample', is used by rbiom to link sample metadata to samples in the abundance table.

biom$metadata$group <- sample(c('A', 'B'), biom$n_samples, TRUE)
biom %<>% mutate(Obese = BMI >= 30, Sex = NULL)
biom %<>% rename('Years Old' = "Age")
biom$metadata

mutate(), rename()

Subset Samples

Removing samples from the metadata will remove those samples from the entire rbiom object.

biom %<>% subset(`Body Site` == "Anterior nares")
biom$metadata
biom

subset(), slice()

Compute

Functions that end in _table or _matrix return calculation results for use outside of rbiom. The _table suffix indicates the returned object will be a tibble data.frame with one computed value per row. Alternatively, _matrix will return a base R matrix.

taxa_table(biom, rank = "Phylum")

taxa_matrix(biom, rank = "Phylum")[1:4, 1:8]

Functions for taxa abundance (like the above taxa_table()) are prefixed by taxa_. Similarly, adiv_ is used for alpha diversity and bdiv_ for beta diversity.

adiv_matrix(), adiv_table(), bdiv_table(), bdiv_ord_table(), taxa_table(), taxa_matrix()

Visualize

The plotting functions in rbiom make it easy to produce informative visualizations about alpha diversity, beta diversity, and taxa abundance, and explore associations between those metrics and sample metadata.

See the Plot Types article for an overview of all the different plot options, and Mapping Metadata to Aesthetics for guidance on using colors, shapes, and patterns to represent metadata values.

For example, to display an ordination colored by body site:

bdiv_ord_plot(biom = hmp50, color.by = "Body Site")

adiv_boxplot(), adiv_corrplot(), bdiv_boxplot(), bdiv_corrplot(), bdiv_heatmap(), bdiv_ord_plot(), rare_stacked(), rare_corrplot(), rare_multiplot(), taxa_stacked(), taxa_boxplot(), taxa_corrplot(), taxa_heatmap()

Analyze

Visualizations are an excellent way to observe differences between sample groups. When color.by is set, boxplots, corrplots, and ord_plots will include the following:

• p-values and brief methodology on the plot itself.
• $stats attribute with a detailed statistics table.
• $stats$code attribute with the R code used to generate the table.

The Statistics article reviews this topic in greater detail.

p <- adiv_boxplot(hmp50, color.by = "Body Site", facet.by = "Sex")
p
p$stats

p$stats$code

adiv_stats(), bdiv_stats(), distmat_stats(), taxa_stats()

Share

Plots are essentially ggplot objects, and ggplot2::ggsave() can save them as PNG, PDF, SVG, and many other image file types.

To save your rbiom object to a BIOM file, use write_biom(). When sharing a dataset with someone who is unfamiliar with accessing BIOM files, write_xlsx() is also an option.

A few additional resources to know about:

• R Markdown - reproducible analyses/reports.
• flexdashboard - interactive dashboards.
• atima - rbiom in website form.

ggplot2::ggsave(), patchwork::wrap_plots(), write_biom(), write_xlsx()

cmmr/rbiom documentation built on April 28, 2024, 6:38 a.m.