📚 课程目录第2章 sctransform算法原理介绍

2. sctransform算法原理：参数化负二项回归

好了，咱们正式开始啃这块硬骨头。

很多做单细胞的朋友，上来就直接跑SCTransform()。但你要问我，不懂原理就去调参，迟早要翻车。我之前帮一个团队排查批次效应，他们跑了SCTransform，结果聚类还是按样本分群。一看日志，参数没调对，正则化强度设得太高，生物信号被当噪声抹掉了。

所以，花半小时搞懂它背后的数学直觉，绝对值。

单细胞数据最让人头疼的特征是什么？过离散。

你想想看，同一个基因在不同细胞里的UMI计数，方差往往远大于均值。这在泊松分布下是解释不了的——泊松要求方差等于均值。现实数据呢？方差经常是均值的几倍甚至几十倍。

我记得第一次处理10X数据时，画了个基因表达量的分布图，拖了个长长的尾巴。当时就知道，泊松肯定没戏。

负二项分布 (NB) 自带一个色散参数，可以同时建模均值和方差。公式长这样：

Var(Y) = μ + μ² / θ

其中 θ 是色散参数。θ 越小，过离散越严重。当 θ 趋近无穷大时，负二项就退化为泊松。

    核心理解： 泊松是负二项的特例。单细胞数据几乎总是过离散的，所以我们必须用负二项。这不是学术洁癖，这是数据本身的要求。

sctransform并没有用传统的广义线性模型那套框架。它的建模思路很直接，说白了就是三步：

但问题来了：如果对两万个基因逐个跑负二项回归，计算量太大了。我当年在集群上跑过一次，等了两个小时还没出结果。

sctransform的聪明之处在于，它用了一个参数化技巧，把拟合过程大大简化了。

sctransform模型假设：对于基因 g 在细胞 c 中的UMI计数 Y_{gc}，其分布为：

Y_{gc} ~ NB(μ_{gc}, θ_g)
其中： log(μ_{gc}) = β_{0g} + β_{1g} * log(UMI_total_c)

这里的 log(UMI_total_c) 就是每个细胞的总UMI测序深度取对数。为什么要用对数？因为测序深度和表达量之间，基本上是幂律关系，取对数之后变成线性，好拟合。

然后重点来了：sctransform假设这些回归系数（β₀和β₁）与基因的均值之间，存在一个平滑的非线性关系。

换句话说，sctransform先用一小部分基因（默认2000个）拟合出完整的回归模型，然后通过核回归，预测出所有基因的回归系数。

避坑指南：我曾经在只有5000个细胞的小数据集上，强行调整了参数，让sctransform用所有基因做完整拟合。结果运行时间暴涨了5倍，但最终的聚类结果跟默认参数没差别。所以，除非你的数据非常特殊，否则默认的2000个基因子集就够用了。

我整理了一张流程图，把整个参数化建模过程串起来：

模型拟合完成后，每个基因都有了完整的参数。接下来计算Pearson残差：

z_{gc} = (y_{gc} - μ̂_{gc}) / √(μ̂_{gc} + μ̂_{gc}² / θ̂_g)

这个残差z_{gc}，就是你最终拿去做降维、聚类的输入值。

为什么不用原始的logCPM、TPM？因为那些方法没能把测序深度的影响完全剥离。sctransform的残差，去除的是"由测序深度带来的表达量变化"，保留的是"超出预期的那部分信号"。

重要警告：Pearson残差的取值范围理论上可以是正无穷或负无穷。但实践中绝大多数的值都在[-3, 3]区间内。如果你看到大量超过这个范围的值，说明你的模型拟合有问题，或者数据中存在严重的异常细胞。我排查过一个案例，是因为部分细胞的双细胞率过高，导致残差爆炸。

我把sctransform与传统标准化方法的核心差异，整理成了一张表：

说到这里，有朋友可能会问：「你花了这么大力气去拟合β₀和β₁，这两个系数到底代表什么？」

我个人习惯这么理解：

我记得在分析肾脏单细胞数据时，β₁特别大的基因大多是核糖体蛋白基因。这说明这些基因的表达量高度依赖于细胞总体转录活性，而不是细胞类型特异的信号。在做降维之前，把这类基因的残差压低，反而让真实的细胞类型标记基因的信号更突出。

    一句话总结： 参数化负二项回归不是玄学。它让每个基因的标准化过程有了一个"专属模型"。你只需要提供一个测序深度，模型就能告诉你：这个表达量，是正常的测序深度效应，还是真实的生物学信号。

嗯，原理部分就讲这么多。下一层，我们会把这点理论拆成代码，一行一行跑通。