Skip to main content

2.2.2 A WFA Example


下面用一个 6 个碱基的 DNA 序列比对 例子，完整走一遍 WFA 的波前扩展和回溯。

1. 选例子​

设两条序列：

R = A C G T A C
Q = A C G T C A

索引：

R: R1=A  R2=C  R3=G  R4=T  R5=A  R6=C
Q: Q1=A  Q2=C  Q3=G  Q4=T  Q5=C  Q6=A

采用最简单的 Levenshtein 距离：

  • 匹配：代价 0
  • 错配：代价 1
  • 插入：代价 1
  • 删除：代价 1

目标：把 R 变成 Q，最少需要几次编辑。

肉眼观察：

R: A C G T A C
Q: A C G T C A

前 4 个 ACGT 完全一样。后面：

  • R 第 5 个是 A，第 6 个是 C
  • Q 第 5 个是 C，第 6 个是 A

可以这样做：

  1. 删除 R5=A
  2. 在末尾插入 Q6=A

所以编辑距离是 2。

下面看 WFA 怎么按分数一层层找到这个结果。

2. WFA 基本记号​

定义 DP 矩阵中的位置 (i, j)：

  • i：已经用了 R 的前 i 个字符
  • j：已经用了 Q 的前 j 个字符

定义对角线：

k=i−jk = i - j

主对角线是 k = 0。

WFA 用波前 X_s(k) 表示：

在分数 s 下，对角线 k 上能到达的最远 i 坐标。

因为同一条对角线上走得越远越好，所以每条对角线只记最远位置。

关键规则：

  • 匹配：留在同一条对角线，免费前进。
  • 错配：留在同一条对角线，分数 +1，位置前进 1。
  • 插入：跳到 k-1 对角线，分数 +1。
  • 删除：跳到 k+1 对角线，分数 +1。

每次编辑操作后，继续沿对角线做匹配延伸，直到不能匹配为止。

3. 分数 s = 0：只靠匹配能走多远​

从 (0,0) 出发，对角线 k = 0。

依次比较：

R1=A  vs  Q1=A   匹配
R2=C  vs  Q2=C   匹配
R3=G  vs  Q3=G   匹配
R4=T  vs  Q4=T   匹配
R5=A  vs  Q5=C   不匹配

所以 s = 0 时，主对角线 k = 0 最远到达：

(4,4)

波前表：

分数 s对角线 k最远位置 (i,j)
00(4,4)

终点是 (6,6)，还没到。

4. 分数 s = 1：允许一次编辑​

从 s = 0 的 (4,4) 出发，考虑一次编辑操作。

4.1 错配：留在 k = 0​

R5=A 和 Q5=C 错配，分数 +1，到达：

(5,5)

然后继续匹配：

R6=C  vs  Q6=A   不匹配

所以 k = 0 上最远到 (5,5)。

4.2 删除：跳到 k = 1​

删除 R5=A，相当于 R 前进、Q 不动：

从 (4,4) 到 (5,4)，对角线：

k=5−4=1k = 5 - 4 = 1

然后继续匹配：

R6=C  vs  Q5=C   匹配

到达 (6,5)。此时 R 结束，Q 还剩 Q6=A，不能继续匹配。

所以 k = 1 上最远到 (6,5)。

4.3 插入：跳到 k = -1​

插入 Q5=C，相当于 Q 前进、R 不动：

从 (4,4) 到 (4,5)，对角线：

k=4−5=−1k = 4 - 5 = -1

然后继续匹配：

R5=A  vs  Q6=A   匹配

到达 (5,6)。此时 Q 结束，R 还剩 R6=C，不能继续匹配。

所以 k = -1 上最远到 (5,6)。

4.4 s = 1 的波前​

分数 s对角线 k最远位置 (i,j)
11(6,5)
10(5,5)
1-1(5,6)

终点 (6,6) 还没到。

5. 分数 s = 2：允许两次编辑，到达终点​

从 s = 1 的波前继续。

5.1 从 k = 1 的 (6,5) 出发​

插入 Q6=A：

从 (6,5) 到 (6,6)，对角线：

k=6−6=0k = 6 - 6 = 0

到达终点 (6,6)，分数从 1 增加到 2。

5.2 从 k = -1 的 (5,6) 出发​

删除 R6=C：

从 (5,6) 到 (6,6)，对角线：

k=6−6=0k = 6 - 6 = 0

同样到达终点 (6,6)，分数从 1 增加到 2。

所以第一次到达终点是在 分数 s = 2。

因此最优 Levenshtein 距离：

d=2d = 2

6. 波前汇总​

分数 s对角线 k最远位置 (i,j)
00(4,4)
11(6,5)
10(5,5)
1-1(5,6)
20(6,6)

可以看到，WFA 不是填满整个 6×6 矩阵，而是只计算到分数 2 的波前。

7. 回溯：找回具体比对​

我们选一条路径：

从终点 (6,6) 回溯。

在 s = 2，到达 (6,6) 的操作是从 s = 1 的 k = 1 的 (6,5) 插入 Q6=A。

所以指针指向：

s=1, k=1, (6,5)

在 s = 1，(6,5) 是怎么来的？

是从 s = 0 的 (4,4) 删除 R5=A 到 (5,4)，然后匹配 R6=C 和 Q5=C 延伸到 (6,5)。

所以指针指向：

s=0, k=0, (4,4)

在 s = 0，从 (4,4) 沿对角线反向匹配：

R4=T  vs  Q4=T
R3=G  vs  Q3=G
R2=C  vs  Q2=C
R1=A  vs  Q1=A

回到 (0,0)。

8. 最终比对​

把回溯路径展开：

R: A C G T A C -
Q: A C G T - C A

对应关系：

位置RQ操作
1AA匹配
2CC匹配
3GG匹配
4TT匹配
5A-删除
6CC匹配
7-A插入

编辑操作：

  • 删除 R5=A
  • 插入 Q6=A

总共 2 次编辑，Levenshtein 距离 = 2。

9. 这个例子说明了什么？​

  1. WFA 按分数 s = 0, 1, 2, … 一层层扩展波前。
  2. 匹配是免费的，所以会沿对角线一直跳。
  3. 错配、插入、删除各让分数 +1，并切换到相邻或同一条对角线。
  4. 每条对角线只保留最远位置。
  5. 第一次到达终点时的分数就是最优编辑距离。
  6. 三角回溯矩阵记录每个格子来自哪条对角线、哪个分数，从而可以重建比对。

这个 6 碱基例子虽然简单，但已经包含了 WFA 的核心机制：波前、对角线、匹配延伸、编辑操作、回溯。TALCO 后面利用的正是这些回溯路径在 tile 边界上的汇聚性质。