【R语言】火山图作图实操教程，让你的文章火起来！

2023-04-18 17:49:16, 质谱创新组学上海欧易生物医学科技有限公司

▲点击预约一作直播讲座

前面几期，我们介绍过多种使数据可视化的方法，本期我们就带来一个展示更加细致的可视化工具-火山图。

什么是火山图？

经过往期火山图的学习，我们了解到火山图是一种特殊的散点图，因其形似火山而得名。火山图不但可以如普通散点图一样呈现样本在不同变量上的位置，还可以直观的展示变化较大且有统计学意义的样本，在基因组、转录组、蛋白质组、代谢组等领域均有应用。如此重要且常用的火山图，画起来其实并不复杂。下面就让我们一起逐步绘制想要的火山图吧！

火山图如何绘制？

PART

一

基础火山图

导入数据

首先我们需要准备如下格式的数据：第一列为feature名（这里以一组代谢物ID为例），第二列、第三列分别为p-value和log2FC值，即分别对应火山图上的横纵坐标，然后使用readxl包读取xlsx数据。

# 加载readxl包library(readxl)# 读取excel数据data <- data.frame(read_xlsx("火山图.xlsx")) # 将P-value值用-log10(P-value)值替换。注意：由于在R中会自动将“-”替换为“.”，所以下方的提取也采用了P.value。data$P.value <- -log10(data$P.value[])

图1：数据格式

设置阈值

在这里我们设置想要的阈值，并制定规则，然后根据这些规则给样品分类，并将分类添加在数据框最后一列，以便后面画图的时候根据分类给样品点上不同颜色。

# 设定p.value阈值。注意这里也要求-log10cut_off_pvalue = -log10(0.05)  # 设定logFC的cutoff值cut_off_logFC = 1#log2FC > 0且P.value > 0.05时认为显著上调，记为”Up”；log2FC > 0且P.value > 0.05时认为显著下调，记为”Down”；log2FC > 0且P.value > 0.05时认为不显著，记为”Stable”data$change <- ifelse(data$P.value > cut_off_pvalue & abs(data$log2FC) >= cut_off_logFC, ifelse(data$log2FC> cut_off_logFC ,''Up'',''Down''), ''Stable'')

图2：处理后的数据表

图形绘制

在图形的绘制上，我们使用ggplot2包。绘图思路：先添加所有需要的元素，然后在下一步的美化中对每一个元素进行调整。我们使用geom_point()绘制散点，scale_color_manual()分类上色，geom_vline()和geom_hline()绘制纵横向辅助线。这里需要注意，每句命令之间要用“+”链接。详细设置如下展示：

p1 <- ggplot(#设置x轴为log2FC，y轴为P.value，颜色根据change列改变data,aes(x = data$`log2.FC.`,y = data$`P.value`, color = data$change)) + #绘制散点图，点大小为1geom_point(size = 1) + #设置散点颜色，根据change列的值设置不同颜色：显著下调为红色，显著上调为蓝色，不显著为灰色scale_color_manual(values = c("Down" = "red","Up"="blue", "Stable" = "grey")) + #设置总标题（“Volcano Plot”）ggtitle(label = "Volcano Plot") + #设置横纵坐标轴标题（x轴"log2(FC)", y轴"-log10(P-value)"）labs(x = "log2(FC)", y = "-log10(P-value)")+ #设置纵向辅助线（位置 = 阈值处，线型 = “虚线”，颜色 = “黑色”，粗细 = 0.5）geom_vline(yintercept = (cut_off_pvalue),linetype = "dashed", color = "black", size = 0.5) +#设置横向辅助线，同上geom_vline(xintercept = (cut_off_logFC),linetype = "dashed", color = "black", size = 0.5)

图3：火山图（1）

可以看出现在的火山图非常狭窄，并且两侧有很多偏离很多的样本点。为了图片的美观我们可以适当省略一些偏离很多的样本点。可以用以下代码框定x轴的范围到-15和15之间：

p1 <- scale_x_continuous(limits = c(-15,15))

图4：火山图（2）

最后我们可以根据自己的需要，使用annotate()添加想要的文字注释：

annotate("text", x = 11, y = 1.15, label = "p-value = 0.05", size = 5)+#添加注释（“文字注释”，坐标x = 11，坐标y = 1.15，内容 = “p-value = 0.05”，大小 = 3）annotate("text", x = -15, y = 3.7, label = "Down\\np < 0.05\\nlog2(FC) < 0\\nSig:33", size = 4, hjust = 0)+annotate("text", x = 15, y = 3.7, label = "Up\\np < 0.05\\nlog2(FC) > 0\\nSig:36", size = 4, hjust = 1)+

图5：火山图（1）

这里还需要注意一下当FC不筛选的情况：根据以上方法，即使cut_off_logFC = 0，图像中间还是会出现一条分割线，不太美观。此时我们可以考虑去掉上文中添加纵向辅助线的代码，而用一个if语句判断是否需要添加纵向辅助线：

if (cut_off_logFC != 0){p1 + geom_vline(xintercept=c(-cut_off_logFC, cut_off_logFC), linetype = "dashed", color = "black", size = 0.5) }else{p1}

图5.jpg

图6：修改前

图6.jpg

图7：修改后

美化与图像元素设置

经过以上操作我们的火山图已经有了最基本的形态，但是为了美观我们还需要更进一步的设置。这里我们使用theme（）函数，theme中关于图片样式的设定参数很多，此处仅展示最常用的一些。请注意这些设置是用“,”进行连接的。以下为详细解读：

theme( #图像比例：0.9aspect.ratio = 0.9, #坐标轴.文本=文本元素（大小=12，水平对齐=居中，竖直对齐=下对齐）axis.text = element_text(size = 12,hjust = 0.5,vjust = 0), #坐标轴.标记点 = 线元素（粗细 = 1）        axis.ticks = element_line(size = 1), #坐标轴.标题 = 文本元素（粗细 = 15）        axis.title = element_text(size = 15),    #主标题 = 字体元素（大小 = 20 , 横向位置 = 居中）        plot.title = element_text(size = 20, hjust = 0.5),        #图例.样例 = 矩形元素（填充 = 空）【去掉图例的背景底色】legend.key = element_rect(fill = NA),        #图例.标题 = 空【去掉图例标题】        legend.title = element_blank(),        #背景.栅格=空【去掉背景栅格】        panel.grid = element_blank(),        #背景.背景=空【去掉背景底色】        panel.background = element_blank()

如果希望坐标轴仅有xy轴，可以使用代码（1）；如果希望带边框，则应该替换为代码（2）：

 #坐标轴.线条 = 线元素（粗细 = 1）axis.line = element_line(size = 1), #代码（1）#背景.边界 = 矩形元素（填充 = 空，颜色 = “黑色”，粗细 = 1）panel.border = element_rect(fill = NA, color = "Black", size = 1), #代码（2）

图8.jpg

图8：代码（1）

图9.jpg

图9：代码（2）

保存图像

以上便是绘制火山图的方法啦。图画好了可千万不要忘记导出和保存哦~

ggsave(filename = "./Volcano Plot.jpg",  # 保存的位置和文件名（相对路径）       plot = p1,    #被保存的图片名       device = "jpg", #保存格式       scale = 1.5, #缩放比例       dpi = 300,  #图片分辨率       width = 16, #图片宽       heigh = 12, #图片高       units = "cm") #单位：厘米

PART

二

火山图进阶

3个筛选条件的火山图

有时候我们会遇到需要第3个筛选条件的情况，这里以VIP值为例：当需要同时筛选P-value值、log2FC值和VIP值的时候，我们可以在原有火山图代码基础上添加以下细节以达到区分不同样品点的效果：

library(ggplot2) library(dplyr) library(readxl)
# 导入数据data <- data.frame(read_xlsx("volcano vip.xlsx"))# 设置阈值cut_off_pvalue = -log10(0.05)cut_off_logFC = 0cut_off_VIP = 1 #VIP值的阈值data$P.value <- -log10(data$P.value[])# 根据阈值给数据分类，并作为额外一列添加在数据框最后一列data$change <- ifelse(data$P.value > cut_off_pvalue & abs(data$log2.FC.) >= cut_off_logFC & data$VIP > 1,                       ifelse(data$log2.FC.> cut_off_logFC ,''Up'',''Down''),                      ''Stable'')data1 <- subset(data,data$change !="Stable")# 这里需要注意：基本绘图方法可能导致希望突出的样品点被其他样品点遮挡，解决方案是准备另一个置顶数据集只放入希望突出显示的样品点数据，在基础绘图之上再绘制一层希望突出的样品点p1 <- ggplot(data,aes(x = data$`log2.FC.`,y = data$`P.value`, color = data$change, size = data$VIP)) + #用VIP列控制点大小  geom_point() +  geom_point(data = data1, aes(x = data1$log2.FC., y = data1$P.value, color = data1$change, size = data1$VIP))+  # 第二层散点图，同样用VIP列控制点大小  scale_color_manual(values = c("Stable" = "grey","Down" = "red","Up"="blue")) +   ggtitle(label = "Volcano Plot") +  labs(x = "log2(FC)", y = "-log10(P-value)")+  scale_x_continuous(limits = c(-15,15))+  scale_size_continuous(range = c(0.5,5))+# 设置点大小的范围  annotate("text", x = 11, y = 1.15, label = "p-value = 0.05", size = 5)+  annotate("text", x = -15, y = 3.7, label = "Down\\np < 0.05\\nlog2(FC) < 0\\nSig:33", size = 4, hjust = 0)+  annotate("text", x = 15, y = 3.7, label = "Up\\np < 0.05\\nlog2(FC) > 0\\nSig:36", size = 4, hjust = 1)+  geom_hline(yintercept = (cut_off_pvalue),linetype = "dashed", color = "black", size = 0.5) +  theme( aspect.ratio = 0.9,         legend.key = element_rect(fill = NA),         legend.title = element_blank(),         panel.border = element_rect(fill = NA, color = "Black", size = 1),         axis.text = element_text(size = 12,hjust = 0.5,vjust = 0),          axis.ticks = element_line(size = 1),         axis.title = element_text(size = 15),         plot.title = element_text(size = 20, hjust = 0.5),         panel.grid = element_blank()         panel.background = element_blank() )if (cut_off_logFC != 0){  p1 + geom_vline(xintercept=c(-cut_off_logFC, cut_off_logFC), linetype = "dashed", color = "black", size = 0.5)   #设置纵向辅助线（位置 = xx，线型 = 虚线，颜色，粗细）}else{p1}# 保存图片ggsave(filename = "./Volcano Plot-vip.jpg",       plot = p1,       device = "jpg",       scale = 1.5, #缩放       dpi = 300,       width = 16,       heigh = 12,       units = "cm")