前往小程序,Get更优阅读体验!
立即前往
首页
学习
活动
专区
圈层
工具
发布
首页
学习
活动
专区
圈层
工具
MCP广场
社区首页 >专栏 >Scanpy 分析 3k PBMCs:寻找 marker 基因

Scanpy 分析 3k PBMCs:寻找 marker 基因

作者头像
数据科学工厂
发布于 2025-06-09 03:54:04
发布于 2025-06-09 03:54:04
7800
代码可运行
举报
运行总次数:0
代码可运行

引言

图片
图片

本系列讲解 使用Scanpy分析单细胞(scRNA-seq)数据教程[1],持续更新,欢迎关注,转发!

寻找 marker 基因

来给每个细胞簇里差异表达明显的基因排个序。通常情况下,要是 AnnData .raw 属性已经提前设置好了,就会拿它来用。最便捷快速的方式就是做** t 检验**。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
sc.tl.rank_genes_groups(adata, "leiden", method="t-test")
sc.pl.rank_genes_groups(adata, n_genes=, sharey=False)
图片
图片
代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
sc.settings.verbosity =   # reduce the verbosity

Wilcoxon 秩和检验得出的结果跟 t 检验差不多。不过,要是要发文章,更建议用 Wilcoxon 秩和检验。另外,还可以考虑用一些更厉害的差异检验工具包,比如 MAST、limma、DESeq2,要是用 python 的话,还有新出的 diffxpy。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
sc.tl.rank_genes_groups(adata, "leiden", method="wilcoxon")
sc.pl.rank_genes_groups(adata, n_genes=, sharey=False)
图片
图片

保存结果:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
adata.write(results_file)

作为一种替代方案,可以用逻辑回归来给基因排个序。主要的不同点是,这次用的是多变量方法,而常规的差异检验是单变量的。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
sc.tl.rank_genes_groups(adata, "leiden", method="logreg", max_iter=)
sc.pl.rank_genes_groups(adata, n_genes=, sharey=False)
图片
图片

除了 IL7R(只有 t 检验能找到)和 FCER1A(只有其他两种方法能找到)之外,其他标记基因在所有方法里都能找出来。

图片
图片

再定义一个标记基因的列表,方便后面用。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
marker_genes = [
    *["IL7R", "CD79A", "MS4A1", "CD8A", "CD8B", "LYZ", "CD14"],
    *["LGALS3", "S100A8", "GNLY", "NKG7", "KLRB1"],
    *["FCGR3A", "MS4A7", "FCER1A", "CST3", "PPBP"],
]

把之前保存有 Wilcoxon 秩和检验结果的那个对象重新加载进来。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
adata = sc.read(results_file)

把每个簇 0、1、……、7 排名前 10 的基因在一个数据框里展示出来。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
pd.DataFrame(adata.uns["rank_genes_groups"]["names"]).head()
图片
图片

弄一个表格,把得分和组别都列出来。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
result = adata.uns["rank_genes_groups"]
groups = result["names"].dtype.names
pd.DataFrame(
    {
        f"{group}_{key[:]}": result[key][group]
        for group in groups
        for key in ["names", "pvals"]
    }
).head()
图片
图片

和单个簇对比一下:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
sc.tl.rank_genes_groups(adata, "leiden", groups=["0"], reference="1", method="wilcoxon")
sc.pl.rank_genes_groups(adata, groups=["0"], n_genes=)
图片
图片

要是想仔细看看某个特定组的情况,可以用 sc.pl.rank_genes_groups_violin

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
sc.pl.rank_genes_groups_violin(adata, groups="0", n_genes=)
图片
图片

重新把带有差异表达计算结果(也就是跟其他组对比出来的差异表达)的对象加载进来:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
adata = sc.read(results_file)

sc.pl.rank_genes_groups_violin(adata, groups="0", n_genes=)
图片
图片

要是想看某个特定基因在不同组之间的变化,就用下面这个方法。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
sc.pl.violin(adata, ["CST3", "NKG7", "PPBP"], groupby="leiden")
图片
图片

给细胞类型做上标记。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
new_cluster_names = [
    "CD4 T",
    "B",
    "FCGR3A+ Monocytes",
    "NK",
    "CD8 T",
    "CD14+ Monocytes",
    "Dendritic",
    "Megakaryocytes",
]
adata.rename_categories("leiden", new_cluster_names)

sc.pl.umap(
    adata, color="leiden", legend_loc="on data", title="", frameon=False, save=".pdf"
)
图片
图片

现在细胞类型已经标记好了,来看看标记基因的可视化效果:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
sc.pl.dotplot(adata, marker_genes, groupby="leiden");
图片
图片

还有一种很简洁的小提琴图:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
sc.pl.stacked_violin(adata, marker_genes, groupby="leiden");
图片
图片

在这个分析过程中,AnnData 产生了一些注释:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
adata
图片
图片
代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
# `compression='gzip'` saves disk space, and slightly slows down writing and subsequent reading
adata.write(results_file, compression="gzip")

可以用 h5ls大致看看文件内容。文件格式以后可能还会进一步优化。不过不用担心,所有的读取功能都会保持向后兼容。

如果你想把文件分享给那些只打算用来做可视化的小伙伴,可以通过删除密集的缩放和校正数据矩阵来减小文件大小。不过别担心,文件里还是有 adata.raw 里用于可视化的原始数据的。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
adata.raw.to_adata().write("./write/pbmc3k_withoutX.h5ad")

Reference

[1] 

Source: https://scanpy.readthedocs.io/en/stable/tutorials/basics/clustering-2017.html

本文参与 腾讯云自媒体同步曝光计划,分享自微信公众号。
原始发表:2025-06-08,如有侵权请联系 cloudcommunity@tencent.com 删除

本文分享自 冷冻工厂 微信公众号,前往查看

如有侵权,请联系 cloudcommunity@tencent.com 删除。

本文参与 腾讯云自媒体同步曝光计划  ,欢迎热爱写作的你一起参与!

评论
登录后参与评论
0 条评论
热度
最新
推荐阅读
目录
  • 引言
  • 寻找 marker 基因
领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档