python版本的差异表达基因火山图绘制

生信技能树

发布于 2025-06-23 12:27:54

6440

文章被收录于专栏：生信技能树生信技能树

matplotlib库是python中使用最广泛的绘图库勒，功能强大，虽然比不上R的绘图，但是还需要学一下~ 本次就使用 matplotlib库来绘制火山图

Matplotlib库官网给你：https://matplotlib.org/

示例

火山图的本质是一个散点图，先看看Matplotlib库官网的案例散点图绘制：

https://matplotlib.org/stable/plot_types/basic/scatter_plot.html#sphx-glr-plot-types-basic-scatter-plot-py

导入python模块：

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
plt.style.use('_mpl-gallery')

生成随机数据：

# make the data
np.random.seed()
x =  + np.random.normal(, , )
y =  + np.random.normal(, , len(x))
# size and color:
sizes = np.random.uniform(, , len(x))
colors = np.random.uniform(, , len(x))
x
y

绘图：

# plot
# 创建画布并设置大小
fig, ax = plt.subplots(figsize=(, ))  # 设置画布大小为 8x6 英寸

ax.scatter(x, y, s=sizes, c=colors, vmin=, vmax=)

ax.set(xlim=(, ), xticks=np.arange(, ),
       ylim=(, ), yticks=np.arange(, ))

plt.show()

绘制火山图

知道了上面的函数使用后，再结合kimi，很快就能得到一个火山图啦！

首先是拿到绘图的数据使用前面做过的一个芯片数据的差异结果吧：2万个基因少一半也不影响最后的差异分析富集结果啊？

数据为：https://ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE17351，拿到之后进行芯片预处理并做差异表达分析得到一个差异结果，或者微信找我发你：Biotree123。

或者百度云盘：链接: https://pan.baidu.com/s/1sAXzlxs4jU24ZcAW4LThqQ?pwd=uavp 提取码: uavp 。

读取差异基因结果：

# 读取差异分析结果
# 记住文件路径改成自己的
# index_col=0：第一列读取为行名
deg = pd.read_csv("./GSE17351/DEG.csv", index_col=)
deg.head()

总共有14080个基因，9列：

# 查看有多少行，有多少列
deg.shape
# (14080, 9)

计算 -log10(pvalue)

# 数据整理
deg['-log10(pval)'] = -np.log10(deg['P.Value'])  # 计算 -log10(pval)
print(deg.head())

# 统计差异个数
deg.g.value_counts()

获取标记的top基因：

# 按 -log10(p-value) 降序排序
deg = deg.sort_values(by='P.Value', ascending=True)

# 提取前10个显著基因
top_10_genes = deg.head()
top_10_genes

top_10_genes['P.Value']

绘图：

# plot
# 创建画布并设置大小
fig, ax = plt.subplots(figsize=(, ))  # 设置画布大小为 8x6 英寸
ax.scatter(x=deg.loc[deg.g=="stable",'logFC'].values, y=deg.loc[deg.g=="stable",'-log10(pval)'].values, s=, c='grey',alpha=0.9)

# 标记显著差异表达的基因
ax.scatter(x=deg.loc[deg.g=="up",'logFC'].values, y=deg.loc[deg.g=="up",'-log10(pval)'].values, c='red', alpha=0.5, s=)
ax.scatter(x=deg.loc[deg.g=="down",'logFC'].values, y=deg.loc[deg.g=="down",'-log10(pval)'].values, c='blue', alpha=0.5, s=)

# 添加阈值线
ax.axhline(y=-np.log10(0.05), color='black', linestyle='--', linewidth=)
ax.axvline(x=, color='black', linestyle='--', linewidth=)
ax.axvline(x=-1, color='black', linestyle='--', linewidth=)

# 在图上标注前10个基因的符号
for i, row in top_10_genes.iterrows():
    ax.text(row['logFC'], row['-log10(pval)'], row['name'], fontsize=,ha='center', va='bottom')

# 添加标题和标签
ax.set_title('Volcano Plot of Differential Gene Expression', fontsize=)
ax.set_xlabel('Log2 Fold Change', fontsize=)
ax.set_ylabel('-Log10(p-value)', fontsize=)

# 添加图例
ax.legend()

plt.show()

结果如下：

基因有一些重叠到一起了，使用 adjustText 分开：

conda activate sc
pip install adjustText -i https://mirrors.tuna.tsinghua.edu.cn/pypi/web/simple

重新绘图：

from adjustText import adjust_text  # 需要安装 adjustText 库

# plot
# 创建画布并设置大小
fig, ax = plt.subplots(figsize=(, ))  # 设置画布大小为 8x6 英寸
ax.scatter(x=deg.loc[deg.g=="stable",'logFC'].values, y=deg.loc[deg.g=="stable",'-log10(pval)'].values, s=, c='grey',alpha=0.5)

# 标记显著差异表达的基因
ax.scatter(x=deg.loc[deg.g=="up",'logFC'].values, y=deg.loc[deg.g=="up",'-log10(pval)'].values, c='red', alpha=0.5, s=)
ax.scatter(x=deg.loc[deg.g=="down",'logFC'].values, y=deg.loc[deg.g=="down",'-log10(pval)'].values, c='blue', alpha=0.5, s=)

# 添加阈值线
ax.axhline(y=-np.log10(0.05), color='black', linestyle='--', linewidth=)
ax.axvline(x=, color='black', linestyle='--', linewidth=)
ax.axvline(x=-1, color='black', linestyle='--', linewidth=)

# 使用 adjustText 自动调整文本位置
texts = []
for i, row in top_10_genes.iterrows():
    text = ax.text(row['logFC'], row['-log10(pval)'], row['name'], fontsize=, ha='center', va='bottom', color='black')
    texts.append(text)

# 调整文本位置以避免重叠
adjust_text(texts, arrowprops=dict(arrowstyle='->', color='gray'))

# 添加标题和标签
ax.set_title('Volcano Plot of Differential Gene Expression', fontsize=)
ax.set_xlabel('Log2 Fold Change', fontsize=)
ax.set_ylabel('-Log10(p-value)', fontsize=)

# 添加图例
ax.legend()

plt.show()