📚 课程目录第2章 RNA velocity的基本原理

2. RNA velocity的基本原理：未剪接与剪接mRNA

好了，这一章咱们来啃核心。RNA velocity 到底在算啥？说白了就是：看一个基因的 mRNA 有多少还没剪接、多少已经剪接好了。这个比例变化，就反映出基因是刚启动转录、正在活跃工作、还是已经准备收工。

你在其它地方看过很多花哨的示意图对吧？我刚开始接触这个概念时，也被那些漂亮的微分方程搞得一头雾水。后来发现，去掉数学外衣，底层逻辑其实很直观。

真核生物的基因转录出来后，首先是 pre-mRNA（未剪接） 形态。这个版本包含内含子和外显子。细胞接着会剪掉内含子、连接外显子，产生 mature mRNA（已剪接）。

在单细胞测序里，我们用 内含子 reads 表示未剪接分子，用 外显子 reads 表示已剪接分子。注意，这里有一个常见的坑——

不要直接拿“内含子区域计数”当作未剪接计数。 因为成熟 mRNA 也可能附带未剪净的内含子。我在某次分析人类胚胎数据时就遇到过这种情况，某基因的“未剪接”信号突然冲高，结果后来发现是内含子保留事件，根本不是 RNA velocity 能解释的。

所以，规范的 RNA velocity 流程会用 无内含子且跨内含子-外显子边界的 read 来定义剪接状态。这是上层算法（velocyto、scVelo）帮你做好的事，但理解了这个原理，你排查问题就快得多。

分子状态	测序信号	生物学含义
未剪接 (unspliced)	内含子区域内含外显子-内含子 junction reads	刚刚转录出来，还没加工
已剪接 (spliced)	外显子区域内不含内含子，且跨外显子-外显子 junction	成熟 mRNA，已做好翻译准备

这张表你背下来也好，但我建议你把它贴在显示器边上——每当你看到某个基因的 unspliced/spliced 比值异常时，先对照检查测序数据的定义。

RNA velocity 的核心模型其实只有两条等式。我写给你看：

du/dt = α - β·u
ds/dt = β·u - γ·s

别害怕。咱们拆开说：

什么意思？u 积累，s 也会积累。 但 s 多了一条降解路径，所以 s 的变化要比 u 慢半拍。这“半拍”就是 velocity 信号的来源。

举个例子，某个基因突然被激活（α 变大），u 立刻上升。但由于剪接需要时间，s 要过一会儿才能跟上。这中间就产生了 u 高 s 低 → u 降低后 s 升高 的轨迹。说白了，就是 u 和 s 的相位差。

个人经验：我在分析小鼠大脑发育数据时，发现一个叫 Neurod1 的基因在神经干细胞阶段 u/s 比值异常高。这就是典型的“刚被激活、正在转录”状态。后续的拟时序分析也证实，这批细胞正准备分化成神经元。

上面说的是动态过程。但静态数据里，每个细胞只有一个时间点的 u 和 s 值。我们怎么算 velocity？答案是：假设在系统层面，大部分基因长时间处于稳态下。

稳态是什么？产生量 == 消耗量。当你把大量稳态细胞画出来，u 和 s 会形成一条理论曲线。偏离这条曲线的细胞，就是“正在变化”的信号。

我这么说你可能有点抽象，但用图就一目了然了。来看这个结构：

看明白了吗？u 是被转录直接驱动的，s 是被 u 的剪接驱动的，而 s 还会被降解。 这三者的时间差，就是我们要的信息。

有了 u 和 s 的单细胞数据后，每个细胞每个基因都能得到一个 residual（残差）。这个残差就是：真实观测的 s 值减去稳态曲线上对应的 s 值。

如果残差 > 0，说明 s 比预期要多，也就是 u 正在被剪接成 s，信号方向是“增加”。

如果残差 < 0，说明 s 比预期要少，也就是 s 在降解，或者 u 的产生在下降。

嗯，这里有一个技术细节：你估出来的 α、β、γ 到底准不准？我踩过这个坑。我记得在某个公共数据上，直接用默认参数跑，结果几个关键神经发育基因居然显示方向相反。后来发现是 剪接速率 β 和降解速率 γ 的初始值设置 不合适。手动调整了参数范围才恢复正常。

避坑指南：当你某个基因的 velocity 方向明显不合常理，先检查它的 u-s 散点图。如果细胞点全部挤在原点附近无变化，那就别期望 velocity 能算出什么。不是所有基因都适合跑 velocity。

每个基因算完之后，我们拿到的是一个超高维的速度向量（基因数量维）。要可视化？不行，得投影到低维空间，比如 PCA 或 UMAP 上。

常用的做法是：
1. 单细胞的低维坐标（比如 UMAP）已经存在。
2. 用每个基因的速度向量拟合出“变化趋势”。
3. 把这些趋势投影到细胞间的欧氏距离上，插值出箭头。

这个过程，说白了就是一个 局部回归 + 方向场插值。scVelo 的 scvelo.tl.velocity_embedding() 就是在干这个活。

# 这是 scVelo 底层伪逻辑概览（不是完整命令）
# 1. 用模型估算 α, β, γ
# 2. 计算每个基因的 velocity = β*u - γ*s
# 3. 将所有基因的 velocity 聚合到 PCA/UMAP 空间
# 4. 在每个细胞邻域内计算平均速度向量 → 绘制箭头

你看，核心步骤其实就这几行。剩下的全是统计学和优化细节。

补充一句：我处理过一些患者肿瘤样本，那里的细胞异质性特别大。速度快慢差异也大。如果你发现大多数箭头都指向同一个方向，先别高兴太早——可能是你选的参数把噪声当信号了。动手调调 n_neighbors 和 min_mass 试试。

任何方法都有假设。RNA velocity 也一样：

嗯，我每次讲到这里都会加一句：工具是在简化中找规律。你用它对就行了，不要过度解读局部细节。

最后，咱们来总结一下这一章的核心要点：