数据预处理：标准化、对数转换、高度可变基因筛选
拿到单细胞数据，第一步就是清洗。原始 AnnData 里的 counts 矩阵，说白了就是每个细胞里检测到多少条 mRNA 分子。但不同细胞的总 UMI 数差异很大，直接拿原始值去做聚类或特征选择，结果会非常难看。
举个例子：某个细胞测序深度特别高，所有基因的计数都偏高，但它并不代表这些基因真的表达更多。所以我们必须先做标准化，然后对数转换，再筛出真正有信息量的基因。这套流程在我做过的十几个单细胞项目里几乎是固定的——没有一次能跳过。
标准化 Normalization
标准化的目标很简单：让每个细胞的表达量具有可比性。我常用的方法是 CPM（Counts Per Million），也就是把每个基因的计数除以该细胞的总 UMI 数，再乘以一百万。
实际上，scran 包或 Scanpy 里更推荐 size factor 标准化。不过我个人习惯直接用 sc.pp.normalize_total，设定 target_sum 为 1e4，这其实等价于 CPM 再乘 10。
import scanpy as sc
# 读取原始 AnnData（假设已经包含 counts 矩阵在 adata.X 中）
adata = sc.read_h5ad('raw_counts.h5ad')
# 标准化：每个细胞的总表达量标准化到 10,000
sc.pp.normalize_total(adata, target_sum=1e4)对数转换 Log Transformation
标准化之后的数据仍然是右偏的——大多数基因表达值在 0 附近，少数基因很高。对数转换可以把数据拉近正态分布，对 PCA、t-SNE 这类假设方差齐性的算法更友好。
标准做法是 log1p，也就是 ln(1 + x)。在 Scanpy 里就是一行代码：
sc.pp.log1p(adata)为什么是 log1p 而不是 log2 或 log10？我个人理解是：log1p 在零值处不产生负无穷，且转换后的方差相对稳定。实际项目中我也试过 log2，效果差别不大，但社区习惯用 log1p，那就跟着主流走，避免后续工具兼容问题。
adata.raw 做两次 log1p！我曾经犯过这个错——第一次标准化+log1p 存到了 adata.X，然后又对 adata.raw 再 log1p 一次，导致后面 HVG 筛选全乱套。建议先在 adata.raw 里保留原始整数计数，只对 adata.X 做标准化和对数转换。
高度可变基因筛选 Highly Variable Genes
对数转换之后，表达矩阵里还是有两万个基因。实际上大部分基因在所有细胞里表达量都很低，或者完全没变化。我们需要挑出那些在细胞间差异最大的基因——也就是 高度可变基因 (HVG)。后面的降维、聚类全指望这批基因。
Scanpy 提供了两种主流方法：seurat_v3 和 seurat（基于方差均值关系，和 Seurat 包一致）。我建议用 flavor='seurat_v3'，它直接对原始计数操作，不依赖之前的标准化步骤，效果更稳健。
# HVG 筛选（推荐版本）
sc.pp.highly_variable_genes(adata, flavor='seurat_v3', n_top_genes=2000)
# 或者第二种 flavor（基于方差）
sc.pp.highly_variable_genes(adata, flavor='seurat', n_top_genes=2000)
# 查看筛选标志
adata.var['highly_variable'].head()筛选之后，通常我们会把数据限定到这批 HVG 上，后续分析量一下从两万降到两千，运算速度提升一个数量级。
# 只保留 HVG 用于后续分析
adata = adata[:, adata.var['highly_variable']].copy()
# 检查新形状
print(adata.shape) # (n_cells, ~2000)我一般取 2000～3000 个。取太少可能丢失稀有亚群特征，取太多易受测序噪声干扰。如果细胞数很少（< 500），可以试 1000 个；如果细胞数很多（> 10 万），3000 到 5000 也可以，但注意内存占用。
你可能会问：为什么不像传统转录组那样用差异表达筛选？其实 HVG 筛选的本质是根据方差/均值关系，而不是组间差异。它不需要先验的细胞分群信息，属于无监督的预处理手段。我在处理肿瘤异质性数据时发现，HVG 筛选后的基因集能很好地保留稀有亚群的特征基因，而传统方差过滤会漏掉那些低表达但变化大的基因。
三种主流标准化方法的对比
为了让你在实际项目里选得更顺手，我整理了一个表格：
| 方法 | 适用场景 | 核心逻辑 | 个人经验 |
|---|---|---|---|
| log1p + CPM | 常规单细胞，快速原型 | CPM 后 ln(1+x) | 最通用，Scanpy 默认流程 |
| SCTransform | 测序深度差异大时 | 负二项回归，同时估计方差 | 效果好但慢，适合小数据集 |
| NormalizeData 法 | 与 Seurat 标准 pipeline 对接 | log1p 前除以 size factor | R 用户更熟悉，Python 里用 scanpy 做类似操作 |
我个人更喜欢 log1p + CPM，因为它简单、可重复、社区接受度高。除非你发现批次效应特别严重，可以考虑 SCTransform 或更复杂的回归标准化。
整体流程：一张图讲清
下面这张 SVG 图总结了本章的核心步骤——从原始计数到最终可用的 AnnData。蓝色节点是“做”，绿色节点是“得到的结果”。
常见问题 QA
- 问：log1p 后还有零值怎么办？ 保留零值。单细胞数据本就稀疏，零值代表了“未检测到”。直接保留即可，后续 PCA 会用，有些模型专门处理零膨胀。
- 问：HVG 筛选会影响批次校正吗？ 一般先做批次校正在做 HVG。批次校正后方差结构变化，先筛 HVG 可能会丢失某些批次相关的真实变异。我曾经在整合多个数据集时先筛了 HVG，结果某些稀有细胞类型被过滤掉了，后来调整为“先除批次效应再选 HVG”。
- 问：所有基因都要做 log1p 吗？ 通常是的。如果后面用
sc.tl.rank_genes_groups做差异表达，它默认输入的是对数转换后的数据，算出来的 fold change 才是正确的。
标准化+log1p 后，打印
adata.X.max()，应该在 4～8 之间（如果 CPM 后 log1p，最大约 ln(1+10000)=9.2）。如果看到几十几百的值，说明还没做 log1p。
嗯，讲到这里，数据预处理的基础就完成了。标准化、对数转换、HVG 筛选，这三步奠定了后续一切分析的质量。我见过不少新手一上来就直接做聚类，结果全聚在测序深度上——那都是血泪教训。
最后提一句：保存中间结果。每做完一步就写一个 .h5ad，万一后面发现偏差还能回退。我自己的习惯是 adata.write('norm_log.h5ad') 和 adata.write('hvg.h5ad')。千万别嫌麻烦。