曼哈顿旅游问题
想象你是一个在曼哈顿旅游的游客,想从(0,0) 走到(n,m),每条街道(边)都有不同数量的景点。你的目标是走过景点总数最多的路径。
- 将序列比对问题转化为网格图中的路径搜寻
- 引入了二维动态规划的状态转移思想
- 通过 DAG(有向无环图)处理非规则网格
- 掌握路径回溯(Backtracking)恢复最优解
1. 是什么
Section titled “1. 是什么”曼哈顿旅游问题(Manhattan Tourist Problem)是一个经典的动态规划问题。它将”在网格中寻找权重最大的路径”这一任务形式化,是理解序列比对算法最直观的数学模型。
问题设定:给定一个 的网格,每条水平边和垂直边都有一个权重(景点数)。游客只能向右或向下移动。
目标:找到一条从左上角 到右下角 的路径,使路径上边权重之和最大。
2. 要解决什么生物信息学问题
Section titled “2. 要解决什么生物信息学问题”曼哈顿旅游问题是序列比对(Sequence Alignment)的抽象模型。在生物信息学中,我们需要比较两条 DNA、RNA 或蛋白质序列,找出它们之间的最优对应关系。这个问题的搜索空间随序列长度呈指数增长,而动态规划将其降为多项式时间。
通过将序列比对问题映射为网格路径问题,曼哈顿旅游问题提供了以下生物信息学核心任务的统一框架:
- 全局比对:对齐两条完整的序列(Needleman-Wunsch 算法)
- 局部比对:寻找两条序列中最相似的子区域(Smith-Waterman 算法)
- 重叠比对:寻找两条序列末端的重叠区域(用于基因组组装)
3. 输入与输出
Section titled “3. 输入与输出”输入:
- 网格的维度 和
- 水平边权重矩阵 :从 向右走到 的边权重
- 垂直边权重矩阵 :从 向下走到 的边权重
输出:
- 从 到 的最大权重路径的权重值
- (可选)最大权重路径本身(通过回溯恢复)
4. 核心思想与数学模型
Section titled “4. 核心思想与数学模型”动态规划的基石是最优性原理(Principle of Optimality):从 到 的最优路径,必然包含从 到其经过的每一个中间点 的最优路径。
这意味着我们可以自底向上地求解:先求出到达每个中间点的最优值,再用这些中间结果推导最终答案。
设 为到达交叉点 时能收集到的最大景点数。
要到达 ,只能从上方 或左方 过来:
- (起点没有收集任何景点)
- 第一行只能从左边推导:,对于
- 第一列只能从上面推导:,对于
5. 关键递推与伪代码
Section titled “5. 关键递推与伪代码”MANHATTANTOURIST(down, right, n, m)1 s[0,0] ← 02 for i ← 1 to n3 s[i,0] ← s[i-1,0] + down[i,0]4 for j ← 1 to m5 s[0,j] ← s[0,j-1] + right[0,j]6 for i ← 1 to n7 for j ← 1 to m8 s[i,j] ← max(s[i-1,j] + down[i,j], s[i,j-1] + right[i,j])9 return s[n,m]6. Worked Example
Section titled “6. Worked Example”考虑一个 的网格(即 ),定义权重如下:
垂直边权重 down(从上方到下方):
| j=0 | j=1 | j=2 | j=3 | |
|---|---|---|---|---|
| i=1 | 1 | 2 | 3 | 1 |
| i=2 | 2 | 1 | 2 | 3 |
| i=3 | 3 | 2 | 1 | 2 |
水平边权重 right(从左方到右方):
| j=1 | j=2 | j=3 | |
|---|---|---|---|
| i=0 | 3 | 2 | 4 |
| i=1 | 1 | 4 | 2 |
| i=2 | 2 | 1 | 3 |
| i=3 | 1 | 3 | 2 |
第一步:初始化
Section titled “第一步:初始化”s[0,0] = 0
第一列(只能从上方来):s[1,0] = 0 + down[1,0] = 0 + 1 = 1s[2,0] = 1 + down[2,0] = 1 + 2 = 3s[3,0] = 3 + down[3,0] = 3 + 3 = 6
第一行(只能从左方来):s[0,1] = 0 + right[0,1] = 0 + 3 = 3s[0,2] = 3 + right[0,2] = 3 + 2 = 5s[0,3] = 5 + right[0,3] = 5 + 4 = 9第二步:填充内部
Section titled “第二步:填充内部”s[1,1] = max(s[0,1] + down[1,1], s[1,0] + right[1,1]) = max(3 + 2, 1 + 1) = max(5, 2) = 5 ← 从上方
s[1,2] = max(s[0,2] + down[1,2], s[1,1] + right[1,2]) = max(5 + 3, 5 + 4) = max(8, 9) = 9 ← 从左方
s[1,3] = max(s[0,3] + down[1,3], s[1,2] + right[1,3]) = max(9 + 1, 9 + 2) = max(10, 11) = 11 ← 从左方
s[2,1] = max(s[1,1] + down[2,1], s[2,0] + right[2,1]) = max(5 + 1, 3 + 2) = max(6, 5) = 6 ← 从上方
s[2,2] = max(s[1,2] + down[2,2], s[2,1] + right[2,2]) = max(9 + 2, 6 + 1) = max(11, 7) = 11 ← 从上方
s[2,3] = max(s[1,3] + down[2,3], s[2,2] + right[2,3]) = max(11 + 3, 11 + 3) = max(14, 14) = 14 ← 两者等价
s[3,1] = max(s[2,1] + down[3,1], s[3,0] + right[3,1]) = max(6 + 2, 6 + 1) = max(8, 7) = 8 ← 从上方
s[3,2] = max(s[2,2] + down[3,2], s[3,1] + right[3,2]) = max(11 + 1, 8 + 3) = max(12, 11) = 12 ← 从上方
s[3,3] = max(s[2,3] + down[3,3], s[3,2] + right[3,3]) = max(14 + 2, 12 + 2) = max(16, 14) = 16 ← 从上方完整的 DP 表:
j=0 j=1 j=2 j=3i=0 0 3 5 9i=1 1 5 9 11i=2 3 6 11 14i=3 6 8 12 16最大权重路径:从 到 的最大权重为 。
通过从终点 反向追踪每一步的选择,可以恢复最优路径:
(3,3) ← 上方(2,3) ← 上方(1,3) ← 左方(1,2) ← 上方(0,2) ← 左方(0,1) ← 左方(0,0)路径:
7. 从网格到有向无环图(DAG)
Section titled “7. 从网格到有向无环图(DAG)”曼哈顿是一个规则网格,但现实中(如百老汇大道纵切曼哈顿)或生物学问题中,路径可能更复杂。我们可以将任何网格或比对图抽象为 DAG(Directed Acyclic Graph, 有向无环图)。
- 入度(Indegree)
- 进入一个顶点的边数。在 DP 中,入度决定了该顶点的前驱数量。
- 出度(Outdegree)
- 从一个顶点出去的边数。
- 前驱(Predecessors)
- 所有能直接到达顶点 $v$ 的顶点集合 $P(v)$。DP 递推中需要遍历所有前驱。
DAG 中的最长路径递推
Section titled “DAG 中的最长路径递推”对于任意顶点 ,其得分 为:
这是一个高度通用的递推公式。在曼哈顿网格中,,每个顶点恰好有 2 个前驱。在通用 DAG 中,前驱的数量可以任意。
拓扑排序(Topological Ordering)
Section titled “拓扑排序(Topological Ordering)”为了确保在计算 时,其所有前驱 已经计算完毕,我们需要按照拓扑排序的顺序遍历顶点。
定义:拓扑排序是 DAG 中顶点的一个线性排列,使得对于每条有向边 , 在排列中出现在 之前。
- 在规则网格中,按行、按列或按对角线遍历都是合法的拓扑排序。
- 在通用 DAG 中,拓扑排序保证了所有边的方向都是从序号小的顶点指向序号大的顶点。
- 拓扑排序可以使用 Kahn 算法(基于入度的 BFS)或 DFS 在 时间内完成。
DAG 拓扑排序伪代码
Section titled “DAG 拓扑排序伪代码”TOPOLOGICALSORT(G)1 计算每个顶点的入度2 将所有入度为 0 的顶点放入队列 Q3 order ← []4 while Q 不为空5 v ← Q.dequeue()6 将 v 追加到 order7 for each 邻接顶点 w of v8 删除边(v, w)9 if w 的入度变为 010 Q.enqueue(w)11 return order8. 复杂度与适用前提
Section titled “8. 复杂度与适用前提”| 步骤 | 时间复杂度 | 空间复杂度 |
|---|---|---|
| 初始化边界 | ||
| 填充 DP 表 | ||
| 路径回溯 | (路径长度) | |
| 总计 |
对于序列比对问题,如果两个序列长度分别为 和 ,则复杂度为 。对于人类基因组级别的数据(),完整的 DP 表需要约 个单元格,这在实践中是不可行的。这正是 BWA、minimap2 等工具使用种子-扩展(Seed-and-Extend)策略来避免完整 DP 的原因。
如果只需要最优路径的权重而不需要路径本身,可以将空间复杂度从 降低到 ,因为计算第 行时只需要第 行的结果。
如果需要恢复路径本身,可以使用 Hirschberg 算法,它利用分治法在 时间和 空间内完成。
动态规划在此问题上有效的前提条件:
- 最优子结构:到 的最优路径包含到其前驱的最优路径。
- 无环性:图中没有环。如果存在环(正向权重环),则不存在”最长路径”(可以无限绕环增加权重)。
- 有向性:边有明确的方向,不允许回头。
9. 与真实工具的连接
Section titled “9. 与真实工具的连接”曼哈顿旅游问题几乎就是**序列比对(Sequence Alignment)**的代名词。通过在网格中引入对角线移动,它可以自然地扩展为经典的比对算法。
从网格到比对
Section titled “从网格到比对”| 网格元素 | 比对含义 | 分数 |
|---|---|---|
| 顶点 | 序列 1 的前 个字符与序列 2 的前 个字符对齐 | — |
| 水平边(向右) | 在序列 2 中插入一个空位(Gap) | (空位罚分) |
| 垂直边(向下) | 在序列 1 中插入一个空位(Gap) | (空位罚分) |
| 对角线移动 | 匹配(Match)或错配(Mismatch) | (匹配)或 (错配) |
与经典算法的关系
Section titled “与经典算法的关系”- Needleman-Wunsch(1970):在曼哈顿网格基础上增加对角线移动,用于全局比对。递推关系变为三个方向的取最大值。
- Smith-Waterman(1981):增加”重新从零开始”的选项( 可以取 ),用于局部比对。
- Gotoh(1982):引入仿射空隙罚分(Affine Gap Penalty),将 DP 矩阵扩展为三个矩阵,分别处理匹配/错配、空位开启和空位延伸。
现代工具中的 DP
Section titled “现代工具中的 DP”现代比对工具并不直接运行完整的 Needleman-Wunsch DP(对长序列来说太慢),而是使用种子-扩展策略:
- 种子阶段:通过索引(如 FM-Index)快速找到精确匹配的短片段(种子)。
- 扩展阶段:在种子周围使用带空隙罚分的 DP 进行局部精确对齐。
- 评分阶段:根据 DP 得分和映射质量(MapQ)评估比对结果的可信度。
这一策略的核心仍然是动态规划,只是将其限制在局部小范围内运行,从而实现了在基因组尺度上的高效比对。