R绘图模板——进化树+柱状堆积图+热图！

 

点击上方

“科研后花园”

关注我们

代码如下：

1、设置工作环境、加载R包

rm(list=ls())#clear Global Environment
setwd('D:\\test')#设置工作路径

#加载R包
library (reshape2)
library(ggplot2) # Create Elegant Data Visualisations Using the Grammar of Graphics
library(plyr) # Tools for Splitting, Applying and Combining Data
library(ggtree) # an R package for visualization of tree and annotation data
library(vegan) # Community Ecology Package
library(ape) # Analyses of Phylogenetics and Evolution
library( RColorBrewer) # not installed on this machine
library(aplot) # Decorate a 'ggplot' with Associated Information

2、加载数据并进行数据处理

#读取数据
df1 <- read.table(file="Genus.txt",sep="\t",header=T,check.names=FALSE)
#查看前6行
head(df1)
##利用循环处理具有重复的数据
#初始化
data<-aggregate(E ~ Tax,data=df1,sum)
#重命名
colnames(data)[2]<-"example"
for (i in colnames(df1)[2:length(colnames(df1))]){
  #计算每列的和
  data1<-aggregate(df1[,i]~Tax,data=df1,sum)
  colnames(data1)[2]<-i  
  #合并列
  data<-merge(data,data1,by="Tax")
}
df2<-data[,-2]
rownames(df2)=df2$Tax#修改行名
df3=df2[,-1]#删除多的列

#计算物种总丰度并降序排列
df3$rowsum <- apply(df3,1,sum)
df4 <- df3[order (df3$rowsum,decreasing=TRUE),]
df5 = df4[,-6]#删除求和列
#求物种相对丰度
df6 <- apply(df5,2,function(x) x/sum(x))
#由于之间已经按照每行的和进行过升序排列，所以可以直接取前10行
df7 <-  df6[1:10,]
df8 <- 1-apply(df7, 2, sum) #计算剩下物种的总丰度
#合并数据
df9 <- rbind(df7,df8)
row.names(df9)[11]="Others"
#变量格式转换,宽数据转化为长数据,方便后续作图
df_genus <- melt(df9)
names(df_genus)[1:2] <- c("Taxonomy","sample")  #修改列名
#颜色
col <-colorRampPalette(brewer.pal(12,"Paired"))(11)

3、柱状堆积图绘制

##柱状堆积图绘制绘图
p1<-ggplot(df_genus, aes( x = sample,y=100 * value,fill = Taxonomy))+#geom_col和geom_bar这两条命令都可以绘制堆叠柱形图
  geom_col(position = 'stack',width = 0.8)+
  #geom_bar(position = "stack", stat = "identity", width = 0.6) 
  scale_y_continuous(expand = c(0,0))+# 调整y轴属性，使柱子与X轴坐标接触
  labs(x=NULL,y="Relative Abundance(%)",#设置X轴和Y轴的名称以及添加标题
       fill="Taxonomy")+
  guides(fill=guide_legend(keywidth = 1, keyheight = 1)) +
  coord_flip()+
  theme(panel.grid = element_blank(),
        panel.background = element_blank())+
  scale_fill_manual(values = col)+
  theme(axis.text.x = element_text(color = "black",size = 10))
p1

4、进化树绘制

####聚类树
data2 <- df3[-6]
#计算距离矩阵
df_dist <- vegdist(t(data2),method = 'bray')#使用bray curtis方法计算距离矩阵
####进行层次聚类,可选择方法有single、complete、median、mcquitty、average 、centroid 、ward 
df_hc <- hclust(df_dist,method="average")#使用类平均法进行聚类
#绘图
plot(as.dendrogram(df_hc),type="rectangle",horiz=T)#实现将垂直的聚类树变成水平聚类树，绘图
# 将聚类结果转成系统发育格式
df_tree <- as.phylo(df_hc)
# 对树分组
gro <- list(group1=c("A","E"),
            group2=c("C","D"),
            group3=c("B"))
# 将分组信息和进化树组合到一起
tree<-groupOTU(df_tree,gro)
# ggtree绘图
p2 <- ggtree(tree,size=1.2)+
  geom_tiplab(aes(color=group),size=5,align = F,
              offset = 0.008,show.legend = F)+
  geom_tippoint(aes(shape=group,color=group),
              show.legend = T,
              size=4)
p2

5、热图绘制

#### 热图绘制
#去排名前20的属
df7 <-  df6[1:19,]
df8 <- 1-apply(df7, 2, sum) #计算剩下物种的总丰度
#合并数据
df9 <- rbind(df7,df8)
row.names(df9)[20]="Others"
#添加分组信息
df9 <- as.data.frame(df9)
df9$group <- factor(c(rep("a",8),rep("b",6),rep("c",4),rep("d",2)))
df9$tax <- rownames(df9)
#变量格式转换,宽数据转化为长数据,方便后续作图
df_heatmap <- melt(df9,id.vars = c("group","tax"))
#指定显示顺序
df_heatmap$tax <- factor(df_heatmap$tax,levels = df9$tax)
#颜色
col2 <-c("#1aafd0","#6a67ce","#ffb900","#fc636b")
#绘图
p2 <- ggplot()+ 
  geom_point(df_heatmap,mapping = aes(x = tax, y = variable, size = value, color=group))+
  scale_color_manual(values = col2)+
  scale_size_continuous(range = c(1, 10))+
  theme(panel.background = element_blank(),
        legend.key = element_blank(),
        axis.text.x = element_text(color = "black",angle=90,size = 10,vjust = 0.5,hjust = 0),
        axis.text.y = element_text(color = "black",size = 10),
        panel.grid.major = element_line(color = "gray"),#网格线条颜色
        panel.border = element_rect(color="black",fill=NA))+#边框色
  labs(x=NULL,y=NULL)+
  annotate("rect", xmin = 0.5, xmax = 8.5, ymin = 0, ymax = 5.8, alpha = 0.2,fill="#1aafd0") +
  annotate("rect", xmin = 8.5, xmax = 14.5, ymin = 0, ymax = 5.8, alpha = 0.2,fill="#6a67ce") +
  annotate("rect", xmin = 14.5, xmax = 18.5, ymin = 0, ymax = 5.8, alpha = 0.2,fill="#ffb900") +
  annotate("rect", xmin = 18.5, xmax = 21, ymin = 0, ymax = 5.8, alpha = 0.2,fill="#fc636b") +
  scale_x_discrete(position = "top")
p2

6、将柱状堆积图与树图和热图进行组合

#####将柱状堆积图与树图进行组合
p1 <- p1+theme(axis.text.y = element_text(color = "black",size=12))
col3 <-colorRampPalette(brewer.pal(9,"Set1"))(3)
p3 <- ggtree(tree,size=1)+
  geom_tippoint(aes(shape=group,color=group),
                show.legend = T,
                size=4)+
  scale_color_manual(values = col3)
p4 <- p2+theme(axis.text.y = element_blank())
p <- p1%>%insert_left(p3,width = 0.4)%>%
  insert_right(p4,width = 0.9)
p

PS: 以上内容是小编个人学习代码笔记分享，仅供参考学习，欢迎大家一起交流学习。

代码及源数据获取：https://github.com/wzsBio/code

温馨提示

如果你喜欢本文，请分享到朋友圈，想要获得更多信息，请关注我。