RNA-seq (4)——差异分析可视化，ggplot2 绘制高颜值火山图

以解牛之法析生信，观微雀之形览科研。

乔粒说：在上篇，我们进行了使用DESeq2对airway数据的进行转录组差异分析，RNA-seq(3):用DESeq2做差异表达分析-以airway数据为例

这篇我们用火山图呈现转录组差异分析的结果，这里使用RStudio中的ggplot包。

火山图怎么看：

以下图示例：

筛选标准：阈值：|log2FC| > 1 且 padj < 0.05

①横轴 X：log₂ Fold Change（对数化表达差异倍数） 横轴代表基因在地塞米松处理组与对照组间的表达差异倍数，并经以 2 为底的对数转换： log 2FC＞0（竖虚线右侧）：上调基因，药物处理后表达量高于对照组； log 2FC＜0（竖虚线左侧）：下调基因，药物处理后表达量低于对照组； -1 < log 2FC < 1：表达变化幅度很小，一般不视为显著差异基因。

②纵轴 Y：-log₁₀(padj) 校正后 P 值负对数 Y 值越大 → 校正后 P 值越小 → 基因表达差异统计学显著性越高； 底部水平虚线 y=-log10(0.05)：显著性分界线； 点在虚线上方：padj < 0.05，差异有统计学意义； 点在虚线下方：padj > 0.05，差异无统计学意义。

③图例标注： -log10(padj)，颜色直接对应显著性高低：青、黄绿、橙：显著性逐步升高；颜色越偏红，该基因差异越可靠。

④显著上调基因：log 2FC>1&padj<0.05

⑤显著下调基因：log 2FC<−1&padj<0.05

⑥底部水平虚线 ：y=-log10(0.05)：显著性分界线；

⑦两条竖直虚线： x= -1、x= 1：常规差异倍数阈值，代表表达变化 2 倍；

火山图图片展示：

下边我们来示例如何绘制这些火山图吧！

四、火山图绘制

1. 提取数据

此步提取数据及绘图，需要上一章节的运行结果，如果是刚刚进来的小伙伴可以参考运行： RNA-seq(3):用DESeq2做差异表达分析-以airway数据为例 得到差异分析的结果再来运行此章节。

#提取数据
#加载包
library(ggplot2)
library(ggrepel)

# 3.1添加基因符号映射
# 从 airway 对象的 rowData 提取基因符号，并命名为与行名对应的向量
gene_symbols <- rowData(airway)$gene_name
names(gene_symbols) <- rownames(airway)

# 在 res_df 中新增一列 gene_symbol
res_df$gene_symbol <- gene_symbols[rownames(res_df)]
res_df$gene_symbol[is.na(res_df$gene_symbol)] <- rownames(res_df)[is.na(res_df$gene_symbol)]

# 3.2 筛选标注基因（各取前5）
# 分别提取上调和下调的显著基因（去除 NA）
up_genes <- res_df[!is.na(res_df$padj) & res_df$change == "Up", ]
down_genes <- res_df[!is.na(res_df$padj) & res_df$change == "Down", ]

# 按 padj 排序，各取前5（不足则全取）
top_up <- up_genes[order(up_genes$padj), ][1:min(5, nrow(up_genes)), ]
top_down <- down_genes[order(down_genes$padj), ][1:min(5, nrow(down_genes)), ]

# 合并标注数据框（包含 gene_symbol 列）
top_label <- rbind(top_up, top_down)
head(top_label)

这里我们可以看到，构建的数据包含了gene_symbal、pvalue、padj、change列。

我们使用的阈值为：|log2FC| > 1 且 padj < 0.05

这里我们使用padj来绘制火山图，可以根据不同的情况使用pvalue或者padj。

2.基础火山图

绘制火山图，展示差异基因。

#基础火山图绘制
p <- ggplot(res_df, aes(x = log2FoldChange, y = –log10(padj), color = change)) +
geom_point(size = 1.5, alpha = 0.6) +
scale_color_manual(values = c("Up" = "#E41A1C",
"Down" = "#377EB8",
"Stable" = "grey70")) +
geom_vline(xintercept = c(–1, 1), linetype = "dashed", color = "grey40") +
geom_hline(yintercept = –log10(0.05), linetype = "dashed", color = "grey40") +
theme_classic() +
theme(plot.title = element_text(hjust = 0.5)) +
labs(x = "log2 Fold Change",
y = "-log10(adjusted p-value)",
color = "Regulation") +
ggtitle("Dexamethasone vs Control (airway dataset)")

print(p)
ggsave("volcano_plot.png", p, width = 8, height = 6,dpi=600)

3. 火山图：替换颜色+显示TOP5 Gene_namen

这个不够美观，换个颜色，显示TOP5的差异基因名字（包含上调下调），

# 换个颜色
p2 <- ggplot(res_df, aes(x = log2FoldChange, y = –log10(padj), color = change)) +
# 点：调大尺寸，增强透明度，增加边框（stroke）使点更清晰
geom_point(size = 3, alpha = 0.7, stroke = 0.2, shape = 16) +

# 颜色：上调用更亮的红色，下调用深蓝，稳定用浅灰
scale_color_manual(values = c("Up" = "#e94234", # 橙色红
"Down" = "#269846", # 深蓝
"Stable" = "#d3d3d3"),
breaks = c("Up", "Down") # 图例只显示 Up 和 Down
) +

# 阈值线：加粗并调色
geom_vline(xintercept = c(–1, 1), linetype = "dashed", color = "grey30", size = 0.5) +
geom_hline(yintercept = –log10(0.05), linetype = "dashed", color = "grey30", size = 0.5) +

# 基因标签：调大字体，加粗，调整颜色为黑色，限制重叠
geom_text_repel(data = top_label,
aes(label = gene_symbol),
size = 4,
fontface = "bold",
color = "black",
box.padding = 0.5,
point.padding = 0.3,
max.overlaps = 30,
segment.color = "grey50",
segment.size = 0.3) +

# 主题：经典白底，调整标题和图例
theme_classic(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5, size = 16, face = "bold"),
axis.title = element_text(size = 14, face = "bold"),
axis.text = element_text(size = 12),
legend.position = "right",
legend.title = element_text(size = 12, face = "bold"),
legend.text = element_text(size = 11),
# 可选：添加轻微网格线
panel.grid.major = element_line(color = "grey95", size = 0.3),
panel.grid.minor = element_blank()
) +
coord_cartesian(xlim = c(–10, 10))+
labs(x = expression("log"[2] ~ "Fold Change"),
y = expression("-log"[10] ~ "(adjusted p-value)"),
color = "Regulation") +
ggtitle("Dexamethasone vs Control (airway dataset)")

print(p2)

ggsave("volcano_plot_2.png", p2, width = 8, height = 6,dpi=600)

4.火山图绘制：显著渐变色

火山图绘制：这个火山图还是不够好看，加上显著渐变色。

注：代码中放置了一些配色，可以根据需要运行（删除对应的colors前的“#”，并且将上部运行的colors前加上“#”即可）

#火山图绘制：+显著渐变色
# 配色：蓝色(不显著) → 黄色(中等) → 红色(极显著)
library(ggplot2)
library(ggrepel)

p3 <- ggplot(res_df, aes(x = log2FoldChange, y = –log10(padj))) +
# 核心：颜色根据显著性(-log10(padj))渐变
geom_point(aes(color = –log10(padj)),
size = 5,
alpha = 1.0) +
geom_point(
data = top_label,
aes(color = –log10(padj)),
size = 4 , # 点的大小，按需调整数值
stroke = 3, # 白色描边宽度
shape = 16, # 实心圆点
alpha = 1
) +

# 自定义渐变色（放置了好几种版本的，使用时删除前边的”#"即可）
scale_color_gradientn(
#colors = c("#2E86AB", "#A23B72", "#F18F01", "#C73E1D"),
#colors = c("#0000FF","#00FFFF","#FFFF00","#FF0000"),# 经典配色
#colors = c("#1f78b8","#a6cee3","#fb9a99","#e31a1c"), #冷色调高级蓝紫
# colors = c("#4A6FE3","#85A0FF","#FF9F43","#E63946"),#莫兰迪低饱和
#colors = c("#3a1c71","#d76d77","#ffaf7b","#ff0000"),
#colors = c("#253494", "#41b6c4", "#fed976", "#fc4e2a", "#b10026") ,
# colors = c( "#2166AC", "#67A9CF", "#D9D9D9", "#F4A582", "#B2182B" ),
colors= c( "#0E2A85", "#52CFC0", "#FFE02E", "#F55E36", "#D12626" ),
name = "-log10(padj)"
) +

# 阈值虚线
geom_vline(xintercept = c(–1, 1), linetype = "dashed", color = "grey30", linewidth = 0.5) +
geom_hline(yintercept = –log10(0.05), linetype = "dashed", color = "grey30", linewidth = 0.5) +

# 标注之前筛选的Top基因（基因符号）
geom_text_repel(
data = top_label,
aes(label = gene_symbol),
size = 4, fontface = "bold", color = "black",
box.padding = 0.5, max.overlaps = 30
) +

# 坐标轴与主题
coord_cartesian(xlim = c(–10, 10)) +
theme_classic(base_size = 14) +
theme(
plot.title = element_text(hjust = 0.5, size = 16, face = "bold"),
legend.position = "right"
) +
labs(
x = expression(log[2]~FoldChange),
y = expression(–log[10]~(padj)),
title = "Dexamethasone vs Control (airway dataset)"
)

# 出图
print(p3)
# 保存高清图片
ggsave("volcano_gradient.png", p3, width = 10, height = 8, dpi = 600)

5.火山图绘制：+上下调基因数目统计，xy轴自适应，细胞大小变化

添加上上调下调基因统计，xy轴自适应，细胞大小颜色随数值大小变化，

#添加上上调下调基因统计，xy轴自适应，细胞大小颜色随数值大小变化
library(ggplot2)
library(ggrepel)
library(dplyr)

# 数据预处理
# 过滤掉padj为NA的基因，避免绘图报错
plot_data <- res_df %>% filter(!is.na(padj))

# 统计上下调基因数量，用于图上标注
up_num <- sum(plot_data$change == "Up")
down_num <- sum(plot_data$change == "Down")

# 定义参考图同款配色：蓝(下调)-灰白-红(上调)
fc_colors <- c("#2166AC", "#67A9CF", "#D9D9D9", "#F4A582", "#B2182B")

# 自适应坐标轴+标注位置计算
# X轴：取最大绝对值，实现左右对称
max_fc <- max(abs(plot_data$log2FoldChange), na.rm = TRUE)

# Y轴：计算数据点的最高位置（所有基因点的最大y值）
max_pval <- max(–log10(plot_data$padj), na.rm = TRUE)
# Y轴上限：在最高点基础上多留10个单位，充足空白放标注
y_upper <- max_pval + 10

# 标注分层（全部在最高点之上，绝对不与细胞点重叠）
#y_arrow <- max_pval + 40 # 箭头在最高点上方40个单位（最靠上）
y_label <- max_pval + 35 # Down/Up文字在最高点上方35个单位
y_num <- max_pval + 25 # 数字在最高点上方25个单位
# ———- 绘图代码 ———-
p_advanced <- ggplot(plot_data, aes(x = log2FoldChange, y = –log10(padj))) +
# 底层：所有基因点 —— 颜色映射logFC，大小映射显著性
geom_point(
aes(color = log2FoldChange, size = –log10(padj)),
alpha = 1.3 #细胞的颜色深浅
) +

# 颜色渐变：蓝-白-红，对称色阶
scale_color_gradientn(
colors = fc_colors,
limits = c(–max(abs(plot_data$log2FoldChange)), max(abs(plot_data$log2FoldChange))),
name = "logFC"
) +

# 点大小渐变
scale_size_continuous(
range = c(0.3, 8),
name = "-log10(padj)"
) +

# 阈值虚线
geom_vline(xintercept = c(–1, 1), linetype = "dashed", color = "black", linewidth = 0.6) +
geom_hline(yintercept = –log10(0.05), linetype = "dashed", color = "black", linewidth = 0.6) +

# Top标注基因
geom_point(
data = top_label,
color = "black",
fill = "#EC5F5F" ,
shape = 21,
size = 7,
stroke = 0.7
) +

# 基因名称标注
geom_text_repel(
data = top_label,
aes(label = gene_symbol),
size = 4,
color = "black",
box.padding = 0.8,
point.padding = 0.4,
max.overlaps = 20,
segment.color = "black",
segment.size = 0.5
) +

# ===================== 修改：自适应上下调标注（无固定坐标） =====================
annotate("text", x = –4, y = y_label, label = "Down", color = "#67A9CF", size = 4, fontface = "bold") +
annotate("text", x = –4, y = y_num, label = down_num, color = "#67A9CF", size = 4, fontface = "bold") +
annotate("text", x = 4, y = y_label, label = "Up", color = "#F4A582", size = 4, fontface = "bold") +
annotate("text", x = 4, y = y_num, label = up_num, color = "#F4A582", size = 4, fontface = "bold") +
# 自适应箭头
# annotate("segment", x = -3.5, xend = -4.2, y = y_arrow, yend = y_arrow, color = "#67A9CF", size = 1, arrow = arrow(length = unit(0.2, "cm"))) +
# annotate("segment", x = 3.5, xend = 4.2, y = y_arrow, yend = y_arrow, color = "#F4A582", size = 1, arrow = arrow(length = unit(0.2, "cm"))) +

# 图例顺序
guides(
size = guide_legend(order = 1),
color = guide_colorbar(order = 2)
) +

# 主题
theme_bw(base_size = 14) +
theme(
panel.grid = element_line(color = "grey90", linewidth = 0.3),
legend.position = "right",
legend.box = "vertical",
axis.title = element_text(face = "bold"),
plot.title = element_text(hjust = 0.5, size = 16, face = "bold")
) +

coord_cartesian(
xlim = c(–max_fc, max_fc), # X轴自适应对称
ylim = c(0, y_upper+30) # Y轴自适应+顶部留白
) +

# 坐标轴标题
labs(
x = "log2 Fold Change",
y = "-log10(padj)",
title = "Dexamethasone vs Control (airway dataset)"
)

# 可视化+保存
print(p_advanced)
ggsave("advanced_volcano_1.png", p_advanced, width = 10, height = 8, dpi = 600)

乔粒说：今天的差异分析绘制火山图展示就到这里啦，欢迎批评指正，我们下期见！

以解牛之法析生信，观微雀之形览科研。