1. 第三篇要解决的核心问题：从跑通到跑明白

scanpy 单细胞分析走到第三步，绝大多数人已经能把官网那份外周血单个核细胞示例的数据读进来、做完质控、挑完高变基因、跑完 PCA 和邻居图、分出十几个 cluster、画出 UMAP 了。但接下来就卡住了：cluster 编号 0 到 15 摆在图上，谁是谁？这些群里哪几个基因能拿去做后续实验验证？这群细胞是从哪个群分化过来的？图想去投稿，配色和排版怎么调？这些问题官方文档里散落在五六个不同的示例页面，每个页面只讲一段，串不起来。

这篇内容的定位很明确：它是 scanpy 单细胞分析官网示例系列的第三篇，专门处理流程中后段那几个"跑是跑出来了，但结果不敢用"的环节，包括差异表达与 marker 基因筛选、细胞类型注释与基因集打分、PAGA 与扩散伪时间的轨迹推断、绘图系统的参数拆解与二次定制，最后附一份报错排查速查表。适合已经能独立跑通前两步、手上有一份自己数据的中级使用者，也适合刚接手别人 notebook、看不懂其中参数为什么这么设的人。我会把默认值、参数背后的统计假设、以及哪些参数改了会翻车，尽量讲透。

1.1 前置状态自检：你的 adata 里应该有什么

在往下走之前，先花两分钟确认一下手里的对象是否达标。这一步省不得，因为后面所有操作都建立在这个状态之上，缺一项就会在某个莫名其妙的报错上耗掉半小时。

打开 notebook 敲几行检查：

import scanpy as sc adata = sc.read_h5ad('pbmc_clustered.h5ad') print(adata) # 看 shape 和 obs/var 的字段 print(adata.obs.columns.tolist()) # 应该有 leiden 或 louvain print(adata.obsm.keys()) # 应该有 X_pca 和 X_umap print(adata.uns.keys()) # 应该有 neighbors、pca、umap print(adata.raw) # 是否保留了 log 化的原始矩阵

一个健康的中间态至少满足几点：adata.X是 log1p 归一化后的表达（稀疏 CSR 更好），adata.raw保存了全基因的 log 化矩阵（因为高变基因筛选后adata只剩下 2000 个左右基因，做 marker 检测时会漏掉大量低表达但有意义的基因），adata.obs里有一个分类变量记录分群结果，adata.obsm['X_pca']和X_umap都存在，adata.uns['neighbors']里的连通性矩阵没被清空。如果adata.raw是 None，后面用use_raw=True就会一路报KeyError，这是新手最常见的翻车点之一。

1.2 本篇四个模块的取舍逻辑

为什么是这四个模块，而不是继续讲降维、讲批次校正？因为按实际项目的推进顺序，聚类完成之后紧接着要回答的科学问题基本就这三类：这个群是什么细胞、这群细胞和其他群有什么表达差异、这群细胞处在发育或激活的哪个阶段。绘图系统单独拎出来讲，是因为它贯穿全程，而且官方示例里每个函数只给了最简调用，真正的调参经验得靠自己踩。

至于批次校正、空间转录组、多组学整合这些，它们更依赖具体实验设计和数据来源，放在单篇里讲只能是泛泛而谈，反而稀释了密度，所以这一篇不展开。四个模块之间不是并列关系，注释和差异表达在实操中是来回迭代的——先差异分析找 marker，再拿 marker 去注释，注释完可能发现某个群拆得太细，回头调分辨率重新聚类，再重新做差异分析。理解这个循环，比记住某个函数怎么写更重要。

提示：把每一步的结果都写回adata.obs或adata.uns，而不是扔进一个临时变量打印一下就完事。三天后回来复现的时候，你会感谢当时的自己。

2. 差异表达分析：rank_genes_groups 的参数与陷阱

sc.tl.rank_genes_groups是这一切的起点。它的官方签名看起来参数不算多，但真正影响结果是否可信的，其实集中在三四个地方。很多人直接无参调用，得到一张漂亮的热图，然后拿着 logFC 排前几名的基因去做 qPCR，最后一无所获——问题往往出在检验方法、背景基因集和表达比例过滤这几处。

2.1 三种检验方法的选择逻辑与耗时对比

method参数有四个可选值：'t-test'、't-test_overestim_var'、'wilcoxon'、'logreg'。默认是't-test'，这一点在文档里写得很小，容易被忽略。

这四种方法不是"谁更先进"，而是针对不同数据假设的不同工具。t-test假设两组方差相同，单细胞数据里一个群 500 个细胞、另一个群 5000 个细胞，方差天然不等，所以官方提供了t-test_overestim_var，它对不等方差做过度估计校正，在群大小悬殊时更保守，也更容易给出可信的显著性。wilcoxon是非参数秩和检验，不假设正态分布，对单细胞这种零膨胀、偏态分布的表达矩阵更友好，代价是慢——两三万个细胞跑起来，耗时可能是 t 检验的十几倍。

logreg是逻辑回归，返回的是每个基因作为分类特征的权重，没有 p 值。它适合的场景是"我想知道哪些基因组合能把这两个群区分开"，而不是"哪个基因在两个群之间显著上调"。用法上有个细节，logreg结果的uns里没有pvals和logfoldchanges，如果你后面的代码统一用sc.get.rank_genes_groups_df取表，会因为字段缺失报错，得单独处理。

我自己的习惯是：探索阶段全部用t-test_overestim_var，速度快、够用；结论定稿前，把关键的那两三个比较换成wilcoxon复核一遍，看 top 基因列表有没有大的重排。如果重排严重，说明这些基因的显著性本身不稳，往下游交付时要打折扣。

2.2 几个决定结果质量的隐藏参数

比method更容易被忽视的是下面这几个。

第一是use_raw。它的默认值是None，含义是"如果adata.raw存在就用 raw"。听起来很贴心，但问题在于adata.raw里的基因集可能和adata.var不一致，而且如果你中途对adata.X做过sc.pp.scale，adata.X里就全是 z-score，直接拿去做 t 检验毫无意义。稳妥做法是显式写use_raw=True，并且确认adata.raw.X是 log1p 归一化后的数值。

第二是pts。默认False，一旦打开为True，会额外计算每个基因在每个群里的表达细胞比例，存在uns['rank_genes_groups']['pts']里。这个字段的最大用途是过滤——一个基因 logFC 很大但只在 3% 的细胞里有表达，大概率是噪声或极少数细胞类型混入。开启它几乎不增加耗时，强烈建议一直开着。

第三是corr_method。默认是'benjamini-hochberg'，做的是 FDR 校正。别改成'bonferroni'，在单细胞这种动辄检测两万个基因的场景下，Bonferroni 会把绝大多数基因的校正后 p 值压到 1.0，等于自废武功。

第四是tie_correct。默认False。单细胞数据里大量基因的表达值是完全相同的（尤其是低表达基因，一大堆 0），秩和检验遇到大量相同秩会算不准。设成True会启用更严格的同分处理，耗时略增但更准确。这个参数只在wilcoxon下有意义。

第五是n_genes。默认值在文档里是None（表示全部基因都参与计算），但可视化时sc.pl.rank_genes_groups有个自己的n_genes=20默认。别混淆这两个。

sc.tl.rank_genes_groups( adata, groupby='leiden', method='t-test_overestim_var', use_raw=True, pts=True, corr_method='benjamini-hochberg', tie_correct=False, key_added='rank_leiden' )

用key_added给结果起个明确的名字，是为了后面可能要对另一个分组变量再做一次差异分析时不覆盖前一次的结果。官方示例里常常省略这一步，导致第二次调用把uns['rank_genes_groups']直接顶掉，回滚都回不来。

2.3 从 uns 里把结果干净地取出来

结果存在adata.uns['rank_genes_groups']里，是一个结构化的记录数组，直接打印很难看。早年的写法是一堆手动的recarray索引，现在有更干净的接口：

result = sc.get.rank_genes_groups_df(adata, group='0', key='rank_leiden') print(result.head(20))

返回的是标准 DataFrame，列为names、scores、logfoldchanges、pvals、pvals_adj、pct_nz_group、pct_nz_reference（后两列在pts=True时才有）。拿到 DataFrame 之后，过滤条件可以自己定，比如：

sig = result[ (result['pvals_adj'] < 0.05) & (result['logfoldchanges'] > 1) & (result['pct_nz_group'] > 0.25) ].sort_values('scores', ascending=False)

这三个条件的含义分别是统计显著、效应量足够、表达广泛。缺任何一个，筛出来的基因都可能不可靠。我在一个真实项目里见过只按logfoldchanges排序、不看表达比例的情况，排第一的基因在目标群里的表达比例只有 2%，后来验证发现是双细胞（doublet）残留。

2.4 五类可视化图各自适合什么场景

scanpy 为差异分析结果配了一整套rank_genes_groups_*系列函数，功能上互相重叠但侧重不同，混用会浪费很多时间。

sc.pl.rank_genes_groups输出的是每个群一列、纵向排列的基因名加数值的正负柱状图。它最大的价值是快速扫一遍所有群的 top 基因，判断分群是否合理。参数sharey=False让各列独立显示前 N 个基因，否则会按全局统一，低分群被压扁。

sc.pl.rank_genes_groups_dotplot是最常被截图引用的那张图。横轴是基因、纵轴是群，点的大小代表表达比例、颜色代表平均表达量的标准化值。n_genes=5表示每个群取前 5 个，standard_scale='var'按基因做标准化。这张图信息密度高，但群一多就挤，超过 15 个群建议拆成两张。

sc.pl.rank_genes_groups_stacked_violin用堆叠小提琴展示分布，比点图更能看出表达的离散程度和是否双峰。缺点是横向占用空间大，一行放不下几个基因。

sc.pl.rank_genes_groups_matrixplot是色块热图，适合基因数量在 20 到 50 之间的中尺度展示。sc.pl.rank_genes_groups_heatmap则是给单个细胞粒度的完整热图，几百个基因乘几千个细胞也画得动，但除非要展示异质性，一般用不上。

2.5 我在这一步踩过的坑

第一个坑：对没有做过滤的基因集做差异分析。有些基因在整个数据集里只在一个细胞里表达，t 检验会给出一个巨大的 logFC 和一个不显著的 p 值，但如果你只看 logFC 排序，它就会排到最前面。所以在rank_genes_groups之前，先用sc.pp.filter_genes(adata, min_cells=10)把极低表达基因剔掉，这一步对结果稳定性的提升比换检验方法还明显。

第二个坑：混淆了 raw 与 X。如果中途执行过sc.pp.scale(adata, max_value=10)，adata.X就变成了带负数的 z-score。此时use_raw=True用的是 raw 里的 log 数据没问题，但如果你后续想画sc.pl.violin看表达，又忘了指定use_raw，画出来的图会有负值，坐标轴看起来很奇怪。养成"凡是看表达量就写use_raw=True"的条件反射。

第三个坑：对同一批细胞做了多轮不同分组变量的差异分析，却只保留了一份uns。key_added一定要用，命名规则建议是rank_加分组变量名。

第四个坑：p 值信仰。单细胞的 p 值受细胞数影响极大，一万个细胞和十万个细胞做同一个比较，前者不显著后者极显著，但生物学差异没变。所以判断一个 marker 好不好，效应量和表达比例永远优先于 p 值。

3. 细胞类型注释：把 cluster 编号变成能读懂的名字

聚类算法给出的是 0 到 15 的数字，注释要做的是把这些数字变成 T 细胞、B 细胞、单核细胞、树突状细胞这类能被生物学同行读懂的标签。这个过程没有一键方案，官方示例里也是半手动完成的，靠的是"差异基因 + 已知 marker 集 + 打分"三路交叉验证。

3.1 手动注释的完整流程与打分策略

我通常的做法是先把每个群的 top 20 差异基因导出来，人工扫一遍有没有熟悉的 marker。这个"人工扫"的动作看起来原始，但它是目前可靠性最高的环节，因为自动化注释工具依赖的参考数据集未必覆盖你研究的组织类型。

扫完之后，进入验证阶段：选一批公认的谱系 marker，比如 T 细胞用 CD3D/CD3E/CD2，B 细胞用 CD79A/MS4A1/CD19，NK 细胞用 NKG7/GNLY，单核细胞用 LYZ/CD14/FCGR3A，树突状细胞用 FCER1A/CST3，血小板用 PPBP。然后做一张点图：

marker_dict = { 'T': ['CD3D', 'CD3E', 'CD2'], 'B': ['CD79A', 'MS4A1', 'CD19'], 'NK': ['NKG7', 'GNLY', 'KLRD1'], 'Mono': ['LYZ', 'CD14', 'FCGR3A'], 'DC': ['FCER1A', 'CST3'], 'Platelet': ['PPBP', 'PF4'] } sc.pl.dotplot(adata, marker_dict, groupby='leiden', standard_scale='var', use_raw=True, cmap='RdBu_r')

点图给出的判断依据有三个层次：点的大小看表达比例，颜色深浅看标准化后的平均表达，两者同时满足才算是这个群的 marker。如果某个群对某组 marker 显示"点大但颜色淡"，说明表达广泛但强度低，可能是污染或低质量细胞；"点小但颜色深"则说明只有极少数细胞高表达，更像是混入的少量其他细胞类型。

3.2 score_genes 的底层逻辑与参数含义

sc.tl.score_genes干的事情不复杂：给定一组目标基因，再随机抽一组表达量分布相近的对照基因，两者平均表达量相减，得到一个打分，存在adata.obs的新列里。它的价值在于把一个群的判定从"看这个 marker 有没有点"升级成"这个群整体上是不是这套程序的转录特征更强"。

关键参数是ctrl_size和n_bins。n_bins=25的含义是把所有基因按平均表达量分成 25 个区间，然后从目标基因所在的每个区间里抽同样数量的对照基因。这么设计是为了校正表达丰度带来的系统性偏差——高表达基因天然数值大，如果对照组随便抽，打分会被丰度主导。ctrl_size=50表示每个区间抽 50 个，总共 25 乘 50 大约 1250 个对照基因。

sc.tl.score_genes( adata, gene_list=marker_dict['T'], score_name='T_score', ctrl_size=50, n_bins=25, random_state=0, use_raw=True )

这里有两点值得强调。random_state必须固定，否则每次运行抽到的对照基因不同，打分会有细微波动，导致你在不同 notebook 里看到的阈值不一样，反复怀疑自己。另外基因列表不要给太长的，几十个足够，给几百个会把信号平均掉，打分反而变得没有区分度。

3.3 细胞周期打分与生长期细胞的识别

sc.tl.score_genes_cell_cycle是一个特化的封装，内置了 S 期和 G2M 期的经典基因列表，跑完之后在obs里加三列：S_score、G2M_score、phase（取值 S、G2M 或 G1）。

sc.tl.score_genes_cell_cycle(adata, s_genes=s_genes, g2m_genes=g2m_genes) sc.pl.umap(adata, color=['S_score', 'G2M_score', 'phase'], ncols=3)

为什么要单独做这个？因为增殖中的细胞转录组会被细胞周期强烈扭曲，同一个细胞类型在 G1 和 G2M 期可能被分到不同 cluster。做分化轨迹分析时如果不处理，伪时间的第一主成分很可能捕捉到的是细胞周期而不是分化方向，结论就全错了。

处理方式有两种思路。保守做法是只把phase作为观察变量叠到 UMAP 上看分布，如果发现某个 cluster 里 S/G2M 细胞高度富集，注释时就要标注"增殖型 XX 细胞"，而不是把它当作独立亚群。激进做法是用sc.pp.regress_out(adata, ['S_score', 'G2M_score'])把周期信号回归掉，但这会引入人工痕迹，也不适合用在稀疏矩阵上，现在更多人倾向于不回归、只在解读时注意。

注意：打分用的基因列表必须是物种对应的。拿人的基因符号去打小鼠数据，score_genes不会报错，只会返回一堆接近 0 的分数，然后你会以为这批细胞没有增殖。

3.4 注释结果回写与可复现性

注释的最终产物是新的一列celltype。别直接覆盖leiden，两者要并存：

cluster2celltype = { '0': 'CD4 T', '1': 'CD14 Mono', '2': 'B', '3': 'CD8 T', '4': 'NK', '5': 'FCGR3A Mono', '6': 'DC', '7': 'Platelet' } adata.obs['celltype'] = adata.obs['leiden'].map(cluster2celltype).astype('category')

.astype('category')这一步别省。scanpy 里所有按分类变量上色的操作都依赖 pandas 的 category dtype，如果是普通 object 列，sc.pl.umap会在图例顺序、配色分配上出各种奇怪问题，比如类别按字符串排序而不是按数值排序，0, 1, 10, 11, 2这种。

更关键的是把这个映射字典也存下来。我习惯写进adata.uns['cluster2celltype']，跟着 h5ad 一起保存：

adata.uns['cluster2celltype'] = cluster2celltype adata.write_h5ad('pbmc_annotated.h5ad', compression='gzip')

这样即使三个月后别人拿到文件，也能追溯每个数字是怎么变成名字的。评审问起来"你的 CD4 T 是怎么定的"，能拿出映射关系和当初的点图，比口头解释有说服力得多。

3.5 注释阶段的常见误判

误判一：把"表达某 marker"等同于"是某细胞"。单核细胞和树突状细胞共享大量基因，巨噬细胞和单核细胞在部分组织里几乎无法用少数几个 marker 区分。遇到边界模糊的群，宁可标成 "Mono/DC" 这种复合标签，也不要硬给一个精确名字。

误判二：忽略组织背景。同样一套 marker，在外周血里分得清清楚楚，到了肿瘤组织里因为应激反应、代谢重编程，表达谱会大幅漂移。判断标准要以你自己的数据分布为准，marker 集只是参考。

误判三：分辨率过高导致过度拆分。一个明确是 B 细胞的群被拆成三个 cluster，如果三个群的差异基因全是核糖体基因、线粒体基因或热休克蛋白，那说明是技术因素而不是生物学亚型。这时候该做的是降低resolution重新聚类，而不是费劲给三个群分别命名。

4. 轨迹推断：PAGA 与扩散伪时间的完整实操

UMAP 图上的距离没有严格含义，这是常识。但很多人看着 UMAP 上两个群靠得近，就认为它们之间有分化关系，这个推断在单细胞领域被反复证伪。轨迹推断要解决的就是这件事：用图结构和随机游走的方式，给出有数学依据的细胞状态转移方向。

4.1 为什么聚类图看不出分化顺序

聚类的本质是把相似的细胞归到一起，它天然抹掉了群内部的连续变化。一个正在分化的细胞群，在聚类结果里可能被切成三个阶段，也可能因为过渡细胞太少而被并入相邻群。更麻烦的是，聚类对"方向"这个概念完全无感，它只知道 A 和 B 相似，不知道是 A 变成 B 还是 B 变成 A。

PAGA 的思路不同。它在已有的 cluster 划分基础上，计算群与群之间的连通性——本质上是问"如果随机选一个 A 群的细胞，它的近邻里有多大比例是 B 群的细胞"。这个比例除以随机期望，就得到连接强度。连接强度高说明两个群之间存在大量中间态细胞，很可能是连续过渡的；连接强度低说明是两个独立的群体。

扩散伪时间（DPT）走的是另一条路：在细胞级别的邻居图上做扩散，计算每个细胞到达指定根细胞需要多少步随机游走，得到一个连续的时间值。它的前提是假设数据存在一个分支结构，而且你得指定起点。

4.2 扩散图与 PAGA 连通性图的分工

两者的输出形式完全不同，用途也不同。扩散图（sc.tl.diffmap）是细胞级别的降维结果，坐标轴代表扩散成分，前几个成分携带了主要的拓扑结构信息。PAGA 图是群级别的抽象图，节点是 cluster，边是连通强度。

sc.tl.diffmap(adata, n_comps=15) sc.pl.diffmap(adata, color='leiden')

n_comps默认 15，意思是计算 15 个扩散成分。这个数字看起来随意，其实有讲究：扩散成分按携带的结构信息量递减排列，前 3 到 5 个通常就够看，但 DPT 计算如果只保留太少成分，会把分支结构压平。默认 15 是个经验值，一般不用改。要注意的是diffmap依赖adata.uns['neighbors']里的连通性矩阵，如果之前为了省内存把它清了，这里会报ValueError。

PAGA 这边：

sc.tl.paga(adata, groups='leiden') sc.pl.paga(adata, plot=False, threshold=0.03)

threshold控制边的显示下限，低于这个值的连接不画。默认在sc.pl.paga里是None（画所有边），实际用的时候几乎总得设一个，否则图上会出现一团浆糊。0.01 到 0.05 之间是个合理的起点，具体看你的群数量和连接密集程度。threshold只影响可视化，不影响tl.paga计算出来的完整连通性矩阵，那个存在adata.uns['paga']['connectivities']里。

4.3 PAGA 初始化 UMAP 的实际效果

这一步是官方示例里容易被当作可选步骤、实际却很有用的一个技巧。普通的sc.tl.umap从随机位置开始优化，结果的全局拓扑结构（哪个群在哪个方向）每次运行都可能不同。用 PAGA 的结果作为初始位置，UMAP 会在保持局部邻域关系的同时，尽量复现群级别的连接结构：

sc.tl.umap(adata, init_pos='paga', random_state=0)

效果上，原来散落在 UMAP 边缘、被切成几块的同一个谱系，会更倾向于聚到一条连续的带上。代价是 UMAP 的局部精细结构会略有损失，如果你主要关心亚群聚类而不是全局走向，可以不用这一步。我通常的做法是两种都跑一遍存成两个 embedding，比较之后选更符合生物学预期的那张，同时把选择理由记录在文档里——这一点很重要，因为 UMAP 布局本身带有随机性和主观性，选图的依据必须写下来。

4.4 根节点选择与 DPT 伪时间计算

DPT 的起点必须手动指定，这是整个流程里最主观的一步，也是最容易出错的一步。根细胞的选法通常是"选一群公认处于最早状态的细胞"，比如造血数据里选造血干/祖细胞群，肠上皮数据里选隐窝底部的细胞群。

import numpy as np adata.uns['iroot'] = np.flatnonzero(adata.obs['leiden'] == '0')[0] sc.tl.dpt(adata) sc.pl.umap(adata, color=['dpt_pseudotime', 'leiden'])

np.flatnonzero返回的是满足条件的细胞在obs里的位置索引，取第一个即可。这里有个陷阱：如果leiden是 category 类型，adata.obs['leiden'] == '0'得到的布尔序列是没问题的；但如果你的群名已经被改成了 'HSC' 这类字符串，就得对应改条件。

根细胞选在群边缘还是群中心，会影响伪时间的起点值（边缘细胞可能接近 0，中心细胞是一个小的正数），但不影响整体的相对顺序。所以更稳妥的做法是多选几个根细胞取平均，或者用n_branchings参数让算法自己找分支点：

sc.tl.dpt(adata, n_branchings=1, allow_kendall_tau_shift=True)

4.5 结果解读与常见错误

解读 DPT 结果有三条纪律。第一，伪时间值本身没有绝对单位，"0.8" 不代表任何具体时间，它只表示在随机游走意义下离根细胞更远。第二，只有落在同一条分支上的细胞之间才能比较伪时间，分叉之后的两条分支共享一段共有路径，过了分叉点的数值不能直接对比早晚。第三，伪时间的方向完全由根节点决定，选错根节点会得到一条看起来合理但完全反向的轨迹。

常见错误一：KeyError: "Could not find key 'iroot'"。原因是没设根细胞或者设置的位置不对，uns里没有这个键。

常见错误二：伪时间取值呈双峰。这通常意味着你的数据里有两个互不连通的群，随机游走需要绕很远才能到达，导致中间出现断崖。可以先检查 PAGA 连通图是不是分裂成两块，如果是，说明这两部分不该放在一起做轨迹，得先确认它们是否属于同一谱系。

常见错误三：伪时间几乎全是 0 或全是 1。检查一下根细胞所在的群是不是占了绝大多数细胞，或者diffmap的n_comps设得太小导致结构被压平。反过来，如果根细胞群特别小，伪时间会迅速饱和，也会看起来没区分度。

5. 绘图系统拆解：从 sc.pl 到出版级图

scanpy 的绘图函数有个统一的设计哲学：底层是 matplotlib，但把常用参数都预先调好了，所以能一行出图；同时又留了return_fig和ax两个口子，让需要精修的人能接着往下改。理解这个双层结构，比记住几十个参数名有用得多。

5.1 sc.pl 的返回值机制与二次定制

原始的sc.pl.umap直接调用的行为是画完就结束，save参数会把图存到sc.settings.figdir指定的目录。想在图的基础上继续加东西，必须用return_fig=True：

fig = sc.pl.umap(adata, color='celltype', return_fig=True) ax = fig.axes[0] ax.set_title('') ax.spines['top'].set_visible(False) fig.set_size_inches(6, 5) fig.savefig('umap_celltype.pdf', bbox_inches='tight', dpi=300)

注意fig.axes[0]取的是第一个子图。当color传入多个变量时，会生成多个子图，axes里有多个元素，得逐个处理。另外bbox_inches='tight'几乎是必加的，否则保存的时候图例容易被裁掉一半，这是投稿阶段最常见的格式问题。

对于 dotplot、matrixplot、stacked_violin 这几个返回绘图对象的函数，定制的接口不同。它们返回的是DotPlot、MatrixPlot、StackedViolin对象，有各自的方法链：

dp = sc.pl.dotplot(adata, marker_dict, groupby='celltype', use_raw=True, return_fig=True) dp.style(cmap='RdBu_r', edge_color='gray').legend() dp.make_figure() dp.savefig('dotplot.pdf', bbox_inches='tight')

.style()里能改颜色映射、边框、字体大小，.legend()用来调整图例位置。这套接口是 scanpy 里设计得相对现代的，值得花时间熟悉。

5.2 高维散点图的参数调优清单

sc.pl.umap和sc.pl.diffmap共享一套散点图参数，其中真正影响可读性的就那么几个。

size控制点半径，默认是根据细胞数自动算的。细胞数上万时默认值往往偏小，一片糊；一两千细胞时默认值又偏大，点叠在一起。手动给size=20到size=80之间试试，具体看画布尺寸。

alpha是透明度，默认 1.0（部分版本）。细胞高度重叠的区域，降低到 0.6 左右能看出密度层次，但低于 0.3 就太淡了。

frameon控制是否画坐标轴框线。发表用的 UMAP 一般不要边框，因为 UMAP 坐标数值本身没有意义，画刻度反而误导读者。设frameon=False。

legend_loc决定图例位置，'on data'把标签直接打在群中心，'right margin'放到右侧，'none'不显示。群多的时候'on data'更清爽，但要配legend_fontoutline=2，否则标签在深色点上完全看不清。

palette指定配色。群数超过 10 个的时候，默认配色循环会开始重复，导致两个不同的群颜色一样。这时候要么换成sc.pl.palettes.default_102，要么干脆换成连续色阶加上标号，避免混淆。

groups和groupby的组合很实用：groupby='celltype', groups=['T', 'NK']只画这两个群，其他群的细胞不显示。做亚群展示时比 subset 数据再画快得多。

vmin和vmax用于基因表达上色。默认取该基因的最小最大值，遇到一个极端高表达的细胞就会把整张图压成一片浅色。设vmax为 95 分位数是常用做法：

import numpy as np expr = adata[:, 'CD3D'].X.toarray().flatten() sc.pl.umap(adata, color='CD3D', use_raw=True, vmax=np.percentile(expr, 99))

5.3 导出、拼图与尺寸控制

sc.settings.figdir设好之后，sc.pl.umap(..., save='_xxx.pdf')会自动存成umap_xxx.pdf。文件名前缀由函数名决定，后缀由save参数决定，格式由扩展名决定。这个机制很方便，但有个坑：它不处理bbox_inches，复杂图例被裁的问题依旧存在，所以重要图我还是走return_fig加手动savefig。

拼图方面，sc.pl.umap的color参数传列表就能生成多子图，ncols控制列数，wspace控制子图间距。三个以上的基因加上几个分类变量，用这个方式一次出图比手动拼快得多。真要精细排版，还是导出各个子图后在排版软件里拼，fig.set_size_inches单独控制每张图的宽高，保证字体大小一致。

5.4 中文字体与坐标轴的几个细节

如果图里要出现中文标签（比如群名用中文），必须先配置 matplotlib 的字体，否则全是一排方块：

import matplotlib matplotlib.rcParams['font.sans-serif'] = ['Noto Sans CJK SC'] # 换成你本机已安装的字体 matplotlib.rcParams['axes.unicode_minus'] = False

字体名必须与本机实际安装的完全一致，写错了不报错，只是静默回退到默认字体然后继续显示方块。axes.unicode_minus那一行是为了让负号正常显示，很多中文字体自带的负号字形有问题。

另一类细节是小提琴图和点图的横向空间。sc.pl.violin的rotation=45让基因名倾斜，避免横向重叠；sc.pl.dotplot的基因名如果太长，var_names可以传一个字典做重命名：

sc.pl.dotplot(adata, {'T cell': ['CD3D', 'CD3E']}, groupby='celltype')

这个用法文档里提了一句但很容易漏，实际上非常实用——既能分组展示，又能把基因名换成好看的显示名。

6. 报错与异常排查速查表

单细胞分析里真正的挫败感不来自算法难，而来自那些"代码看着没错但结果就是不对"的时刻。下面这张表是我这几年攒下来的高频问题，按报错信息、根本原因、处理办法三列整理。

6.1 高频报错对照表

6.2 结果看起来不对劲的排查思路

有一类问题不报错，但结果明显不合逻辑，这时候排查顺序比排查工具重要。

第一步永远是核对数据来源。我把adata.raw的 shape、adata.X的数值范围打印出来看，log1p 之后最大值通常不超过 10，如果有两位数说明没归一化；如果出现负数说明被 scale 过。这一步能解决三成以上的"结果不对"。

第二步看分群数量是否合理。外周血单个核细胞的数据，两三千细胞分出 8 到 12 个群是常态，如果只分出 3 个，说明resolution太低或者 PCA 成分太少；如果分出 30 个，没有任何一群超过 10 个细胞，那就是过拟合，把噪声当成了结构。

第三步把关键基因的表达叠到 UMAP 上看。如果一个据说在 T 细胞里高表达的基因在 UMAP 上呈现随机分布，说明要么基因名拼错了（大小写、别名问题），要么use_raw指向的矩阵不对。基因名大小写是个隐蔽的坑，CD3D和Cd3d在不同物种规范里都有，但对不上就是取不到。

第四步做一次独立验证。挑几个群的 marker，用另一套数据（同项目的另一个样本）重复差异分析，看 top 基因是否重合。重合度低说明前面的结果不稳定，可能是批次效应或者细胞数太少。

7. 管线化与性能：让分析能重复跑第二遍

一份分析的价值，一半在于结论，另一半在于别人能不能复现。我第一次做的单细胞分析，三个月后自己想重跑一遍，发现中间隔了几十个手工修改的单元格，改了什么完全记不清。后来我把这些习惯固定下来，重跑时间从一整天缩短到一个小时。

7.1 内存与稀疏矩阵的处理习惯

单细胞矩阵稀疏度通常在 90% 以上，绝大多数元素是 0。一旦某一步把它转成稠密数组，内存占用会瞬间涨十倍以上，几万细胞的数据直接吃满三十二 G。所以矩阵格式必须一直保持 CSR 稀疏：

from scipy.sparse import issparse print(type(adata.X), issparse(adata.X))

如果发现变成了numpy.ndarray，通常是被sc.pp.scale或者某个手写的 numpy 操作转过去的。sc.pp.scale有个zero_center参数，默认True会稠密化；设成False只做缩放不做中心化，能保持稀疏。真要中心化的时候，只对高变基因做，别对全基因做。

另一个习惯是及时的切片顺序。sc.pp.highly_variable_genes之后，先adata.raw = adata，再adata = adata[:, adata.var.highly_variable]，把全基因信息留在 raw 里，主矩阵瘦身。反过来先切片再设 raw，raw 里也就只剩高变基因了，后面做 marker 分析会漏掉一大片。

adata.write_h5ad('xxx.h5ad', compression='gzip')这个压缩参数值得加，文件体积通常能减掉一半以上，读取速度损失可以接受。几个 G 的文件反复读写，日积月累省的时间很可观。

7.2 中间结果的保存与版本管理

我现在的规矩是每个大阶段结束存一个文件：01_raw.h5ad、02_qc.h5ad、03_clustered.h5ad、04_annotated.h5ad。文件名带序号，一眼能看出顺序。每个文件存完立刻用一次干净的内核重新读进来验证一遍，确保没有存坏。

参数记录不要靠脑子。我的做法是在adata.uns['params']里存一个字典，记录关键参数的取值：

adata.uns['params'] = { 'n_top_genes': 2000, 'n_pcs': 30, 'n_neighbors': 10, 'resolution': 0.8, 'seed': 0, 'qc_min_genes': 200, 'qc_max_pct_mt': 20 }

这个字典跟着 h5ad 走，任何人打开文件都能看到当时是怎么设的。比在 notebook 里写注释可靠得多，注释会被删，uns不会。

版本管理上，scanpy、anndata、numpy、scipy四个包的版本要一起记。anndata在 0.8 到 0.10 之间对uns里嵌套结构的支持改过，跨版本读写 h5ad 偶尔会出兼容问题。最简单的办法是在项目根目录放一个requirements.txt，锁死版本号。

7.3 参数记录与随机种子

随机种子这件事，说小很小，说大能让人怀疑人生。sc.tl.umap、sc.tl.leiden、sc.tl.score_genes都有随机成分，不固定种子的话，同一份数据跑两次，UMAP 的旋转角度可能完全不同，score_genes的打分值会有百分之几的浮动。如果你的分群边界本来就模糊，浮动可能直接把几个细胞换到另一个群去。

固定方式有两种。全局设sc.settings.seed = 0，会影响所有接了种子的函数；或者在每个函数调用里显式写random_state=0。我倾向于后者，因为它是局部的、可读的，看代码就知道哪一步用了什么种子。特别是leiden，它的random_state直接影响分群结果，同一份数据不同种子可能得到不同群数，这个必须固定并记录。

n_jobs参数在pp.pca、pp.neighbors里都能设，默认 1。设成 4 在多数场景下能快两三倍，但设得比机器核数还高反而变慢，因为线程调度开销。调的时候实测一次耗时，别照抄别人的值。

还有个小技巧：跑耗时的步骤之前，先用adata[:2000, :].copy()抽一个子集把整条路走通，确认没有语法错误和结构错误，再上全量数据。子集上跑一遍通常几十秒，全量可能是几十分钟，这个时间差足够你发现三个 bug。

说到最后，分享一个我自己反复验证过的体会：单细胞分析里最耗时间的从来不是写代码，而是判断"这个结果能不能信"。工具本身给的是一堆数字和图形，把它们变成结论的那一步，靠的是对数据来源的熟悉、对生物学背景的了解，以及对异常现象的敏感度。我见过太多分析卡在注释那一步，不是因为不会用dotplot，而是因为不知道该拿哪些 marker 去问。

所以我的建议是，每跑完一个模块，都强迫自己写下三句话：这一步的输入是什么状态、我做了什么假设、输出结果如果和预期不符最可能是什么原因。这三句话写下来，下次重跑或者交接给别人的时候，价值远超一段漂亮的代码。