本系列讲解 使用Scanpy
分析单细胞(scRNA-seq)数据教程[1],持续更新,欢迎关注,转发!
来给每个细胞簇里差异表达明显的基因排个序。通常情况下,要是 AnnData
的 .raw
属性已经提前设置好了,就会拿它来用。最便捷快速的方式就是做** t 检验**。
sc.tl.rank_genes_groups(adata, "leiden", method="t-test")
sc.pl.rank_genes_groups(adata, n_genes=, sharey=False)
sc.settings.verbosity = # reduce the verbosity
Wilcoxon 秩和检验得出的结果跟 t 检验差不多。不过,要是要发文章,更建议用 Wilcoxon 秩和检验。另外,还可以考虑用一些更厉害的差异检验工具包,比如 MAST、limma、DESeq2,要是用 python 的话,还有新出的 diffxpy。
sc.tl.rank_genes_groups(adata, "leiden", method="wilcoxon")
sc.pl.rank_genes_groups(adata, n_genes=, sharey=False)
保存结果:
adata.write(results_file)
作为一种替代方案,可以用逻辑回归来给基因排个序。主要的不同点是,这次用的是多变量方法,而常规的差异检验是单变量的。
sc.tl.rank_genes_groups(adata, "leiden", method="logreg", max_iter=)
sc.pl.rank_genes_groups(adata, n_genes=, sharey=False)
除了 IL7R(只有 t 检验能找到)和 FCER1A(只有其他两种方法能找到)之外,其他标记基因在所有方法里都能找出来。
再定义一个标记基因的列表,方便后面用。
marker_genes = [
*["IL7R", "CD79A", "MS4A1", "CD8A", "CD8B", "LYZ", "CD14"],
*["LGALS3", "S100A8", "GNLY", "NKG7", "KLRB1"],
*["FCGR3A", "MS4A7", "FCER1A", "CST3", "PPBP"],
]
把之前保存有 Wilcoxon 秩和检验结果的那个对象重新加载进来。
adata = sc.read(results_file)
把每个簇 0、1、……、7 排名前 10 的基因在一个数据框里展示出来。
pd.DataFrame(adata.uns["rank_genes_groups"]["names"]).head()
弄一个表格,把得分和组别都列出来。
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()
和单个簇对比一下:
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
。
sc.pl.rank_genes_groups_violin(adata, groups="0", n_genes=)
重新把带有差异表达计算结果(也就是跟其他组对比出来的差异表达)的对象加载进来:
adata = sc.read(results_file)
sc.pl.rank_genes_groups_violin(adata, groups="0", n_genes=)
要是想看某个特定基因在不同组之间的变化,就用下面这个方法。
sc.pl.violin(adata, ["CST3", "NKG7", "PPBP"], groupby="leiden")
给细胞类型做上标记。
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"
)
现在细胞类型已经标记好了,来看看标记基因的可视化效果:
sc.pl.dotplot(adata, marker_genes, groupby="leiden");
还有一种很简洁的小提琴图:
sc.pl.stacked_violin(adata, marker_genes, groupby="leiden");
在这个分析过程中,AnnData 产生了一些注释:
adata
# `compression='gzip'` saves disk space, and slightly slows down writing and subsequent reading
adata.write(results_file, compression="gzip")
可以用 h5ls
大致看看文件内容。文件格式以后可能还会进一步优化。不过不用担心,所有的读取功能都会保持向后兼容。
如果你想把文件分享给那些只打算用来做可视化的小伙伴,可以通过删除密集的缩放和校正数据矩阵来减小文件大小。不过别担心,文件里还是有 adata.raw
里用于可视化的原始数据的。
adata.raw.to_adata().write("./write/pbmc3k_withoutX.h5ad")
Reference
[1]
Source: https://scanpy.readthedocs.io/en/stable/tutorials/basics/clustering-2017.html
扫码关注腾讯云开发者
领取腾讯云代金券
Copyright © 2013 - 2025 Tencent Cloud. All Rights Reserved. 腾讯云 版权所有
深圳市腾讯计算机系统有限公司 ICP备案/许可证号:粤B2-20090059 深公网安备号 44030502008569
腾讯云计算(北京)有限责任公司 京ICP证150476号 | 京ICP备11018762号 | 京公网安备号11010802020287
Copyright © 2013 - 2025 Tencent Cloud.
All Rights Reserved. 腾讯云 版权所有