From 47fe4de8891277aaa2ffafec00631c66602bb58f Mon Sep 17 00:00:00 2001 From: Louis Date: Mon, 13 Oct 2025 10:26:20 +0200 Subject: [PATCH] Created a script to generate nice visuals for slides --- generate-figures.R | 82 ++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 82 insertions(+) create mode 100644 generate-figures.R diff --git a/generate-figures.R b/generate-figures.R new file mode 100644 index 0000000..7c6d2e2 --- /dev/null +++ b/generate-figures.R @@ -0,0 +1,82 @@ +source("utils.R") +library(stringr) + +if (!file.exists("data/hmc-phyloseq-europe.Rds")) { + source("load-europe.R") +} else { + cpd_phyloseq_europe <- readRDS("data/hmc-phyloseq-europe.Rds") +} + +df <- otu_table(cpd_phyloseq_europe) +colnames(df) <- paste0("$S_{", seq(1, ncol(df)), "}$") + + +nr <- 6 +nc <- 6 + +sub_df <- ellipse_table(df, nr = nr, nc = nc) + +taxo_df <- + rownames(sub_df)[-(ceiling(nr / 2) + 1)] |> + str_replace_all(" ", "") |> + strsplit("\\.") |> + do.call(what = "rbind") |> + as.data.frame() + +colnames(taxo_df) <- c("Kingdom", "Phylum", "Class", "Order", "Family", "Genus") + +rownames(sub_df) <- c(taxo_df[seq(1, ceiling((nr / 2))), "Genus"], rownames(sub_df)[(ceiling(nr / 2) + 1)], taxo_df[seq(ceiling((nr / 2)) + 1, nrow(taxo_df)), "Genus"]) + + +kable(sub_df, format = "latex", escape = FALSE, align = "c") |> + writeLines(con = file.path("outfiles", "otu-tabular.tex")) + +library(data.tree) +library(ape) +library(ggtree) +library(ggplot2) +library(dplyr) +# Create a path string column from your hierarchy +levels_vec <- c("Kingdom", "Phylum", "Class", "Order", "Family", "Genus") +taxo_df$pathString <- apply( + taxo_df[, levels_vec], + 1, + function(x) paste(x, collapse = "/") +) + +# Now build the tree +taxo_tree <- as.Node(taxo_df) + +taxo_phylo <- as.phylo(taxo_tree) + +# 3. Extract node data +p <- ggtree(taxo_phylo, layout = "rectangular") +node_data <- p$data + +# 4. Build lookup table: each label → its rank +lookup <- taxo_df %>% + tidyr::pivot_longer(cols = all_of(levels_vec), names_to = "rank_name", values_to = "label") + +# 5. Join node data with taxonomy ranks +node_data2 <- node_data %>% + left_join(lookup, by = "label") %>% + mutate(rank_name = factor(rank_name, levels = levels_vec)) # enforce order + +internal_nodes <- node_data2 + +library(ggokabeito) +# 6. Plot with colored nodes +(phylo_plot <- ggtree(taxo_phylo, layout = "slanted") %<+% node_data2 + + geom_point(aes(color = rank_name), size = 3, na.rm = TRUE) + + geom_text2(aes(label = label, color = rank_name), + hjust = 1.1, vjust = -1, size = 4, na.rm = TRUE + ) + + scale_colour_okabe_ito(name = "Taxonomic Rank", order = c(1:3, 5:7)) + + theme_void() + + theme(legend.position = "none")) +library(tikzDevice) +# Set up tikzDevice to use standalone document class +options(tikzDocumentDeclaration = "\\documentclass[10pt]{standalone}") +tikz(file = file.path("outfiles", "otu-phylo.tex"), height = 9, , width = 15, standAlone = TRUE) +phylo_plot +dev.off()