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

```{css, echo=FALSE}

.sourceCode { display: none; } math[display="block"] { font-size: 22px; }

img { background-color: #fff; border: 2px gray solid; }

.center { margin-left: auto; margin-right: auto; }

.otu_matrix { width: 80%; border-width: 0; border-collapse: separate; padding: 10px; font-size: 14pt; margin-left: auto; margin-right: auto; }

.def_table { font-size: 12pt; margin-left: auto; margin-right: auto; } .def_table td { padding: 10px; }

.level2 { margin-top: 40px; }

.level1 { margin-top: 50px; }

```r
knitr::opts_chunk$set(
  collapse = TRUE,
  comment  = "#>",
  dev      = 'ragg_png' )

Introduction

The different UniFrac algorithms are listed below, along with examples for calculating them.

Input Data

json <- '{"id":"","comment":"","date":"2025-01-29T22:14:00Z","format":"1.0.0","type":"OTU table","format_url":"http://biom-format.org","generated_by":"rbiom 2.0.13","matrix_type":"sparse","matrix_element_type":"int","shape":[5,2],"phylogeny":"(((OTU_1:0.8,OTU_2:0.5):0.4,OTU_3:0.9):0.2,(OTU_4:0.7,OTU_5:0.3):0.6);","rows":{"1":{"id":"OTU_1"},"2":{"id":"OTU_2"},"3":{"id":"OTU_3"},"4":{"id":"OTU_4"},"5":{"id":"OTU_5"}},"columns":[{"id":"Sample_1"},{"id":"Sample_2"}],"data":[[2,0,9],[3,0,3],[4,0,3],[0,1,1],[1,1,4],[2,1,2],[3,1,8]]}'

biom <- rbiom::as_rbiom(json, underscores = TRUE)
mtx  <- t(as.matrix(biom))
phy  <- biom$tree

L <- phy$edge.length
A <- c(9,0,0,0,9,6,3,3)
B <- c(7,5,1,4,2,8,8,0)

• An OTU matrix with two samples and five OTUs.
• A phylogenetic tree for those five OTUs.

knitr::kable(mtx, format="html", table.attr='class="otu_matrix" cellspacing="0"', align='c')

par(xpd = NA)
plot(
  x          = phy, 
  direction  = 'downwards', 
  srt        = 90, 
  adj        = 0.5, 
  no.margin  = TRUE,
  underscore = TRUE,
  x.lim      = c(0.5, 5.5) )

ape::edgelabels(phy$edge.length, bg = 'white', frame = 'none', adj = -0.4)

Definitions

The branch indices (green circles) are used for ordering the $L$, $A$, and $B$ arrays. Values for $L$ are drawn from the input phylogenetic tree. Values for $A$ and $B$ are the total number of OTU observations descending from that branch; $A$ for Sample_1, and $B$ for Sample_2.

local({

  phy$edge.length <- c(1, 1, 1, 1, 2, 1, 2, 2)

  par(xpd = NA)
  plot(
    x               = phy, 
    direction       = 'downwards', 
    srt             = 90, 
    adj             = 0.5, 
    no.margin       = TRUE,
    underscore      = TRUE,
    x.lim           = c(.8, 6) )

  ape::edgelabels(1:8, frame = 'circle')

  ape::edgelabels(paste('A =', A), bg = 'white', frame = 'none', adj = c(-0.4, -1.2))
  ape::edgelabels(paste('B =', B), bg = 'white', frame = 'none', adj = c(-0.4,  0.0))
  ape::edgelabels(paste('L =', L), bg = 'white', frame = 'none', adj = c(-0.3,  1.2))
})

$n = 8$	Number of branches
$A = \{`r A`\}$	Branch **weights** for Sample_1.
$B = \{`r B`\}$	Branch **weights** for Sample_2.
$A_T = 15$	Total OTU counts for Sample_1.
$B_T = 15$	Total OTU counts for Sample_2.
$L = \{`r L`\}$	The branch **lengths**.

Unweighted

• Lozupone et al, 2005: Unweighted UniFrac
• R Package rbiom: bdiv_matrix(bdiv = "unifrac", weighted=FALSE)
• R Package phyloseq: UniFrac(weighted=FALSE)
• R Package abdiv: unweighted_unifrac()
• qiime2 qiime diversity beta-phylogenetic --p-metric unweighted_unifrac
• mothur: unifrac.unweighted()

First, transform A and B into presence (1) and absence (0) indicators.

\begin{align*} A &= \{`r A`\} \\ A' &= \{`r as.numeric(A > 0)`\} \end{align*}

\begin{align*} B &= \{`r B`\} \\ B' &= \{`r as.numeric(B > 0)`\} \end{align*}

Then apply the formula:

\begin{align} U &= \displaystyle \frac{\sum_{i = 1}^{n} L_i(|A'i - B'_i|)}{\sum{i = 1}^{n} L_i(max(A'_i,B'_i))} \ \ U &= \displaystyle \frac{L_1(|A'_1-B'_1|) + L_2(|A'_2-B'_2|) + \cdots + L_n(|A'_n-B'_n|)}{L_1(max(A'_1,B'_1)) + L_2(max(A'_2,B'_2)) + \cdots + L_n(max(A'_n,B'_n))} \ \ U &= \displaystyle \frac{0.2(|1-1|) + 0.4(|0-1|) + \cdots + 0.3(|1-0|)}{0.2(max(1,1)) + 0.4(max(0,1)) + \cdots + 0.3(max(1,0))} \ \ U &= \displaystyle \frac{0.2(0) + 0.4(1) + 0.8(1) + 0.5(1) + 0.9(0) + 0.6(0) + 0.7(0) + 0.3(1)}{0.2(1) + 0.4(1) + 0.8(1) + 0.5(1) + 0.9(1) + 0.6(1) + 0.7(1) + 0.3(1)} \ \ U &= \displaystyle \frac{0.4 + 0.8 + 0.5 + 0.3}{0.2 + 0.4 + 0.8 + 0.5 + 0.9 + 0.6 + 0.7 + 0.3} \ \ U &= \displaystyle \frac{2}{4.4} \ \ U &= 0.4545455 \end{align}

Weighted

• Lozupone et al, 2007: Raw Weighted UniFrac
• R Package rbiom: bdiv_matrix(bdiv = "unifrac", weighted=TRUE, normalized=FALSE)
• R Package phyloseq: UniFrac(weighted=TRUE, normalized=FALSE)
• R Package abdiv: weighted_unifrac()
• qiime2 qiime diversity beta-phylogenetic --p-metric weighted_unifrac

\begin{align} W &= \sum_{i = 1}^{n} L_i|\frac{A_i}{A_T} - \frac{B_i}{B_T}| \ \ W &= L_1|\frac{A_1}{A_T} - \frac{B_1}{B_T}| + L_2|\frac{A_2}{A_T} - \frac{B_2}{B_T}| + \cdots + L_n|\frac{A_n}{A_T} - \frac{B_n}{B_T}| \ \ W &= 0.2|\frac{9}{15} - \frac{7}{15}| + 0.4|\frac{0}{15} - \frac{5}{15}| + \cdots + 0.3|\frac{3}{15} - \frac{0}{15}| \ \ W &= 0.02\overline{6} + 0.1\overline{3} + 0.05\overline{3} + 0.1\overline{3} + 0.42 + 0.08 + 0.2\overline{3} + 0.06 \ \ W &= 1.14 \end{align}

Normalized

• Lozupone et al, 2007: Normalized Weighted UniFrac
• R Package rbiom: bdiv_matrix(bdiv = "unifrac", weighted=TRUE, normalized=TRUE)
• R Package phyloseq: UniFrac(weighted=TRUE, normalized=TRUE)
• R Package abdiv: weighted_normalized_unifrac()
• qiime2 qiime diversity beta-phylogenetic --p-metric weighted_normalized_unifrac
• mothur: unifrac.weighted()

\begin{align*}

N &= \displaystyle \frac {\sum_{i = 1}^{n} L_i|\frac{A_i}{A_T} - \frac{B_i}{B_T}|} {\sum_{i = 1}^{n} L_i(\frac{A_i}{A_T} + \frac{B_i}{B_T})} \ \

N &= \displaystyle \frac {L_1|\frac{A_1}{A_T} - \frac{B_1}{B_T}| + L_2|\frac{A_2}{A_T} - \frac{B_2}{B_T}| + \cdots + L_n|\frac{A_n}{A_T} - \frac{B_n}{B_T}|} {L_1(\frac{A_1}{A_T} + \frac{B_1}{B_T}) + L_2(\frac{A_2}{A_T} + \frac{B_2}{B_T}) + \cdots + L_n(\frac{A_n}{A_T} + \frac{B_n}{B_T})} \ \

N &= \displaystyle \frac {0.2|\frac{9}{15} - \frac{7}{15}| + 0.4|\frac{0}{15} - \frac{5}{15}| + \cdots + 0.3|\frac{3}{15} - \frac{0}{15}|} {0.2(\frac{9}{15} + \frac{7}{15}) + 0.4(\frac{0}{15} + \frac{5}{15}) + \cdots + 0.3(\frac{3}{15} + \frac{0}{15})} \ \

N &= \displaystyle \frac {0.02\overline{6} + 0.1\overline{3} + 0.05\overline{3} + 0.1\overline{3} + 0.42 + 0.08 + 0.2\overline{3} + 0.06} {0.21\overline{3} + 0.1\overline{3} + 0.05\overline{3} + 0.1\overline{3} + 0.66 + 0.56 + 0.51\overline{3} + 0.06} \ \

N &= \displaystyle \frac{1.14}{2.326667} \ \ N &= 0.4899713 \end{align*}

Generalized

• Chen et al. 2012: Generalized UniFrac
• R Package GUniFrac: GUniFrac(alpha=0.5)
• R Package abdiv: generalized_unifrac(alpha = 0.5)
• qiime2 qiime diversity beta-phylogenetic --p-metric generalized_unifrac -a 0.5

\begin{align*}

G &= \displaystyle \frac {\sum_{i = 1}^{n} L_i(\frac{A_i}{A_T} + \frac{B_i}{B_T})^{0.5} |\displaystyle \frac {\frac{A_i}{A_T} - \frac{B_i}{B_T}} {\frac{A_i}{A_T} + \frac{B_i}{B_T}} |} {\sum_{i = 1}^{n} L_i(\frac{A_i}{A_T} + \frac{B_i}{B_T})^{0.5}} \ \

G &= \displaystyle \frac { L_1(\frac{A_1}{A_T} + \frac{B_1}{B_T})^{0.5} |\displaystyle \frac {\frac{A_1}{A_T} - \frac{B_1}{B_T}} {\frac{A_1}{A_T} + \frac{B_1}{B_T}}| + \cdots + L_n(\frac{A_n}{A_T} + \frac{B_n}{B_T})^{0.5} |\displaystyle \frac {\frac{A_n}{A_T} - \frac{B_n}{B_T}} {\frac{A_n}{A_T} + \frac{B_n}{B_T}}| }{ L_1(\frac{A_1}{A_T} + \frac{B_1}{B_T})^{0.5} + \cdots + L_n(\frac{A_n}{A_T} + \frac{B_n}{B_T})^{0.5} }
\ \

G &= \displaystyle \frac { 0.2(\frac{9}{15} + \frac{7}{15})^{0.5} |\displaystyle \frac {\frac{9}{15} - \frac{7}{15}} {\frac{9}{15} + \frac{7}{15}}| + \cdots + 0.3(\frac{3}{15} + \frac{0}{15})^{0.5} |\displaystyle \frac {\frac{3}{15} - \frac{0}{15}} {\frac{3}{15} + \frac{0}{15}}| }{ 0.2(\frac{9}{15} + \frac{7}{15})^{0.5} + \cdots + 0.3(\frac{3}{15} + \frac{0}{15})^{0.5} }
\ \

G &\approx \displaystyle \frac {0.03 + 0.23 + 0.21 + 0.26 + 0.49 + 0.08 + 0.27 + 0.13} {0.21 + 0.23 + 0.21 + 0.26 + 0.77 + 0.58 + 0.60+ 0.13} \ \

G &= \displaystyle \frac{1.701419}{2.986235} \ \ G &= 0.569754

\end{align*}

Variance Adjusted

• Chang et al, 2011: Variance Adjusted Weighted (VAW) UniFrac
• R Package abdiv: variance_adjusted_unifrac()
• qiime2 qiime diversity beta-phylogenetic --p-metric weighted_normalized_unifrac --p-variance-adjusted

\begin{align*}

V &= \displaystyle \frac {\sum_{i = 1}^{n} L_i\displaystyle \frac {|\frac{A_i}{A_T} - \frac{B_i}{B_T}|} {(A_T + B_T)(A_T + B_T - A_i - B_i)} } {\sum_{i = 1}^{n} L_i\displaystyle \frac {\frac{A_i}{A_T} + \frac{B_i}{B_T}} {(A_T + B_T)(A_T + B_T - A_i - B_i)} } \ \

\end{align*}

cmmr/rbiom documentation built on June 10, 2025, 9:24 p.m.