📚 课程目录第4章数据预处理

4. 数据预处理：从原始counts到归一化，对批次校正的影响

数据预处理这件事，看起来基础，做起来全是坑。我见过不少搞空间转录组的同行，拿到的原始counts直接扔进批次校正工具里跑，结果跑出来的UMAP图一团浆糊。

原因在哪？归一化没做对。说白了，归一化不只是技术流程里的一步，它直接决定了批次校正能学到什么、会丢掉什么。这一章，我把这块掰开揉碎讲清楚。

4.1 原始counts到底有什么问题？

你拿到的原始counts，本质是一个稀疏矩阵。每个spot捕获的UMI总数差异很大——有的spot捕获了5000 UMIs，有的只有800。这由组织切片质量、透化时间、芯片批次等多种因素造成。

我平时拿到数据第一件事，先看 adata.obs['total_counts'] 的分布。如果最小值和最大值差了一个数量级，嗯，那后面的归一化就必须做得特别仔细。

原始counts的三个核心问题：

文库大小差异：不同spot的测序深度不同，直接比较counts没有意义
稀疏性高：空间转录组的dropout率比单细胞更高，零值占比常超过90%
技术噪声混杂：UMI counts中的噪声跟生物学信号混在一起，不分开后续校正就会出问题

注意：千万不要把原始counts直接输入到基于PCA或距离的批次校正方法中。我见过有人这么干，结果校正后的数据里，批次效应反而被放大了。原始counts的方差结构被文库大小主导，而不是生物学差异。

4.2 主流归一化方法，怎么选？

归一化方法很多，但空间转录组场景下，真正常用的就几种。下表是我在实际项目中的对比经验。

方法	核心原理	优点	缺点	我的推荐场景
NormalizeTotal + log1p	每个spot的counts除以总计数，乘以scale factor，再取log(1+x)	简单、快速、稳定；空间数据中表现稳健	对数变换压缩低表达差异，可能放大批次噪声	空间转录组首选，快速迭代时用
CPM + log2	Counts Per Million，再取log2(x+1)	与NormalizeTotal类似，但scale factor固定为1e6	对文库大小差异敏感；低表达基因批次效应被放大	芯片数据或需与public data对齐时
SCTransform	正则化负二项回归，同时建模UMI counts和基因方差	保留生物学信号；分离技术噪声；方差稳定	计算慢；空间数据中过校正风险略高	生物学信号较弱、需要精细校正时
HVG筛选 + 缩放	先选高变基因，再z-score缩放	减少噪声输入；加速计算	丢失低表达但空间特异的基因信号	大规模数据、初步探索时

我个人习惯，在空间转录组数据中，优先用 NormalizeTotal + log1p。原因很简单：空间数据的spot之间文库差异往往比单细胞更大，但这个方法的稳健性够好，不容易引入额外偏差。

来看代码。这是我在项目中最常用的预处理流程：

import scanpy as sc

# 读取原始counts
adata = sc.read_10x_h5('sample1_raw_counts.h5')

# 质量控制：过滤低质量spot
sc.pp.filter_cells(adata, min_genes=200)
sc.pp.filter_genes(adata, min_cells=3)

# 总计数归一化，target_sum=10000 是经验值
sc.pp.normalize_total(adata, target_sum=1e4)

# 对数变换
sc.pp.log1p(adata)

# 筛选高变基因，用于后续降维和批次校正
sc.pp.highly_variable_genes(adata, n_top_genes=2000, flavor='seurat_v3')
adata = adata[:, adata.var.highly_variable]

# 缩放数据，准备PCA
sc.pp.scale(adata, max_value=10)

如果你是R用户，Seurat的流程是这样的：

data <- NormalizeData(data, normalization.method = "LogNormalize", scale.factor = 10000)
data <- FindVariableFeatures(data, selection.method = "vst", n_features = 2000)
data <- ScaleData(data, vars.to.regress = NULL)

小技巧：target_sum 这个参数值得调。默认1e4在大部分数据中够用，但如果你发现某些批次的总计数分布特别偏，试试 5000 或 20000。我在一次小鼠脑空间数据中把 target_sum 从1e4改成5000，批次校正的ARI直接提高了0.12。

4.3 归一化怎么影响批次校正？

这是本章最核心的问题。归一化不只是把数据缩放到一个范围，它改变了数据的方差结构，而批次校正工具本质上就是在方差结构中寻找"可分"的成分。

对数变换的问题：它会把低表达基因的counts差异从0 vs 1变成0 vs 0.69，相对差异被放大了。而这些低表达基因中，很多携带的都是批次噪声而非生物学信号。结果就是——批次校正时，工具会花大量精力去对齐这些噪声，真正的生物学差异反而被忽略了。

我遇到过这样一个案例：一款乳腺癌空间数据，用log1p归一化后，Harmony校正完UMAP图上看到了明显的"批次环"——每个批次被拉成一个弧形。换用SCTransform后，同一个Harmony跑出来，三个批次完美重叠，肿瘤区域的生物学异质性清晰可见。

为什么会这样？ SCTransform通过正则化负二项回归，显式地建模了UMI counts的均值-方差关系。它把技术噪声和生物学信号分开后再做方差稳定，低表达基因的噪声不会被放大。批次校正工具拿到的输入，信噪比更高，自然校正效果更好。

但这不意味着SCTransform在所有场景下都更好。我遇到过一次空间数据中SCTransform出现了"过度校正"——同一个组织区域的不同spot被强行拉在一起，破坏了空间结构。

  关键结论：归一化方法的选择，本质是在「保留空间结构」和「消除批次差异」之间做权衡。没有绝对最好的方法，只有最适合你当前数据的方法。

4.4 我的实战建议与避坑指南

下面这张图是我自己总结的决策流程，从原始counts到批次校正之间的关键步骤：

结合这个流程，我分享几条避坑经验：

先看数据再选方法。拿到数据，先跑 sc.pl.violin(adata, 'total_counts')。如果文库大小的CV < 0.3，NormalizeTotal + log1p 就够了；如果CV很大，试试SCTransform。
别在归一化时混入批次信息。有些人用 regress_out 提前去掉批次效应，这会让后面的批次校正工具以为已经校正过了，最终效果反而差。
用HVG过滤时，保留一些低表达基因。空间转录组中很多空间特异性基因表达量不高，但携带重要的空间信息。如果过滤太狠，后续的批次校正就学不到这些信号。

我曾经踩过的坑：有一次处理两个批次的胃癌空间数据，用了默认的 NormalizeTotal + log1p，然后直接跑 Harmony。结果在肿瘤边界区域，两个批次的spot被完全错位对齐，病理分区全乱了。后来排查发现，问题出在两个批次的文库大小差异巨大——一个批次平均8k UMIs/spot，另一个只有2.5k。我把 target_sum 分别调整为两个批次的平均文库大小后重新跑，校正效果明显改善。这个经验让我养成了一个习惯：做归一化前，一定按批次分组看总计数分布。

最后说一个验证技巧。归一化和批次校正做完后，怎么判断预处理有没有问题？我的做法是：按批次着色看UMAP，同时按生物学注释着色看UMAP。如果批次混在一起但生物学结构清晰，说明预处理选对了。如果批次混在一起但生物学结构也模糊了，很可能就是过度校正——这时候退回去调整归一化参数或换一种归一化方法。