使用ggplot2绘制菌群组成堆叠图
2020-07-14 本文已影响0人
你猜我菜不菜
堆叠图主要看不同菌群在该样本中的相对比例，把每个样本的情况或者每个处理组的情况展示出来，可以大致看出样本的重复性好坏和不同处理间的某些菌群变化趋势。
1.门水平的堆叠图
###设置工作目录
getwd()
path1 = ("C:/Users/Administrator/Desktop/R_work/03_irradiation_gutMicro/00_visualization")
setwd(path1)
list.files()
############门水平组成堆叠图#####################################
####读取门水平数据
phylum <- read.delim('phylum_table.txt', row.names = 1, sep = '\t', stringsAsFactors = FALSE, check.names = FALSE)
#将门水平进行排序
phylum$sum <- rowSums(phylum)
phylum <- phylum[order(phylum$sum, decreasing = TRUE), ]
phylum
#选取相对丰度
phylum_top10 <- phylum[1:3, -ncol(phylum)]
phylum_top10['Others', ] <- 100 - colSums(phylum_top10)
phylum_top10
####ggplot2堆叠图
library(reshape2)
library(ggplot2)
#调整数据布局
phylum_top10$Taxonomy <- factor(rownames(phylum_top10), levels = rev(rownames(phylum_top10)))
phylum_top10 <- melt(phylum_top10, id = 'Taxonomy')
#添加分组信息
group <- read.delim('group.txt', sep = '\t', stringsAsFactors = FALSE)
names(group)[1] <- 'variable'
phylum_top10 <- merge(phylum_top10, group, by = 'variable')
#因子排序，对样品排序，对分组排序
phylum_top10$variable <- factor(phylum_top10$variable,
levels = c("UN1_1", "UN1_2", "UN1_3", "UN1_4",
"UN1_5", "UN1_6", "IR1_1", "IR1_2",
"IR1_3", "IR1_4", "IR1_5", "IR1_6",
"UN7_1", "UN7_2", "UN7_3", "UN7_4",
"UN7_5", "UN7_6", "IR7_1", "IR7_2",
"IR7_3", "IR7_4", "IR7_5", "IR7_6",
"UN14_1", "UN14_2", "UN14_3", "UN14_4",
"UN14_5", "UN14_6", "IR14_1", "IR14_2",
"IR14_3", "IR14_4", "IR14_5", "IR14_6"), ordered = TRUE)
phylum_top10$group <- factor(phylum_top10$group,
levels = c("Unirradiated","Irradiated"),ordered = TRUE)
phylum_top10
#作图
#颜色不够用，设置更多梯度变化的颜色
#colourCount = length(unique(phylum_top10$Taxonomy))
#getPalette = colorRampPalette(brewer.pal(9, "Set1"))
p <- ggplot(phylum_top10, aes(variable, value, fill = Taxonomy)) +
geom_col(position = "stack", width = 0.8) +
facet_wrap(~group, scales = 'free_x', ncol = 2) +
scale_fill_manual(values = rev(c( "#E64B35B2", "#4DBBD5B2",
"#3C5488B2","yellow"))) +
labs(x = '', y = 'Relative Abundance(%)', fill = "Phylum") +
theme(panel.grid = element_blank(),
panel.background = element_rect(color = '#636363',
fill = 'transparent'),
strip.text = element_text(size = 20),
strip.text.x = element_text(size = 20, colour = "black",
face = "bold"),
strip.background = element_rect(color="grey95",
fill="grey95",
size=1.5,
linetype="solid"),
legend.title = element_text(size = 15,face = "bold",
vjust = 1, hjust = 0),
legend.text = element_text(size = 12),
legend.key.size=unit(.5,'cm'),
axis.title.x = element_text(size = 16, vjust = 0.5,
hjust = 0.5),
axis.title.y = element_text(size = 15, vjust = 0.5,
hjust = 0.5),
axis.text.x = element_text(angle = 90, size = 10,
vjust = 0.5, hjust = 0.5),
axis.text.y = element_text(size = 15,vjust = 0.5,
hjust = 0.5),
legend.position = "bottom") +
guides(fill=guide_legend(nrow=1,reverse=TRUE))
p
2.科水平的堆叠图
和门水平的代码基本一样，这样重复的代码越来越多，经常是在满眼的代码中寻找参数去修改调整，眼要瞎，逼迫我使用function()和一些循环函数如for，lapply等，这次还是重复代码跑一遍，下次一定要简化工作！
############科水平组成堆叠图#####################################
####读取科水平数据
family <- read.delim('family_table.txt', row.names = 1, sep = '\t', stringsAsFactors = FALSE, check.names = FALSE)
#将门水平进行排序
family$sum <- rowSums(family)
family <- family[order(family$sum, decreasing = TRUE), ]
family
#选取相对丰度
family_top20 <- family[1:21, -ncol(family)]
family_top20
#排序后数据导出成csv文件
write.csv(family_top20, 'family_top20.csv', quote = FALSE)
####ggplot2堆叠图
library(reshape2)
library(ggplot2)
#调整数据布局
family_top20$Taxonomy <- factor(rownames(family_top20),
levels = rev(rownames(family_top20)))
family_top20 <- melt(family_top20, id = 'Taxonomy')
#添加分组信息
group <- read.delim('group.txt', sep = '\t', stringsAsFactors = FALSE)
names(group)[1] <- 'variable'
family_top20 <- merge(family_top20, group, by = 'variable')
#因子排序，对样品排序，对分组排序
family_top20$variable <- factor(family_top20$variable,
levels = c("UN1_1", "UN1_2", "UN1_3", "UN1_4",
"UN1_5", "UN1_6", "IR1_1", "IR1_2",
"IR1_3", "IR1_4", "IR1_5", "IR1_6",
"UN7_1", "UN7_2", "UN7_3", "UN7_4",
"UN7_5", "UN7_6", "IR7_1", "IR7_2",
"IR7_3", "IR7_4", "IR7_5", "IR7_6",
"UN14_1", "UN14_2", "UN14_3", "UN14_4",
"UN14_5", "UN14_6", "IR14_1", "IR14_2",
"IR14_3", "IR14_4", "IR14_5", "IR14_6"), ordered = TRUE)
family_top20$group <- factor(family_top20$group,
levels = c("Unirradiated","Irradiated"),ordered = TRUE)
family_top20
#作图
#颜色不够用，设置更多梯度变化的颜色
#colourCount = length(unique(family_top20$Taxonomy))
#getPalette = colorRampPalette(brewer.pal(9, "Set1"))
p1 <- ggplot(family_top20, aes(variable, value, fill = Taxonomy)) +
geom_col(position = 'stack', width = 0.8) +
facet_wrap(~group, scales = 'free_x', ncol = 2) +
scale_fill_manual(values = rev(c( "#3C5488B2","#00A087B2",
"#F39B7FB2","#91D1C2B2",
"#8491B4B2", "#DC0000B2",
"#7E6148B2","yellow",
"darkolivegreen1", "lightskyblue",
"darkgreen", "deeppink", "khaki2",
"firebrick", "brown1", "darkorange1",
"cyan1", "royalblue4", "darksalmon",
"darkgoldenrod1", "darkseagreen", "darkorchid"))) +
labs(x = '', y = 'Relative Abundance(%)', fill = "Family") +
theme(panel.grid = element_blank(),
panel.background = element_rect(color = '#636363',
fill = 'transparent'),
strip.text = element_text(size = 20),
strip.text.x = element_text(size = 20, colour = "black", face = "bold"),
strip.background = element_rect(color="grey95",
fill="grey95",
size=1.5,
linetype="solid"),
legend.title = element_text(size = 15,face = "bold",
vjust = 1, hjust = 0),
legend.text = element_text(size = 11),
legend.key.size=unit(.4,'cm'),
axis.title.x = element_text(size = 16, vjust = 0.5,
hjust = 0.5),
axis.title.y = element_text(size = 15, vjust = 0.5,
hjust = 0.5),
axis.text.x = element_text(angle = 90, size = 10,
vjust = 0.5, hjust = 0.5),
axis.text.y = element_text(size = 15,vjust = 0.5,
hjust = 0.5),
legend.position="bottom") +
guides(fill=guide_legend(nrow=6,reverse=TRUE))
p1
3.拼图
#拼图
library(patchwork)
p/p1
堆叠图的配色比较重要，颜色不好看起来很乱很躁，谷歌出来的一套配色觉得还行。
