引言
本系列讲解 使用Scanpy分析单细胞（scRNA-seq）数据教程，持续更新，欢迎关注，转发！
寻找 marker 基因
来给每个细胞簇里差异表达明显的基因排个序。通常情况下，要是 AnnData的.raw 属性已经提前设置好了，就会拿它来用。最便捷快速的方式就是做 t 检验。
sc.tl.rank_genes_groups(adata, "leiden", method="t-test")
sc.pl.rank_genes_groups(adata, n_genes=25, sharey=False)
sc.settings.verbosity = 2 # 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=25, sharey=False)
保存结果：
adata.write(results_file)
作为一种替代方案，可以用逻辑回归来给基因排个序。主要的不同点是，这次用的是多变量方法，而常规的差异检验是单变量的。
sc.tl.rank_genes_groups(adata, "leiden", method="logreg", max_iter=1000)
sc.pl.rank_genes_groups(adata, n_genes=25, 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(5)
弄一个表格，把得分和组别都列出来。
result = adata.uns["rank_genes_groups"]
groups = result["names"].dtype.names
pd.DataFrame(
{
f"{group}_{key[:1]}": result[key][group]
for group in groups
for key in ["names", "pvals"]
}
).head(5)
和单个簇对比一下：
sc.tl.rank_genes_groups(adata, "leiden", groups=["0"], reference="1", method="wilcoxon")
sc.pl.rank_genes_groups(adata, groups=["0"], n_genes=20)
要是想仔细看看某个特定组的情况，可以用 sc.pl.rank_genes_groups_violin。
sc.pl.rank_genes_groups_violin(adata, groups="0", n_genes=8)
重新把带有差异表达计算结果（也就是跟其他组对比出来的差异表达）的对象加载进来：
adata = sc.read(results_file)
sc.pl.rank_genes_groups_violin(adata, groups="0", n_genes=8)
要是想看某个特定基因在不同组之间的变化，就用下面这个方法。
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")
