BWA-mem Smith-Waterman 算法

时间:2024-10-14 21:25:11

Smith-Waterman 算法简介

Smith-Waterman 算法 是一种用于局部序列比对的动态规划算法,它可以在两个序列中找到最优的局部比对,允许处理序列中的插入缺失(indels)错配。该算法最早由 Temple F. Smith 和 Michael S. Waterman 于 1981 年提出,因此得名 Smith-Waterman

核心思想

Smith-Waterman 算法通过构建一个得分矩阵,从两个序列中找到得分最高的局部比对,即找到序列中最相似的子序列,哪怕这些序列之间存在错配或 indels(插入和缺失)。

应用场景

该算法通常用于:

  • 局部比对:寻找两个序列中最相似的片段。
  • 比对短序列和长参考序列:如在基因组比对中,BWA-MEM 就使用 Smith-Waterman 算法扩展比对。
  • 序列相似性分析:用来识别基因或蛋白质序列中的相似区域。

与全局比对(如 Needleman-Wunsch 算法)的区别

  • 局部比对:Smith-Waterman 专注于找到两个序列中最相似的一部分,而不是从头到尾的全局比对。这特别有用,当你关心的是局部的相似性,而不是整个序列的相似性。
  • 允许 indels 和错配:在比对过程中,Smith-Waterman 允许插入、缺失和错配,但会根据这些事件给出相应的惩罚分数。

算法步骤

  1. 初始化得分矩阵:构建一个二维得分矩阵,序列 A 作为行,序列 B 作为列。矩阵的每个元素表示当前子序列的比对得分。
  2. 填充得分矩阵:根据相似性(匹配得分)和惩罚(错配或 indels)来填充得分矩阵。计算方式如下:
    • 匹配:如果两个字符匹配,增加匹配得分。
    • 错配:如果两个字符不匹配,扣除错配罚分。
    • 插入或缺失:如果在一个序列中存在插入或缺失,扣除 indel 罚分。
  3. 寻找最大得分:找到得分矩阵中的最大值,这表示局部比对的最优解。
  4. 回溯路径:从最大值开始,沿着得分矩阵回溯,找到局部比对的路径(即序列片段)。

示例:

假设有两个序列:

  • 序列 A:GACTTAC
  • 序列 B:CGTGAATTCAT

构建一个得分矩阵,按照上述步骤进行比对。Smith-Waterman 算法将会找到最相似的局部子序列,并返回局部比对结果。

优点:

  • 局部比对:能够准确找到最相似的片段,即使序列中有较长的不匹配部分。
  • 灵活处理 indels:适合处理带有插入和缺失的序列。

缺点:

  • 时间复杂度较高:Smith-Waterman 的时间复杂度为 O(m*n),其中 m 和 n 分别是两个序列的长度,因此对长序列的比对效率较低。

总结

Smith-Waterman 算法是一种强大的局部序列比对工具,它通过动态规划方法找到最优的局部比对区域,允许插入、缺失和错配。这使得它在基因比对、蛋白质序列分析等领域非常有用。

让我们通过这个矩阵来演示 Smith-Waterman 算法的整个比对过程。我们会逐步填写这个矩阵,并通过这个例子理解算法的流程。

输入序列
序列 A:GACTTAC
序列 B:CGTGAATTCAT
我们将通过填充比对网格,找到它们之间的最优局部比对。

步骤 1:初始化矩阵
首先,我们创建一个大小为 8 x 12 的矩阵(包括初始行和列,行列代表序列的字符)。矩阵的每个格子将存储两个字符之间的得分。

在 Smith-Waterman 算法中,矩阵的第一行和第一列初始化为 0。这是因为比对可以从任何位置开始,不必从头比到尾。
整个 Smith-Waterman 算法的比对过程。接下来,我们将继续基于序列 A (GACTTAC) 和 B (CGTGAATTCAT) 填充矩阵,找出最优的局部比对。

1. 初始化矩阵

首先,初始化矩阵第一行和第一列为 0:

      C   G   T   G   A   A   T   T   C   A   T
    +------------------------------------------
G   0   0   0   0   0   0   0   0   0   0   0
A   0
C   0
T   0
T   0
A   0
C   0

2. 定义打分规则

  • 匹配得分:+1
  • 错配罚分:-1
  • 插入/缺失(indel)罚分:-1

3. 填充矩阵

我们继续按照之前的示例填充剩下的部分。

第 1 行(序列 G 与 C, G, T, G, A, A, T, T, C, A, T)
  • (G, C):错配,得分 0。
  • (G, G):匹配,得分 1。
  • (G, T):错配,得分 0。
  • (G, G):匹配,得分 1。
  • (G, A):错配,得分 0。
  • (G, A):错配,得分 0。
  • (G, T):错配,得分 0。
  • (G, T):错配,得分 0。
  • (G, C):错配,得分 0。
  • (G, A):错配,得分 0。
  • (G, T):错配,得分 0。

因此,第 1 行填充如下:

      C   G   T   G   A   A   T   T   C   A   T
    +------------------------------------------
G   0   0   1   0   1   0   0   0   0   0   0
A   0
C   0
T   0
T   0
A   0
C   0
第 2 行(序列 A 与 C, G, T, G, A, A, T, T, C, A, T)
  • (A, C):错配,得分 0。
  • (A, G):错配,得分 0。
  • (A, T):错配,得分 0。
  • (A, G):错配,得分 0。
  • (A, A):匹配,得分 1。
  • (A, A):匹配,得分 1。
  • (A, T):错配,得分 0。
  • (A, T):错配,得分 0。
  • (A, C):错配,得分 0。
  • (A, A):匹配,得分 1。
  • (A, T):错配,得分 0。

因此,第 2 行填充如下:

      C   G   T   G   A   A   T   T   C   A   T
    +------------------------------------------
G   0   0   1   0   1   0   0   0   0   0   0
A   0   0   0   0   0   1   1   0   0   0   1
C   0
T   0
T   0
A   0
C   0
第 3 行(序列 C 与 C, G, T, G, A, A, T, T, C, A, T)
  • (C, C):匹配,得分 1。
  • (C, G):错配,得分 0。
  • (C, T):错配,得分 0。
  • (C, G):错配,得分 0。
  • (C, A):错配,得分 0。
  • (C, A):错配,得分 0。
  • (C, T):错配,得分 0。
  • (C, T):错配,得分 0。
  • (C, C):匹配,得分 1。
  • (C, A):错配,得分 0。
  • (C, T):错配,得分 0。

因此,第 3 行填充如下:

      C   G   T   G   A   A   T   T   C   A   T
    +------------------------------------------
G   0   0   1   0   1   0   0   0   0   0   0
A   0   0   0   0   0   1   1   0   0   0   1
C   0   1   0   0   0   0   0   0   1   0   0
T   0
T   0
A   0
C   0
继续填充剩余行

通过类似的步骤,我们可以逐步完成整个矩阵的填充。

4. 回溯路径

通过矩阵得分最高的区域,我们可以回溯找到局部最优的比对区域。

这就是 Smith-Waterman 算法的完整流程。你对这个过程有疑问或需要进一步的详细说明吗?