Fork me on GitHub

序列比对(一) Needleman-Wunsch算法

引言

序列比对是生物信息学的经典问题。为了解决生物序列比对的问题,人们已经开发了很多算法,其中包含精确比对,以及为了优化比对速度而诞生的启发式比对算法等等。Needleman-Wunsch算法是一种全局比对算法,也就是会考虑两条序列的所有碱基的比对算法,如下:

  ATCGGGGGGGGG---GGGGGGCGCG
AT-GGGGGGG--CCC------CGCG

算法介绍

(1)建立得分矩阵

所谓得分矩阵,就是用于序列比对时计算比对得分,比对得分越高,就说明该比对越接近最优比对。根据比对的不同情况可以分为:完全匹配,错误匹配,插入空位三种情况。针对这三种不同的比对情况,我们可以根据得分矩阵进行打分处理。

下表是一个得分矩阵的例子:

0 A T C G
0 0 -2 -2 -2 -2
A -2 1 -1 -1 -1
T -2 -1 1 -1 -1
C -2 -1 -1 1 -1
G -2 -1 -1 -1 1

以上只是一个比较简单的例子,没有根据生物学序列的具体问题进行设计。比如:根据不同碱基的错误匹配,可以给出不同的打分,当然这个针对氨基酸序列的比对更有意义;插入空位的情况在生物学序列中的出现的概率要比错误匹配的概率低很多,那么根据这个生物学特性,可以将插入空位的罚分设的更低。

(2)初始化相似性分数矩阵

假设有两条待比对的序列,根据得分矩阵计算每一个碱基之间的相似性分数,从而构建相似性得分矩阵,以(D序列)ACGCG和(Q序列)CCGACC两条序列为例:

- A C G C G
- 0 -2 -2 -2 -2 -2
C -2
C -2
G -2
A -2
C -2
C -2

根据上面的得分矩阵,可以得到碱基与空位的相似性分数为-2,所以把对应位置的相似性分数初始化为-2.

(3)根据下面的逻辑计算各个位置的相似性得分

为了计算某个位置的相似性分数,需要考虑三个方向,横向,纵向和对角线。三个方向,其实是针对某个比对位置的三种比对方式,选择分数最大的一种比对方式。比较好理解的是对角线方向,也就是直接匹配。横向和纵向都是插入空格。

对角线: V(i,j) = v(i-1,j-1) + F(D[i],Q[j])

横向: V(i,j) = v(i-1,j) + F(D[i],-)

纵向:V(i,j) = v(i,j-1) + F(-,Q[j])

选择以上三种计算方式得到的最大相似性分数作为该位置的相似性分数。以此方法填充得到完整的相似性得分矩阵。

(4)通过回溯法,找到最佳比对,并得到最高得分

从相似性得分矩阵的右下角向上回溯,每次取三个方向的最大值作为走向,找到最佳比对,输出比对。

以上是Needleman-Wunsch算法的整体思路,它主要思路是基于动态规划法。