Gotoh 算法

快速概览

Gotoh 算法是处理 affine gap penalty 的标准动态规划算法：它通过引入三个状态矩阵，正确区分 gap 开启和 gap 延长，使比对模型更符合生物学直觉。

• 核心是维护三个独立的 DP 矩阵：匹配状态、X 方向 gap、Y 方向 gap
• 它是所有现代比对工具中 affine gap 模型的理论基础

前置知识

所属板块 核心方法

索引、比对、组装与概率模型构成的核心算法主轴。

阅读目标 帮助建立阅读上下文

先判断这页与你当前问题的关系，再决定是否深入展开。

建议前置 先建立相关基础对象与方法直觉

建议先建立相关基础对象与方法直觉，再进入本页。

是什么

Gotoh 算法是 1982 年提出的用于处理 affine gap penalty 的动态规划算法。它解决了简单线性 gap 模型的问题：

在保持动态规划框架的同时，正确处理 affine gap penalty（gap 开启和 gap 延长代价不同）。

affine gap penalty 更符合生物学直觉：一个长 indel 往往比多个零散小 indel 更合理。

要解决什么生物信息学问题

在简单的 Needleman-Wunsch 算法中，每个 gap 字符的代价是相同的。这导致：

• 一个长度为 5 的连续缺失：5 × g
• 五个分散的 1 bp 缺失：5 × g

两者代价相同，但生物学上前者更常见。

Gotoh 算法引入 affine gap penalty：

$$\text\{GapCost\}(k) = g_o + k \cdot g_e$$

其中：

• g_o：gap opening penalty（开启 gap 的代价）
• g_e：gap extension penalty（延长 gap 的代价）
• k：gap 长度

这样，连续 gap 的代价增长更慢，更符合真实生物学过程。

数学模型

三状态设计

Gotoh 算法维护三个矩阵：

M[i][j]

以匹配或替换（非 gap）结束的最优得分

Ix[i][j]

以在序列 X 中插入 gap（即删除 x_i）结束的最优得分

Iy[i][j]

以在序列 Y 中插入 gap（即插入 y_j）结束的最优得分

边界条件

\begin\{aligned\} M[0][0] &= 0 \\ M[i][0] &= -\infty, \quad i > 0 \\ M[0][j] &= -\infty, \quad j > 0 \\ Ix[0][0] &= -\infty \\ Ix[i][0] &= g_o + i \cdot g_e, \quad i > 0 \\ Ix[0][j] &= -\infty, \quad j > 0 \\ Iy[0][0] &= -\infty \\ Iy[i][0] &= -\infty, \quad i > 0 \\ Iy[0][j] &= g_o + j \cdot g_e, \quad j > 0 \end\{aligned\}

边界条件的含义：

• M 状态不能以 gap 开始，所以当任一序列为空时得分为 -∞
• Ix 状态表示在 X 中插入 gap，所以当 Y 为空时可以开启 gap
• Iy 状态表示在 Y 中插入 gap，所以当 X 为空时可以开启 gap

核心递推公式

当 i > 0 且 j > 0 时：

匹配状态 M：

$$M[i][j] = \max \begin{cases} M[i-1][j-1] + s(x_i, y_j) \\ Ix[i-1][j-1] + s(x_i, y_j) \\ Iy[i-1][j-1] + s(x_i, y_j) \end{cases}$$

X 方向 gap 状态 Ix：

$$Ix[i][j] = \max \begin{cases} M[i-1][j] + g_o + g_e & \text{开启新 gap} \\ Ix[i-1][j] + g_e & \text{延长已有 gap} \end{cases}$$

Y 方向 gap 状态 Iy：

$$Iy[i][j] = \max \begin{cases} M[i][j-1] + g_o + g_e & \text{开启新 gap} \\ Iy[i][j-1] + g_e & \text{延长已有 gap} \end{cases}$$

最终得分

$$F[i][j] = \max \\{ M[i][j], Ix[i][j], Iy[i][j] \\}$$

算法步骤

步骤 1：初始化三个矩阵

创建 (m+1) × (n+1) 的三个矩阵 M、Ix、Iy，并填充边界条件。

步骤 2：填充矩阵

按照行优先或列优先的顺序，使用递推公式填充三个矩阵。

步骤 3：得到最终得分

最终得分为 F[m][n] = max(M[m][n], Ix[m][n], Iy[m][n])。

步骤 4：回溯得到比对

从最终得分所在的状态开始回溯：

1. 如果最终来自 M[i][j]：
   • 将 x_i 和 y_j 对齐
   • 根据来源回溯到 M[i-1][j-1]、Ix[i-1][j-1] 或 Iy[i-1][j-1]
2. 如果最终来自 Ix[i][j]：
   • 将 x_i 与 gap 对齐
   • 根据来源回溯到 M[i-1][j]（开启 gap）或 Ix[i-1][j]（延长 gap）
3. 如果最终来自 Iy[i][j]：
   • 将 gap 与 y_j 对齐
   • 根据来源回溯到 M[i][j-1]（开启 gap）或 Iy[i][j-1]（延长 gap）

示例

考虑两条序列：

• X = GATT（长度 m = 4）
• Y = GCAT（长度 n = 4）

使用以下打分体系：

• 匹配得分：+2
• 错配得分：-1
• Gap opening penalty：g_o = -3
• Gap extension penalty：g_e = -1

步骤 1：初始化

M 矩阵：

	ε	G	C	A	T
ε	0	-∞	-∞	-∞	-∞
G	-∞
A	-∞
T	-∞
T	-∞

Ix 矩阵：

	ε	G	C	A	T
ε	-∞	-∞	-∞	-∞	-∞
G	-4
A	-5
T	-6
T	-7

Iy 矩阵：

	ε	G	C	A	T
ε	-∞	-4	-5	-6	-7
G	-∞
A	-∞
T	-∞
T	-∞

步骤 2：填充矩阵

以 M[1][1] 为例（G 对 G）：

• 来自 M[0][0] + s(G,G) = 0 + 2 = 2
• 来自 Ix[0][0] + s(G,G) = -∞ + 2 = -∞
• 来自 Iy[0][0] + s(G,G) = -∞ + 2 = -∞
• M[1][1] = 2

以 Ix[1][1] 为例：

• 开启 gap：M[0][1] + g_o + g_e = -∞ + (-3) + (-1) = -∞
• 延长 gap：Ix[0][1] + g_e = -∞ + (-1) = -∞
• Ix[1][1] = -∞

以 Iy[1][1] 为例：

• 开启 gap：M[1][0] + g_o + g_e = -∞ + (-3) + (-1) = -∞
• 延长 gap：Iy[1][0] + g_e = -∞ + (-1) = -∞
• Iy[1][1] = -∞

完整矩阵（为简洁起见，用 -9 表示 -∞）：

M 矩阵：

	ε	G	C	A	T
ε	0	-9	-9	-9	-9
G	-9	2	-2	-1	-3
A	-9	-2	1	4	2
T	-9	-3	0	2	6
T	-9	-4	-1	1	5

Ix 矩阵：

	ε	G	C	A	T
ε	-9	-9	-9	-9	-9
G	-4	-9	-9	-9	-9
A	-5	-6	-9	-9	-9
T	-6	-7	-8	-9	-9
T	-7	-8	-9	-10	-9

Iy 矩阵：

	ε	G	C	A	T
ε	-9	-4	-5	-6	-7
G	-9	-9	-6	-7	-8
A	-9	-9	-9	-6	-7
T	-9	-9	-9	-9	-6
T	-9	-9	-9	-9	-9

步骤 3：得到最终得分

F[4][4] = max(M[4][4], Ix[4][4], Iy[4][4]) = max(5, -9, -9) = 5

步骤 4：回溯

从 M[4][4] = 5 开始回溯。根据 M 矩阵数据：

• 第3行（T）：-9 | -3 | 0 | 2 | 6
• 第4行（T）：-9 | -4 | -1 | 1 | 5

回溯过程（选择得分递增的对角线路径）：

1. M[4][4] = 5，可能来自 M[3][3] + s(T,T) = 2 + 2 = 4（不精确匹配，可能有其他路径）
2. 继续向前追溯：检查 M[3][3] = 2 的来源
3. 直至到达 M[0][0] = 0

说明：由于示例矩阵中部分值被简化表示（-9 代表 -∞），完整回溯需要仔细验证每个单元格的实际计算值。

一个可能的比对结果：

GATT
| ||
GCAT

得分计算（假设无 gap）：

• G-G: +2
• A-C: -1（错配）
• T-A: -1（错配）
• T-T: +2 总计：+2

复杂度分析

时间复杂度

• 填充三个矩阵：O(mn)，每个单元格需要常数时间计算
• 回溯：O(m + n)
• 总时间复杂度：O(mn)

空间复杂度

• 存储三个矩阵：O(mn)
• 如果只需要得分不需要回溯：可以优化到 O(min(m, n))（每个矩阵只保留两行）

与简单线性 gap 的对比

维度	线性 gap	Affine gap (Gotoh)
Gap 代价	`k × g`	`g_o + k × g_e`
状态数	1	3
空间复杂度	O(mn)	3 × O(mn)
生物学合理性	较低	较高

与真实工具或流程的连接

Gotoh 算法是：

• 所有现代比对工具中 affine gap 模型的理论基础
• Needleman-Wunsch 的自然扩展，增加了更真实的 gap 模型
• 蛋白质比对中特别重要，因为蛋白序列中的 indel 往往是连续的

现代工具如 BWA、minimap2、BLAST 等都使用 affine gap penalty，其核心思想来自 Gotoh 算法。

常见误区

• Affine gap 只是参数调整：
• 不是。它改变了动态规划的状态设计和递推关系，是算法结构的根本变化。
• Gotoh 算法总是优于线性 gap：
• 不一定。如果数据特征不符合 affine gap 假设，或者参数设置不当，结果可能不如线性 gap。
• 三个状态矩阵可以简化：
• 不能。三个状态是处理 affine gap 的必要条件，简化会失去正确性。

算法变体

局部比对版本

将 Gotoh 算法扩展到局部比对，类似于 Smith-Waterman：

在每个递推公式中加入 0 选项，允许从任意位置重新开始比对。

线性空间优化

类似 Hirschberg 算法，可以将 Gotoh 的空间复杂度从 O(mn) 优化到 O(min(m, n))，但实现更复杂。

参考资料

• Gotoh, O. (1982). An improved algorithm for matching biological sequences. Journal of molecular biology, 162(3), 705-708.
• Durbin, R., Eddy, S., Krogh, A., & Mitchison, G. (1998). Biological sequence analysis: probabilistic models of proteins and nucleic acids. Cambridge university press.

相关页面

• 编辑距离
• 全局比对与局部比对
• 打分矩阵与 gap 罚分
• Affine gap penalty
• Needleman-Wunsch 算法

Needleman-Wunsch 算法

Gotoh 算法的基础框架，处理线性 gap 罚分的全局比对

延伸阅读

Smith-Waterman 算法

Gotoh 思想同样可扩展到局部比对场景

延伸阅读

Affine gap penalty

Gotoh 算法所处理的 affine gap 模型的详细说明

延伸阅读

打分矩阵与 gap 罚分

比对打分体系的全面介绍

延伸阅读

编辑距离

理解序列比对 DP 递推关系的最简起点

延伸阅读