dN/dS (ω) 是一个在进化生物学领域常出现的词，也是我最近分析的时候着重关注的一个对象。这几天翻了一些文献和书籍，对它的计算和设计理念有了一个全新的了解，所以有了这篇博客。

在阅读这篇文章前，你可能需要一定的相关知识基础，例如，你应该已经知道以下概念：

  • 密码子、简并位点、同义突变、非同义突变。
  • 中性理论、正选择、纯化选择。

dN/dS 是什么

这里的 S 即 Synonymous，N 即 Nonsynonymous。

dN 指的是非同义替换率，dS 指的是同义替换率。但我认为另一个知乎博主文章里的说法会更加直观一些：

  • dN: 每个非同义位点上发生的非同义突变的数量。
  • dS: 每个同义位点上发生的同义突变的数量。

计算原理

首先，简单来说，dN 和 dS 的计算涉及到四个变量，分别为同义突变数量（S）、非同义突变数量（N）以及同义位点数量（Ssite）和非同义位点数量（Nsite）。

前两个变量 N 和 S 非常好理解，所以这里需要重点关注的是后两个变量 Nsite 和 Ssite，即每种类型的位点数量是什么？它们是如何得到的？但在这之前，我们先简单回顾一下一些常见的突变类型。

突变类型

我们都知道，一个密码子在发生突变时，其性质可能发生以下几种变化，以密码子 TTT 为例（其编码苯丙氨酸 Phe）：

  • 错义突变：由编码一个氨基酸变为编码另一个氨基酸，例如 TTT(苯丙氨酸)>ATT(异亮氨酸)。
  • 同义突变：突变过后，仍编码相同氨基酸，例如 TTT(苯丙氨酸)>TTC(苯丙氨酸)。同义突变通常发生在密码子的第三个碱基上（有时第一个碱基也可以，例如 AGA>CGA 皆编码精氨酸），此时我们称这个位点具有简并性。需要注意的是，在某些情况下，终止密码子突变为终止密码子（例如 TAA>TAG）也被视为一种同义突变。

除了这一些可能的突变以外，还有一些其他的情况，例如，TAT(酪氨酸) 可能变为一个终止密码子（TAA），这时我们称其为无义突变，在某些情景下，这也可以被称为截断突变 (protein-truncating variant，PTV)。

以上突变类型并不涵盖全部情况，但这些突变是大多 dN/dS 算法的关注对象。

部分情况下，错义突变和无义突变会被统一视为非同义突变，当然在一些算法中，无义突变并不被考虑。这也取决于具体的研究目的和需要。

位点类型

现在我们知道，一个位点的突变可能带来的结果也不同，这也是计算 Ssite 和 Nsite 的关键所在。

首先对于任何情况下的 Nsite 和 Ssite 计算结果，我们都有：

例如，对于一个长度为 3 的密码子，它所计算出来的 Nsite 和 Ssite 总和也必定为 3。

接下来就是计算过程，最简单的计算即遍历每个位点的每种突变情况可能导致的突变类型，对于 TTT，我们可以列举出它的所有可能突变情况，这里直接用杨子恒《计算分子进化》中的一张图作为示例：

我们先仅考虑 k=1 一栏，k=2 在后文中我们会拓展讨论。

其中，对于一个密码子的每一个位点，它都可以突变成除它自身之外的碱基，因此每个位点都有三个 Target codon，其中同义突变的情况仅有 TTT>TTC，而其他都为非同义突变。

因此，这里共有 9 个情况，其中同义突变占 1 个，非同义突变占 8 个。因此，同义位点数量为 1/9 * 3（3 即密码子长度），非同义位点数量为 8/9 * 3。最后则分别为表中的 1/3 和 8/3。

突变计数及 dN/dS 计算

知道位点数量后，我们就可以计算了，依然以一个最简单的情况作为示例，假设我们有一个包含九个 TTT 的 CDS：

  • TTTTTTTTTTTTTTTTTTTTTTTTTTT

通过计算我们可以知道，它具有 3 个同义位点和 24 个非同义位点，假设这里在两个密码子上我们观察到了一个同义突变和一个非同义突变：

  • TTCTTTTTTTTTTTTTTTTTTTTTTAT

那么此时，其 dN = 1/24，dS = 1/3，计算得到 dN/dS = 1/8。

当然，在突变计数步骤中也可能出现更加复杂的情况，例如同一个密码子上出现了两个突变，例如 TTT 变为了 TCA，这种情况就要考虑多条进化途径：

  • TTT > TTA > TCA（苯丙氨酸 > 亮氨酸 > 丝氨酸），这条途径涉及到两个非同义突变。
  • TTT > TCT > TCA（苯丙氨酸 > 丝氨酸 > 丝氨酸），这条途径涉及到一个非同义突变和一个同义突变。

对于这种情况，常见的做法有对不同的途径赋予合适的权重等，更多细节如果看客感兴趣可自行深入了解。

dN/dS 应该是多少

在最简单的 dN/dS 计算中（前文中的 k=1），其中有一个隐性的假设：我们认为每一个位点上每一个突变情况发生的概率都是相等的。

那么，我们继续以 TTT 作为例子，如果每一个突变情况发生的概率都是相等的，理想情况下（无选择压力）我们会获得的 dN/dS 值会是多少？

很简单，发生同义突变的概率是 1/9，发生非同义突变的概率是 8/9，假设发生突变的数量为 x，则通过计算和约分后我们可以发现，dN/dS 值应该是 1。

当然，在现实条件中，每一个突变情况发生的概率并不一定相等，例如，转换（transition）发生的概率就要比颠换（transversion）更高，而很多同义突变都是转换事件，因此如果我们忽略转换和颠换之间的发生概率差异，就会导致 Ssite 低估和 Nsite 的高估，从而导致计算的 dN/dS 偏低。

对于上面这种情况，一种可行的方法就是在计算转换事件时，给它提供的位点数量进行相应加权，也就有了这个表格中 k=2 的情况：

相似的事情还发生在密码子偏好等因素上，而很多方法也针对这些点进行了改进。例如一篇 17 年的 cell 文章中，作者使用 192 个 rate parameters (包括不同的替换类型和链属性，有些类似 3-mer 的突变率模型) 来计算 dN/dS，并发现其他方法要么高估要么低估了 dN/dS。

所以，dN/dS 有一个很关键的设计理念，就是在完全没有选择的情况下，这个值应当等于 1，即同义替代率和非同义替代率应当是相等的。这个概念对于后续的一系列引申也很重要。

dN/dS 有什么用

现在我们能够知道，如果完全考虑了不同替换发生概率的差异并针对它们矫正 Nsite 和 Ssite 的计算，那么 dN/dS 在没有选择的情况下应该是等于 1 的。

在很多地方我们也可以看到这么一个结论：

  • dN/dS < 1，表明存在纯化选择。
  • dN/dS = 1，表明处在中性选择下。
  • dN/dS > 1，表明存在正选择。

也就是说，dN/dS 还可以帮助我们确认一个基因所处的状态，并据此推断出其选择的机制，为我们理解其进化过程提供信息。

那么为什么说 dN/dS = 1 时，基因处于中性选择下呢？这里首先我们要明确一个概念，中性选择并不意味着不存在选择，只是它对生物体的适应性和生存几率没有正面或负面的影响，因此它并不因自然选择而频率增加或下降，在种群内的频率变化全由遗传漂变决定。当然为了方便理解，可以将其视为和 “没有选择” 等效。

同样，中性（neutral）的概念也是如此，一个中性突变就是对生物既无益处也无害处的突变（亦或是带来的益害相抵消）。

然后要明确一个概念，dN/dS = 1 推导出中性选择的前提假设是 —— 同义突变都是中性的。

不妨设想一下，如果同义突变和非同义突变都受到相同强度的强纯化选择，导致其突变数量减少 50%，那么最后算出来的 dN/dS 依然等于 1。

那么这个假设到底成不成立呢？很显然不成立。这个假设主要的根基在于同义突变并不会改变蛋白质的组成和结构，因此从蛋白的角度出发它不会对生物体的表型产生影响。但是在转录翻译阶段，已经有很多研究发现，不同密码子就算编码同一氨基酸，也会存在效率和质量上的不同（这也是密码子偏好的原因之一），所以大多数同义突变其实归根结底是有害的（可见 Shen et al. 2022, Nature）。

但是这影不影响我们使用同义突变作为中性的 Proxy？其实这是个很有争议的问题，因为即使同义突变有害，也很难去量化同义突变到底 “多有害”。今年 4 月份，针对得出 “大多数同义突变都有害” 的文章有人产生了质疑。同样，很多其他的研究表明同义突变可能产生影响的概率比非同义小得多。

违反一个假设的严重性，一要看违反的程度，二要看理论对此的鲁棒性。就以我个人目前进行的分析来看，同义突变受到的影响应该较为轻微。从各种使用 dN/dS 进行选择推断的分析来看，似乎其结果都 “合乎情理”。当然，具体情况如何也有待更详尽的探讨。

总结

其实这篇博客最重要的点还是理清为什么会以 1 这个值作为 dN/dS 的 “分界线”，其实理解了 dN/dS 的计算目的，也就能明白这么做的理由。

这里再总结一下该文章中的一些重点：

①、计算 Nsite 和 Ssite 时，对于突变率高的突变类型，需要对其进行一定的加权以平衡其对 dN 和 dS 计算的贡献。以往的计算中，转换-颠换比率 和 密码子偏好 等都是基于该点考虑对计算进行修正。

②、基于突变发生概率矫正后计算的 dN/dS 可以作为选择情况的参考。如果同义突变的中性假设是正确的或者几近正确的，则当 dN/dS 等于 1 时我们可以认为该蛋白编码基因处于中性选择下。

如果你对代码更感兴趣，并且想知道具体的计算步骤，可以尝试看一看 Biopython 中相应脚本的源码：

https://github.com/biopython/biopython/blob/master/Bio/codonalign/codonseq.py

参考资料

  1. Ziheng Yang, Computational Molecular Evolution.
  2. Martincorena et al., Universal Patterns of Selection in Cancer and Somatic Tissues.
  3. Shen et al., Synonymous mutations in representative yeast genes are mostly strongly non-neutral.
  4. Dhindsa et al., A minimal role for synonymous variation in human disease.