超简单的比较基因组学Ka/Ks进化速率图绘制

看文献学绘图

近日，武汉大学李家儒老师课题组在国际著名杂志The Plant Journal上发表了题为“Integrated evolutionary pattern analyses reveal multiple origins of steroidal saponins in plants”的研究论文，揭示了植物甾体皂苷起源与进化机制。

其中文章中值得对我们分析基因组学比较实用的是图二：基因表达、氨基酸使用频率及基因进化速率的相关性分析。

这个图片中的b图就是参与甾体皂苷生物合成的基因kaks选择压力图的一个绘制（我们一般可以绘制一个基因家族不同物种间或者物种内的进化速率分析）。

代码实操绘图

利用R语言的ggplot2包就可以比较丝滑的进行绘制美美的kaks选择压力进化图了，如果是对R不熟悉或者刚入门的建议看看Y叔（余光创教授）的R实战书这本书《R实战：系统发育树的数据集成操作及可视化》，刚推出时候，卖断货N多次，现在价格只要半价，非常推荐。

代码实操部分

01

代码第一部分rm(list=ls())

library(ggplot2)

library(reshape2)

data <- read.table("all.ks",header = T,sep="\t")

data = na.omit(data)

head(data)

示例文件在生信漫谈交流群，可以直接下载，具体结果如下：

p1 <- ggplot(data,aes(x=KaKs ,y=..density.., color = Species, fill=Species)) +

geom_histogram(position = "identity" ,alpha = 0.5) +

geom_density(alpha = 0) +

xlim(0, 5) +

xlab("KaKs value") +

theme_classic() +

#scale_fill_brewer(palette = 'Paired') +

theme(legend.title=element_blank(),legend.position = c(0.9,0.9))

pdf(file=paste("Kaks_density.pdf", sep = "."), width = 10, height = 7)

p1

可以直接得到如下图所示kaks的选择压力图：

如果你对柱形图不感兴趣，想得到文章类似的图可以改下代码部分，如下所示：#如果只要线条型而不是柱状图

p2 <- ggplot(data,aes(KaKs,color = Species, fill = Species)) +

geom_line(stat="density",alpha=0.8) + xlim(0,5) + xlab("KaKs") +theme_classic() +

theme(axis.title = element_text(size=16),axis.text=element_text(size=16))

pdf(file=paste("newKaks_density.pdf", sep = "."), width = 10, height = 7)

p2

dev.off()

02

代码调整部分如果小伙伴还想调整主题风格，theme_grey()------灰色背景，不带x、y轴线，带白色网格线等

p1 <- ggplot(data,aes(x=KaKs ,y=..density.., color = Species, fill=Species)) +

geom_histogram(position = "identity" ,alpha = 0.5) +

geom_density(alpha = 0) +

xlim(0, 5) +

xlab("KaKs value") +

theme_grey()+

#scale_fill_brewer(palette = 'Paired') +

theme(legend.title=element_blank(),legend.position = c(0.9,0.9))

pdf(file=paste("Kaks_density.pdf", sep = "."), width = 10, height = 7)

p1

theme_bw() 主题-----变色背景、带全方位x、y轴线（上下左右），带灰色网格线。

p1 <- ggplot(data,aes(x=KaKs ,y=..density.., color = Species, fill=Species)) +

geom_histogram(position = "identity" ,alpha = 0.5) +

geom_density(alpha = 0) +

xlim(0, 5) +

xlab("KaKs value") +

theme_bw()+

#scale_fill_brewer(palette = 'Paired') +

theme(legend.title=element_blank(),legend.position = c(0.9,0.9))

pdf(file=paste("Kaks_density.pdf", sep = "."), width = 10, height = 7)

p1