最近做细菌比较基因组学分析，跑了一下Roary软件，出来了一些结果和图，看了看Roary软件github上的作图代码，使用的python语言对系统进化树和gene presence/absence数据的可视化，想要调整图片的一些细节逼着自己去看python代码，谷歌了 Biopython中的phylo模块的使用文档，依然是懵逼的，我想让系统进化树显示support值，显示并调节树上的菌株名字的大小，还是决定用自己稍微熟悉的R语言去解决问题。
就是类似这张图要修改。
发现Y叔的ggtree可以实现类似的图。
但我学习的时候没有找到这个图的完整数据，搞不清数据结构是怎么样的，所以我决定用ggtree展示系统进化树，ggplot2展示矩阵数据，然后使用patchwork两个图拼在一起。
1. 系统进化树可视化
#########系统进化树可视化#######
rm(list = ls())
getwd()
setwd("C:\\Users\\Administrator\\Desktop\\R_work\\07_BD177_genome_work")
list.files()
#载入包
library("ggtree")
library("treeio")
library("tidyverse")
#树文件导入
BD177_tree <- read.newick("accessory_binary_genes.fa.newick",node.label = "support")
#系统进化树可视化
p1 <- ggtree(BD177_tree, branch.length = "none") +
geom_tiplab(size=2, align=TRUE, linesize=.5) +
xlim(0,30) +
geom_nodelab(aes(subset = !isTip,label=round((node/232), 2)),hjust=0.5,color="black", size = 2.5)
p1
2. gene presence/absence矩阵数据的可视化
#矩阵数据导入
gene_pre_ab <- read.csv("gene_presence_absence_heatmap_tidydata.csv", sep = ",")
#宽数据变长数据
gene_pre_ab_long <- gene_pre_ab %>%
pivot_longer(cols = -strain, names_to = "gene", values_to = "value")
#将ggtree生成的系统进化树所有信息存入
p1b = ggplot_build(p1)
head(p1b$data)
#将gene_pre_ab_long数据中的strain与系统进化树数据中的label进行一一对应
gene_pre_ab_long = gene_pre_ab_long %>%
mutate(strain = factor(strain, levels=p1b$data[[3]] %>% arrange(y) %>% pull(label)))
#矩阵可视化
p2 <- ggplot(gene_pre_ab_long, aes(x = fct_reorder(gene, value,.desc = T), y = strain)) +
geom_raster(aes(fill=value)) +
scale_fill_gradient(low="#f6fafe", high = "dodgerblue3") +
labs(x="", y = "", title="") +
theme_bw() +
theme(axis.text.x=element_blank(),
axis.text.y=element_blank(),
plot.title=element_text(size=5)) +
guides(fill = FALSE )
p2
3. 拼图
#拼图与导出
library("patchwork")
library("export")
p1 + p2 + plot_layout(widths = c(1.5, 1))
graph2ppt(file="BD177_tree_heatmap.pptx", width=7, height=10)
Y叔的ggtree包还需要进一步探索，看ggtree的帮助文档发现功能十分强大。
