动态规划算法

快速概览

动态规划是生物信息学最核心的算法范式之一，通过将复杂问题分解为重叠子问题，系统地寻找全局最优解。

理解动态规划的关键在于设计递推公式
状态定义和边界条件比具体算法更重要

所属板块 分析方向与案例

把基础对象与算法方法重新放回真实分析任务与工作流。

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

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

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

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

是什么

动态规划（Dynamic Programming, DP）是一种将复杂问题分解为重叠子问题进行求解的算法范式，其核心思想包括：

最优子结构：问题的最优解包含其子问题的最优解
重叠子问题：子问题会被重复计算，可以存储结果避免重复
无后效性：当前状态只依赖于之前的状态，与如何到达当前状态无关

在生物信息学中，动态规划主要用于：

序列比对（全局、局部、半全局）
编辑距离计算
RNA 二级结构预测
HMM 中的 Viterbi、Forward、Backward 算法

为什么适合生物信息学

生物序列分析问题天然适合动态规划：

序列的线性结构：DNA、RNA、蛋白质序列均为线性排列，易于定义子问题
局部相似性：同源序列往往在局部区域相似，适合逐步比对
明确的递推关系：比对得分可从子比对得分递推得到
可解释性：DP 矩阵可直观展示比对过程

核心要素

1. 状态定义

状态定义是 DP 的首要步骤，必须满足：

完整性：状态必须包含求解问题所需的全部信息
无冗余：状态不应包含无关信息
可递推：状态之间必须有明确的转移关系

示例：序列比对

F[i][j] = 序列 X 的前 i 个字符与序列 Y 的前 j 个字符的最优比对得分

2. 递推公式

递推公式描述如何从子问题状态计算当前状态。

通用形式：

\text{当前状态} = \text{转移函数}(\text{前驱状态})

序列比对递推公式:

F[i][j] = \max \begin{cases} F[i-1][j-1] + s(x_i, y_j) & \text{匹配或替换} \\ F[i-1][j] + g & \text{在 X 中插入 gap} \\ F[i][j-1] + g & \text{在 Y 中插入 gap} \end{cases}

3. 边界条件

边界条件用于初始化 DP 表格的边缘，通常对应空序列或初始状态。

序列比对边界：

\begin{aligned} F[0][0] &= 0 \\ F[i][0] &= i \cdot g \quad \text{（前 i 个字符与空串对齐）} \\ F[0][j] &= j \cdot g \quad \text{（空串与前 j 个字符对齐）} \end{aligned}

4. 回溯策略

DP 表格填充完成后，需从最终状态回溯至初始状态以重构解。

回溯规则：

记录每个状态的最优前驱状态
从最终状态开始，按照前驱指针回溯
回溯路径即为最优解

经典应用

1. 全局序列比对（Needleman-Wunsch）

问题：找到两条序列之间得分最高的全局比对。

Worked Example：

序列 X = GATTACA，Y = GCATGCU

打分：匹配 +1，错配 -1，gap -2

步骤 1：初始化

	ε	G	C	A	T	G	C	U
ε	0	-2	-4	-6	-8	-10	-12	-14
G	-2
A	-4
T	-6
T	-8
A	-10
C	-12
A	-14

步骤 2：填充矩阵

计算 F[1][1]（G vs G）：

匹配：F[0][0] + s(G,G) = 0 + 1 = 1
X 插入 gap：F[0][1] + g = -2 + (-2) = -4
Y 插入 gap：F[1][0] + g = -2 + (-2) = -4
取最大值：F[1][1] = 1

完整矩阵：

	ε	G	C	A	T	G	C	U
ε	0	-2	-4	-6	-8	-10	-12	-14
G	-2	1	-1	-3	-5	-3	-5	-7
A	-4	-1	0	0	-2	-4	-6	-8
T	-6	-3	-2	-1	1	-1	-3	-5
T	-8	-5	-4	-3	-1	0	-2	-4
A	-10	-7	-6	-2	-3	-2	-1	-3
C	-12	-9	-5	-4	-4	-3	-1	-3
A	-14	-11	-7	-3	-5	-4	-2	0

步骤 3：回溯

从 F[7][7] = 0 开始回溯，得到最优比对：

G-A-T-T-A-C-A
| | |   | | |
G-C-A-T-G-C-U

得分：1 -1 -1 -2 -1 -1 -1 = -6

矩阵填充与回溯图解逻辑：

初始化：第一行和第一列根据空位罚分填充。
填充过程：每个格子 $F(i,j)$ 观察左、上、左上三个邻居。
视觉流向：
- 对角线箭头：表示匹配（Match）或错配（Mismatch），序列同时推进。
- 垂直箭头：表示在水平序列中引入空位（Insertion）。
- 水平箭头：表示在垂直序列中引入空位（Deletion）。
回溯路径：从右下角开始，逆着最优来源方向追溯，直到左上角起点。路径上的每一步决定了比对字符串的具体形式。

详见：Needleman-Wunsch 算法详解

2. 局部序列比对（Smith-Waterman）

问题：找到两条序列中相似度最高的局部片段。

与全局比对的关键区别:

边界条件不同：F[i][0] = F[0][j] = 0
递推增加一个选项：F[i][j] = 0（从空开始）
回溯从矩阵中的最大值开始，而非 F[m][n]

递推公式：

F[i][j] = \max \begin{cases} 0 & \text{重新开始} \\ F[i-1][j-1] + s(x_i, y_j) & \text{匹配或替换} \\ F[i-1][j] + g & \text{在 X 中插入 gap} \\ F[i][j-1] + g & \text{在 Y 中插入 gap} \end{cases}

3. 编辑距离（Levenshtein Distance）

问题：计算将一个序列转换为另一个序列所需的最少编辑操作数。

操作类型：插入、删除、替换（每次代价为 1）

递推公式：

D[i][j] = \min \begin{cases} D[i-1][j-1] & \text{如果 } x_i = y_j \\ D[i-1][j-1] + 1 & \text{如果 } x_i \neq y_j \text{（替换）} \\ D[i-1][j] + 1 & \text{删除 } x_i \\ D[i][j-1] + 1 & \text{插入 } y_j \end{cases}

Worked Example：

序列 X = kitten，Y = sitting

	ε	s	i	t	t	i	n	g
ε	0	1	2	3	4	5	6	7
k	1	1	2	3	4	5	6	7
i	2	2	1	2	3	4	5	6
t	3	3	2	1	2	3	4	5
t	4	4	3	2	1	2	3	4
e	5	5	4	3	2	2	3	4
n	6	6	5	4	3	3	2	3

编辑距离 = D[7][7] = 3

编辑路径：

k → s（替换）
e → i（替换）
插入 g

4. RNA 二级结构预测（Nussinov 算法）

问题：给定 RNA 序列，预测其二级结构中的最大碱基配对数。

配对规则：A-U、G-C、G-U（wobble 配对）

状态定义：

N[i][j] = 序列从位置 i 到 j 的最大配对数

递推公式：

N[i][j] = \max \begin{cases} N[i+1][j] & \text{i 不参与配对} \\ N[i][j-1] & \text{j 不参与配对} \\ N[i+1][j-1] + s(i,j) & \text{i 与 j 配对} \\ \max_{i<k<j} \{N[i][k] + N[k+1][j]\} & \text{分叉} \end{cases}

其中 s(i,j) = 1（若 i 和 j 可配对），否则为 0。

Worked Example：

序列 = GGAC

	1(G)	2(G)	3(A)	4(C)
1(G)	0	0	0	1
2(G)		0	0	1
3(A)			0	0
4(C)				0

最优配对：G(1)-C(4)，配对数为 1

5. HMM Viterbi 算法

问题：给定观测序列与 HMM 模型，寻找最可能的隐状态序列。

状态定义：

V[t][i] = 在时刻 t 处于状态 i 的最大概率路径得分

递推公式：

V[t][i] = \max_j \\{V[t-1][j] \cdot a_\{j,i\}\\} \cdot b_i(o_t)

其中：

$a_{j,i}$ 为从状态 j 到状态 i 的转移概率
$b_i(o_t)$ 为状态 i 产生观测 $o_t$ 的发射概率

初始化：

V[1][i] = \pi_i \cdot b_i(o_1)

其中 $\pi_i$ 是初始状态概率。

复杂度分析

时间复杂度

算法	时间复杂度	说明
全局比对	O(mn)	m, n 为序列长度
局部比对	O(mn)	需要搜索最大值
编辑距离	O(mn)	与比对类似
RNA 结构	O(n³)	bifurcation 增加一维
HMM Viterbi	O(T·K²)	T 为时间步，K 为状态数

空间复杂度

策略	空间复杂度	说明
完整存储	O(mn)	存储整个 DP 表格
线性空间	O(min(m,n))	只存储当前行/列
Hirschberg	O(min(m,n))	分治策略，可回溯

线性空间优化：

若仅需得分，可只保留两行或两列
若需回溯，可使用 Hirschberg 算法

算法优化

1. 带状动态规划（Banded DP）

思想：若已知两条序列相似度较高，可限制搜索带宽。

复杂度：从 O(mn) 降至 O(b·min(m,n))，b 为带宽。

适用场景：

高度相似序列的比对
重测序数据分析
已知大致比对位置的精化

2. 自适应带宽

思想：根据局部相似性动态调整带宽。

优势：

相似区域使用窄带宽
差异区域使用宽带宽
平衡精度与效率

3. 剪枝策略

思想：提前剪除不可能成为最优解的分支。

示例：

Ukkonen 算法用于编辑距离
限制最大 gap 长度
设置得分阈值

4. 并行化

策略：

行/列并行填充
分块计算
GPU 加速

挑战：

依赖关系限制了并行度
回溯阶段难以并行

与真实工具的连接

BWA-MEM

BWA-MEM 虽主要采用 seed-and-extend 策略，但在扩展阶段使用局部动态规划进行精确比对。

GATK HaplotypeCaller

在局部区域使用类 Smith-Waterman 算法进行 haplotype 比对。

RNAfold

使用动态规划（Zuker 算法）预测 RNA 二级结构，考虑热力学稳定性而非仅最大化配对数。

BLAST

虽然 BLAST 主要使用启发式方法，但其核心思想可以追溯到动态规划。

扩展阅读

高级主题

仿射 gap penalty：gap 开放和延伸使用不同罚分
双序列比对到多序列比对：从 2D DP 到高维 DP 的挑战
近似算法：当精确 DP 不可行时的近似策略
概率 DP：与 HMM、贝叶斯推断的结合

实现技巧

使用滚动数组减少内存
位运算优化
缓存友好性设计
数值稳定性处理

练习

手动计算 AGCT 和 ACGT 的全局比对得分（匹配 +1，错配 -1，gap -2）
证明编辑距离满足三角不等式
设计一个带状 DP 算法，带宽为 3
比较 Needleman-Wunsch 和 Smith-Waterman 在不同场景下的适用性

参考资料

Needleman, S. B., & Wunsch, C. D. (1970). A general method applicable to the search for similarities in the amino acid sequence of two proteins. Journal of molecular biology, 48(3), 443-453.
Smith, T. F., & Waterman, M. S. (1981). Identification of common molecular subsequences. Journal of molecular biology, 147(1), 195-197.
Durbin, R., Eddy, S., Krogh, A., & Mitchison, G. (1998). Biological sequence analysis: probabilistic models of proteins and nucleic acids. Cambridge university press.

	ε	G	C	A	T	G	C	U
ε	0	-2	-4	-6	-8	-10	-12	-14
G	-2	1	-1	-3	-5	-3	-5	-7
A	-4	-1	0	0	-2	-4	-6	-8
T	-6	-3	-2	-1	1	-1	-3	-5
T	-8	-5	-4	-3	-1	0	-2	-4
A	-10	-7	-6	-2	-3	-2	-1	-3
C	-12	-9	-5	-4	-4	-3	-1	-3
A	-14	-11	-7	-3	-5	-4	-2	0

	ε	s	i	t	t	i	n	g
ε	0	1	2	3	4	5	6	7
k	1	1	2	3	4	5	6	7
i	2	2	1	2	3	4	5	6
t	3	3	2	1	2	3	4	5
t	4	4	3	2	1	2	3	4
e	5	5	4	3	2	2	3	4
n	6	6	5	4	3	3	2	3

	ε	G	C	A	T	G	C	U
ε	0	-2	-4	-6	-8	-10	-12	-14
G	-2	1	-1	-3	-5	-3	-5	-7
A	-4	-1	0	0	-2	-4	-6	-8
T	-6	-3	-2	-1	1	-1	-3	-5
T	-8	-5	-4	-3	-1	0	-2	-4
A	-10	-7	-6	-2	-3	-2	-1	-3
C	-12	-9	-5	-4	-4	-3	-1	-3
A	-14	-11	-7	-3	-5	-4	-2	0

	ε	s	i	t	t	i	n	g
ε	0	1	2	3	4	5	6	7
k	1	1	2	3	4	5	6	7
i	2	2	1	2	3	4	5	6
t	3	3	2	1	2	3	4	5
t	4	4	3	2	1	2	3	4
e	5	5	4	3	2	2	3	4
n	6	6	5	4	3	3	2	3

动态规划算法

是什么

为什么适合生物信息学

核心要素

1. 状态定义

2. 递推公式

3. 边界条件

4. 回溯策略

经典应用

1. 全局序列比对（Needleman-Wunsch）

2. 局部序列比对（Smith-Waterman）

3. 编辑距离（Levenshtein Distance）

4. RNA 二级结构预测（Nussinov 算法）

5. HMM Viterbi 算法

复杂度分析

时间复杂度

空间复杂度

算法优化

1. 带状动态规划（Banded DP）

2. 自适应带宽

3. 剪枝策略

4. 并行化

与真实工具的连接

BWA-MEM

GATK HaplotypeCaller

RNAfold

BLAST

扩展阅读

高级主题

实现技巧

练习

参考资料

相关页面

Needleman-Wunsch 算法

Smith-Waterman 算法

编辑距离

HMM 与 Viterbi 算法

算法总览

	ε	G	C	A	T	G	C	U
ε	0	-2	-4	-6	-8	-10	-12	-14
G	-2	1	-1	-3	-5	-3	-5	-7
A	-4	-1	0	0	-2	-4	-6	-8
T	-6	-3	-2	-1	1	-1	-3	-5
T	-8	-5	-4	-3	-1	0	-2	-4
A	-10	-7	-6	-2	-3	-2	-1	-3
C	-12	-9	-5	-4	-4	-3	-1	-3
A	-14	-11	-7	-3	-5	-4	-2	0

	ε	s	i	t	t	i	n	g
ε	0	1	2	3	4	5	6	7
k	1	1	2	3	4	5	6	7
i	2	2	1	2	3	4	5	6
t	3	3	2	1	2	3	4	5
t	4	4	3	2	1	2	3	4
e	5	5	4	3	2	2	3	4
n	6	6	5	4	3	3	2	3