1.R语言物种组成图，物种堆叠图

      我们测16s多样性等会得到细菌群落的OTU表和注释表，这时候难免要比较不同样本的群落结构的差异。这个时候一般看门或者属，种水平一般是达不到的。于是我们需要画物种组成图类似于下面这张（网上下载的，图1）。那我们需要一个OTU表和注释表合并后的表格，我这个实例里叫做“筛选otu注释1.xlsx”（图2），然后根据这个表格利用我们的代码就能顺利出图啦。

              图1                             

                                                                  图2

下面就是R代码，输入文件：筛选otu注释1.xlsx，A1-C9是样本名，需要你自己修改，假如你的样本第一个是Q1,最后一个是M1，代码里的A1:C9也应该全部改成Q1：M1。改成自己的文件时，确保格式和我的例子一样，第一列是OTU名，第二列是注释，后面的列是样本，第一，二列列名也应该一样。

#install.packages("readxl")
library(readxl)
# 替换以下路径为您的Excel文件路径
data <- read_excel("筛选otu注释1.xlsx")

# 假设第二列的列名为 "Column2"
# 下面的代码将 "Column2" 分割成7个新的列
data_separated <- data %>%
  separate(Taxon, into = c("界", "门", "纲", "目", "科", "属", "种"), sep = ";", remove = TRUE, extra = "merge")

library(dplyr)

# 计算每列的总和
column_sums <- data_separated %>%
  summarise(across(A1:C9, sum, na.rm = TRUE))

# 转换为丰度值
data_separated <- data_separated %>%
  mutate(across(A1:C9, ~ .x / column_sums[[cur_column()]]))


data_aggregated <- data_separated %>%
  group_by(门) %>%
  summarise(across(A1:C9, sum, na.rm = TRUE))

library(ggplot2)

# 1. 计算每个门的总丰度
total_abundance_by_Phylum <- data_aggregated %>%
  rowwise() %>%
  mutate(TotalAbundance = sum(c_across(A1:C9), na.rm = TRUE)) %>%
  ungroup()

# 2. 选择除了NA之外的丰度最高的前9个门
top_10_phyla <- total_abundance_by_Phylum %>%
  filter(!is.na(门)) %>%
  arrange(desc(TotalAbundance)) %>%
  slice_head(n = 9) %>%
  pull(门)

# 3. 标记其余的门和NA为 "Others"
data_aggregated <- data_aggregated %>%
  mutate(门 = ifelse(门 %in% top_10_phyla, 门 , "Others"))


data_aggregated <- data_aggregated %>%
  group_by(门) %>%
  summarise(across(A1:C9, sum, na.rm = TRUE))



# 5. 转换数据为长格式并创建堆叠图
data_long <- data_aggregated %>%
  pivot_longer(-门, names_to = "Sample", values_to = "Abundance")

library(viridis)  # 加载viridis包
library(scales) 
# 计算需要的颜色数量
num_colors <- length(unique(data_long$门))

# 生成HSV颜色
hsv_colors <- hsv(h = seq(0, 1, length = num_colors), s = 0.8, v = 0.8)

# 手动指定样本的顺序 调整X轴排序
Phylumed_samples <- c(paste0("A", 1:12), paste0("B", 1:12), paste0("C", 1:12))

# 将样本名转换为因子，并指定水平为我们手动指定的顺序
data_long$Sample <- factor(data_long$Sample, levels = Phylumed_samples)




# 获取门的唯一值并排序
phylum_order <- unique(data_aggregated$门)

# 绘制堆叠图
ggplot(data_long, aes(x = Sample, y = Abundance, fill = 门)) +
  geom_bar(stat = "identity", position = "fill") +
  theme(panel.background = element_blank(),
        axis.line = element_line(),
        axis.text.x = element_text(angle = 45, hjust = 1)) +
  scale_fill_manual(values = hsv_colors, breaks = phylum_order) +  # 设置颜色和图例顺序
  scale_y_continuous(expand = expansion(mult=c(0.01,0.01)),labels = scales::label_percent(suffix = "")) +
  labs(x = "Sample", y = "Abundance (%)", fill = "Phylum")

最后我们就可以得到物种组成图啦

关注我，以后的内容更精彩哦。
b站号：羽球最强生信人

闲鱼号：小秋家的小卖铺
淘宝号：小秋家的小卖铺