R语言用ggplot2绘制堆积柱状图遇到的问题(y轴顺序、facet_wrap分面顺序)- 需要将名称因子化
1、用ggplot2绘制菌的堆积柱状图时,y轴的顺序是按照菌属的首字母排序的,而不是按照我表格中顺序,代码如下(以门水平的数据示例):
这里,我们先将原始数据处理,筛选出排名前4的门,其它的用Others代替。
#import data and modify it
phylum <- read.table('phylumAbundance.txt', sep="\t", header=T, row.names=1)
phylum.ave <- apply(phylum, 1, FUN=mean)
phylum.2 <- cbind(phylum, phylum.ave)[order(-phylum.ave),]
#sum(phylum.2$phylum.ave[1:4]) #check the summary of top 4 taxonomies
phylum.2 <- subset(phylum.2[1:4,], select=-phylum.ave)
phylum.2 <- rbind(phylum.2, others=apply(phylum.2, 2, function(x){100-sum(x)}))
#apply(phylum.2, 2, function(x){sum(x)}) #check the summary of each column
phylum.2 <- cbind(PhylumID=row.names(phylum.2), phylum.2)
#convert the data format
library(reshape2)
phylum.gg <- melt(phylum.2, id.vars="PhylumID", variable.name="SampleID", value.name="Abundance")
#sum(phypum.gg$Abundance) #check
#plot bar chart
library(ggplot2)
ggplot(phylum.gg, aes(SampleID, Abundance, fill=PhylumID)) +
geom_bar(stat="identity") +
guides(fill=guide_legend(reverse=F)) +
scale_y_continuous(expand=c(0,0))
得到的图如下
经过查找资料后发现,需要将菌的名称因子化(factor),于是在将数据框由“宽”变“长”前,先将菌的名称因子化:(tips:不可在转为长数据后再因子化,这样会报错:显示名称有重复。)
phylum.2$PhylumID <- factor(phylum.2$PhylumID, levels = rev(phylum.2$PhylumID))
重新绘图得到正确的图:
2、为了让图片更美观,将上述数据分成三组“M,L,C”,用“facet_wrap”绘制分面图。一开始,我没有将组别名称因子化,得到的分面结果是按照“C,L,M”排序的。
在绘图前加入分组信息:
phylum.gg$group <- c(rep('M', 30), rep('L', 40), rep('C', 40))
ggplot(phylum.gg, aes(SampleID, Abundance, fill = PhylumID)) +
geom_col(position = 'stack', width = 0.6) +
facet_wrap(~group, scales = 'free_x', ncol = 3) +
scale_fill_manual(values = rev(c('blue', 'orange', 'green', 'yellow', 'gray'))) +
labs(x = '', y = 'Relative Abundance(%)') +
theme(panel.grid = element_blank(), panel.background = element_rect(color = 'black', fill = 'transparent'), strip.text = element_text(size = 12)) +
theme(axis.text = element_text(size = 11), axis.title = element_text(size = 13), legend.text = element_text(size = 11, face='italic'))
最开始我没搞明白为什么会这样,但我发现了是按字母排序的。
方法1(可行):于是我将组别名称改为“A,B,C”,再修改组别名称,从而得到想要的结果。
参考方法:https://ggplot2.tidyverse.org/reference/as_labeller.html
方法2(不可行):后面我知道需要将组别名称因子化才行,于是想将组别名称先做到一张表里,因子化后再与数据合并。
groupData <- read.delim('group.txt', sep = '\t', stringsAsFactors = FALSE)
groupData$group <- factor(groupData$group, levels = rev(groupData$group))
这时会报错:
Error in `levels<-`(`*tmp*`, value = as.character(levels)) :
factor level [2] is duplicated
和处理y轴菌的名称时一样,会提示有重复。
这样做是参考了http://www.360doc.com/content/17/0912/21/46931810_686600613.shtml,但这个里面和我的数据不太一样。
其它:在查找过程中看到的别人绘制堆积柱状图的代码,数据处理过程要比我的简单,供参考。
来源:http://blog.sciencenet.cn/blog-3406804-1166733.html
####读取并挑选 top4丰度门
#读取数据
phylum <- read.table('phylumAbundance.txt', row.names = 1, sep = '\t', header=T)
#求各类群的丰度总和,并排序
phylum$sum <- rowSums(phylum)
phylum <- phylum[order(phylum$sum, decreasing = TRUE), ]
#挑选 top4门类群,并将 top4外的类群合并为“Others”
phylum_top4 <- phylum[1:4, -ncol(phylum)]
phylum_top4['Others', ] <- 100 - colSums(phylum_top4)
#可选输出(如 csv 格式)
#write.table(phylum_top4, 'phylum_top4.txt', sep="\t", quote = FALSE)
####ggplot2 作图(一个简单示例)
library(reshape2) #用于排列数据
library(ggplot2) #ggplot2 作图
#调整数据布局
phylum_top4$Phylum <- factor(rownames(phylum_top4), levels = rev(rownames(phylum_top4)))
phylum_top4 <- melt(phylum_top4, id = 'Phylum')
#添加分组
groupData <- read.delim('group.txt', sep = '\t', stringsAsFactors = FALSE)
names(groupData)[1] <- 'variable'
phylum_top4 <- merge(phylum_top4, groupData, by = 'variable', sort=F)
#绘制带分面的柱状图
real_group <- as_labeller(c(`A`="Model", `B`="LHE", `Control`="Control"))
p <- ggplot(phylum_top4, aes(variable, value, fill = Phylum)) +
geom_col(position = 'stack', width = 0.6) +
facet_wrap(~group, scales = 'free_x', ncol = 3, labeller=real_group) +
scale_fill_manual(values = rev(c('blue', 'orange', 'green', 'yellow', 'gray'))) +
labs(x = '', y = 'Relative Abundance(%)') +
theme(panel.grid = element_blank(), panel.background = element_rect(color = 'black', fill = 'transparent'), strip.text = element_text(size = 12)) +
theme(axis.text = element_text(size = 11), axis.title = element_text(size = 13), legend.text = element_text(size = 11, face='italic'))
p
ggsave('ggplot2_plot.png', p, width = 8, height = 6)