526互联

R : 生成一个堆叠图用于展示OTU在不同分类水平上的相对丰度，并结合一个聚类树进行可视化

发布时间 2023-05-24 21:11:01作者: 王哲MGG_AI

setwd("E:\\中国农业科学院\\20220927宏基因组教学\\02后期分析\\01堆叠图")
rm(list = ls())
library(tidyverse)
library(ggplot2)
library(ggtree)
library(treeio)
library(ggsci)
library(cowplot)
otu = read.table('top10.2.txt', row.names = 1, sep = '\t', head = TRUE, check.names = FALSE)
tree = hclust(vegan::vegdist(t(otu), method = 'bray'),
method = 'average') %>%
as.phylo()

tree = groupClade(tree, .node=c(41))
p1 = ggtree(tree,branch.length="none") + geom_tiplab(color="black",size=3) + xlim(NA, 70)
p1

otu$mean = apply(otu, 1, mean)
otu = otu[order(otu$mean),]
otu$phylum = factor(rownames(otu), levels = rownames(otu))
otu = otu[, -c(ncol(otu)-1)]

p2 = otu %>%
reshape2::melt(id.vars = 'phylum') %>%
ggplot(aes(variable, value, fill = phylum))+
geom_bar(stat = 'identity', position = 'fill', colour= "black",size=.3)+
scale_x_discrete(limits = c("RE185","SE185","SE339","SE54","SE455","GE331","GE533","GE455","GE165","GE185","GE200","SE200","GE339","GE54","RE455","RE533","RE54","SE533","SE165","SE331","BS200","RS200","RE331","RS339","RS54","BS165","RS185","RE200","BS54","RE165","BS533","RE339","BS331","BS185","RS331","RS455","BS339","BS455","RS165","RS533"))+
scale_fill_igv()+
scale_fill_manual(values = c("#A6CEE3","#FFFF99","#B2DF8A","#33A02C","#B15928","#E31A1C","#FDBF6F","#FF7F00","#6A3D9A","#1F78B4","#FB9A99")) +
scale_y_continuous(expand = c(0,0))+
scale_y_continuous(labels = scales::percent)+
coord_flip()+
theme_classic()+
theme(axis.ticks.y = element_blank(),
axis.title.y = element_blank(),
axis.text.y = element_blank(),
axis.line = element_blank())+
labs(y = 'Percentage')

ggdraw()+
draw_plot(p1, 0, 0.055, 0.8, 0.945)+
draw_plot(p2, 0.25, 0, 0.7, 1)

#################################################################################################

这段R语言代码的功能是生成一个堆叠图（stacked plot）用于展示OTU（Operational Taxonomic Units）在不同分类水平上的相对丰度，并结合一个聚类树（cluster tree）进行可视化。

下面是代码的解释：

  1. 设置工作目录和清空变量：

    • setwd("E:\\中国农业科学院\\20220927宏基因组教学\\02后期分析\\01堆叠图")：设置工作目录为指定路径。
    • rm(list = ls())：清空当前环境中的所有变量。

  2. 载入所需的R包：

    • library(tidyverse)：加载Tidyverse包，包含了多个常用的数据处理和可视化包。
    • library(ggplot2)：加载ggplot2包，用于创建图形。
    • library(ggtree)：加载ggtree包，用于创建聚类树的可视化。
    • library(treeio)：加载treeio包，用于读取和处理聚类树数据。
    • library(ggsci)：加载ggsci包，用于设置颜色主题。
    • library(cowplot)：加载cowplot包，用于组合多个图形。

  3. 读取数据：

    • otu = read.table('top10.2.txt', row.names = 1, sep = '\t', head = TRUE, check.names = FALSE)：从名为top10.2.txt的文件中读取OTU数据，使用第一列作为行名，以制表符分隔，文件包含列名。

  4. 构建聚类树：

    • tree = hclust(vegan::vegdist(t(otu), method = 'bray'), method = 'average') %>% as.phylo()：使用OTU数据计算Bray-Curtis距离并进行聚类分析，得到聚类树对象。

  5. 对聚类树进行分组：

    • tree = groupClade(tree, .node=c(41))：将聚类树的节点41进行分组。

  6. 创建第一个图形p1（聚类树）：

    • p1 = ggtree(tree,branch.length="none") + geom_tiplab(color="black",size=3) + xlim(NA, 70)：使用ggtree包绘制聚类树，禁用分支长度显示，设置叶节点标签颜色为黑色、大小为3，设置x轴范围为自动缩放。

  7. 对OTU数据进行处理：

    • otu$mean = apply(otu, 1, mean)：计算OTU各行的平均值，并将结果保存在新列"mean"中。
    • otu = otu[order(otu$mean),]：按照平均值对OTU数据进行升序排序。
    • otu$phylum = factor(rownames(otu), levels = rownames(otu))：将OTU数据的行名设置为因子变量"phylum"，以保持排序顺序。
    • otu = otu[, -c(ncol(otu)-1)]：删除倒数第二列（平均值所在列）。

  8. 创建第二个图形p2（堆叠图）：

    • p2 = otu %>% reshape2::melt(id.vars = 'phylum') %>% ggplot(aes(variable, value, fill = phylum)) + ...：使用reshape2包将OTU数据进行重塑，以便于堆叠图的绘制。设置x轴为变量，y轴为值，颜色填充为"phylum"因子。接下来是一系列设置图形样式和主题的代码。

  9. 组合图形并绘制：

    • ggdraw() + draw_plot(p1, 0, 0.055, 0.8, 0.945) + draw_plot(p2, 0.25, 0, 0.7, 1)：使用cowplot包中的函数将图形p1和p2组合在一起，并设置它们的位置和大小。

最后，代码生成了一个堆叠图和聚类树的组合图形，用于展示OTU的相对丰度和它们在分类水平上的关系。

################################################################################################

逐行解释您提供的R代码：

  1. 设置工作目录：

    • setwd("D:\\Desktop")：将工作目录设置为 "D:\Desktop"。这将是R语言执行文件读取和写入操作的默认位置。

  2. 清空变量：

    • rm(list = ls())：清空当前R环境中的所有对象，以确保开始时没有任何变量存在。

  3. 载入包：

    • library(tidyverse)：加载tidyverse包，该包提供了一套用于数据处理和可视化的R包。
    • library(ggplot2)：加载ggplot2包，该包用于绘制数据图形。
    • library(ggtree)：加载ggtree包，该包用于绘制漂亮的进化树。
    • library(treeio)：加载treeio包，该包用于处理和操作进化树。
    • library(ggsci)：加载ggsci包，该包提供了一套科学颜色主题。
    • library(cowplot)：加载cowplot包，该包用于组合多个图形。

  4. 读取数据：

    • otu = read.table('top10.2.txt', row.names = 1, sep = '\t', head = TRUE, check.names = FALSE)：从'top10.2.txt'文件中读取数据，并将第一列作为行名。文件使用制表符作为列分隔符。

  5. 构建聚类树：

    • tree = hclust(vegan::vegdist(t(otu), method = 'bray'), method = 'average') %>% as.phylo()：根据OTU数据计算样本间的Bray-Curtis距离，然后使用聚类方法构建聚类树。最后，将结果转换为phylo对象。

  6. 对聚类树进行分组：

    • tree = groupClade(tree, .node=c(41))：在聚类树中根据给定的节点（41号节点）进行分组。

  7. 绘制第一个图形p1（聚类树）：

    • p1 = ggtree(tree, branch.length = "none") + geom_tiplab(color = "black", size = 3) + xlim(NA, 70)：使用ggtree包绘制聚类树。设置分支长度为"none"，给叶节点添加标签，并设置标签颜色为黑色、大小为3。同时，限制x轴范围在0到70之间。

  8. 数据处理和绘制第二个图形p2（堆叠图）：

    • otu$mean = apply(otu, 1, mean)：计算每行数据的平均值，并将结果存储在新的列"mean"中。
    • otu = otu[order(otu$mean),]：按照平均值的升序对数据进行排序。
    • otu$phylum = factor(rownames(otu), levels = rownames(otu))：将行名作为新的列"phylum"，并将其设置为有序因子，按照原始数据的行顺序设置水平顺序。
    • otu = otu[, -c(ncol(otu)-1)]：移除倒数第二列。

  9. 使用ggplot2绘制堆叠图：

    • p2 = otu %>% reshape2::melt(id.vars = 'phylum') %>% ggplot(aes(variable, value, fill = phylum)) + geom_bar(stat = 'identity', position = 'fill', colour= "black", size = .3)：使用reshape2包的melt函数将OTU数据转换为长格式，并利用ggplot2包绘制堆叠柱状图。设置x轴为变量名，y轴为数值，堆叠填充颜色按照"phylum"列进行分组。
    • scale_x_discrete(limits = c(...))：设置x轴的离散值限制，按照指定顺序显示。
    • scale_fill_igv()：使用ggsci包中的配色方案设置填充颜色。
    • scale_fill_manual(values = c(...))：手动设置填充颜色的值。
    • scale_y_continuous(expand = c(0,0))：设置y轴的连续值范围，去除空白。
    • scale_y_continuous(labels = scales::percent)：将y轴刻度标签转换为百分比格式。
    • coord_flip()：交换x轴和y轴的位置，使图形为横向堆叠柱状图。
    • theme_classic()：使用经典主题样式。
    • theme(...)：对图形的轴线、刻度、文本等进行自定义的主题设置。
    • labs(y = 'Percentage')：设置y轴标签为"Percentage"。

  10. 组合图形并输出：

    • ggdraw() + draw_plot(p1, 0, 0.055, 0.8, 0.945) + draw_plot(p2, 0.25, 0, 0.7, 1)：使用cowplot包的ggdraw函数创建画布，并利用draw_plot函数将p1和p2两个图形绘制在画布上的指定位置。最后，输出组合图形。

这段代码的目的是读取OTU数据，构建聚类树并绘制聚类树图形，然后对OTU数据进行处理并绘制堆叠图，最后将聚类树图形和堆叠图形组合在一起输出。