In grunwaldlab/metacoder: Tools for Parsing, Manipulating, and Graphing Taxonomic Abundance Data

old way

source('http://bioconductor.org/biocLite.R')
biocLite('phyloseq')

phyloseq_to_taxmap <- function(phyloseq_obj) {
  tax_data <- data.frame(
    taxonomy = vapply(seq_len(nrow(phyloseq_obj@tax_table)), FUN.VALUE = character(1),
                      function(i) {
                        taxon_names <- as.vector(phyloseq_obj@tax_table[i, ])
                        rank_names <- colnames(phyloseq_obj@tax_table)
                        rank_names <- rank_names[!is.na(taxon_names)]
                        taxon_names <- taxon_names[!is.na(taxon_names)]
                        paste(rank_names, taxon_names, sep = "::", collapse = ";;")
                      }),
    phyloseq_id = rownames(phyloseq_obj@tax_table)
  )

  otu_data  <- as.data.frame(t(phyloseq_obj@otu_table))
  otu_data$phyloseq_id <- rownames(otu_data)
  tax_data <- merge(tax_data, otu_data, by.x = "phyloseq_id", by.y = "phyloseq_id")

  parse_taxonomy_table(tax_data, taxon_col = c("class" = 2), other_col_type = "obs_info",
                       class_regex = "^(.*)::(.*)$",
                       class_key = c(rank = "taxon_info", "name"),
                       class_sep = ";;")
}

library(phyloseq)
phy <- readRDS("example.phyloseq")

data <- phyloseq_to_taxmap(phy)

New way with taxa taxmap

# devtools::install_github("grunwaldlab/metacoder")
# devtools::install_github("ropensci/taxa")

library(metacoder)
library(taxa)

phy <- readRDS("example.phyloseq")
obj <- parse_phyloseq(phy)

obj$filter_taxa(taxon_names == "Bacteria", subtaxa = TRUE)

healthy_samples <- obj$data$sam_data$sample_id[obj$data$sam_data$status == "healthy"]
healthy_counts <- obs_apply(obj, "otu_table", function(i) sum(obj$data$otu_table[i, healthy_samples]), simplify = TRUE)

obj %>%
  heat_tree(node_size = healthy_counts, 
            node_color = healthy_counts,
            node_label = taxon_names)

Even newer way

devtools::install_github("grunwaldlab/metacoder")
devtools::install_github("ropensci/taxa")

library(metacoder)

phy <- readRDS("example.phyloseq")
obj <- parse_phyloseq(phy)

# Convert counts to proportions
obj$data$otu_table <- calc_obs_props(obj,
                                     dataset = "otu_table",
                                     cols = obj$data$sam_data$sample_ids,
                                     other_cols = TRUE)

# Calculate per-taxon proportions 
obj$data$tax_table <- calc_taxon_abund(obj, 
                                       dataset = "otu_table", 
                                       cols = obj$data$sam_data$sample_ids)

# Calculate difference between treatments
obj$data$diff_table <- compare_treatments(obj,
                                          dataset = "tax_table",
                                          sample_ids = obj$data$sam_data$sample_ids,
                                          treatments = obj$data$sam_data$status,
                                          other_cols = TRUE)

# Plot differneces
color_interval <- c(-5, 5) # The range of values (log 2 ratio of median proportion) to display
obj %>%
  taxa::filter_taxa(taxon_names == "Bacteria", subtaxa = TRUE) %>%
  heat_tree(node_size_axis_label = "Number of OTUs",
            node_size = n_obs,
            node_color_axis_label = "Log 2 ratio of median proportions",
            node_color = log2_median_ratio,
            node_color_range = diverging_palette(),
            node_color_trans = "linear",
            node_color_interval = color_interval,
            edge_color_interval = color_interval,
            node_label = taxon_names,
            node_label_max = 150)

grunwaldlab/metacoder documentation built on Feb. 22, 2024, 3:47 a.m.