实现watermanegert算法

Question

我正在尝试实现Waterman-Eggert算法来寻找次优的局部序列比对，但是我很难理解如何“分解”不同的对齐组。我有基本的史密斯-沃特曼算法运行良好。在

一个简单的测试，将以下序列与自身对齐：

'HEAGHEAGHEAG'
'HEAGHEAGHEAG'

生成fMatrix，如下所示：

^{pr2}$

为了找到次优路线，例如

'HEAGHEAGHEAG    '
'    HEAGHEAGHEAG'

您必须首先移除最佳对齐（即沿主对角线）并重新计算fMatrix；这称为“去集总”，其中路线的“束”定义为路径相交/共享一对或多对对齐残基的任何路线。除了fMatrix之外，还有一个次级矩阵，它包含关于fMatrix构建方向的信息。在

构建fMatrix和回溯矩阵的代码片段如下：

# Generates fMatrix.
for i in range(1, length):
    for j in range(1, length):
        matchScore = fMatrix[i-1][j-1] + simMatrixDict[seq[i-1]+seq[j-1]]
        insScore = fMatrix[i][j-1] + gap
        delScore = fMatrix[i-1][j] + gap
        fMatrix[i][j] = max(0, matchScore, insScore, delScore)

        # Generates matrix for backtracking.
        if fMatrix[i][j] == matchScore:
            backMatrix[i][j] = 2      
        elif fMatrix[i][j] == insScore:
            backMatrix[i][j] = 3        # INSERTION in seq - Horizontal
        elif fMatrix[i][j] == delScore:
            backMatrix[i][j] = 1        # DELETION in seq - Vertical

        if fMatrix[i][j] >= backtrackStart: 
            backtrackStart = fMatrix[i][j]
            endCoords = i, j                
return fMatrix, backMatrix, endCoords

为了删除这个最佳对齐方式，我尝试使用这个回溯矩阵来回溯fMatrix（按照最初的Smith Waterman算法），并在我前进的过程中设置fMatrix[i][j] = 0，但这并没有删除整个块，只删除了该块中的精确对齐。在

对于一些背景信息，Smith Waterman算法的Wikipedia页面解释了fMatrix是如何构造的，here上有关于回溯如何工作的解释。对Waterman-Eggert算法作了粗略的解释here。在

谢谢。在

Answer 1

好吧。这里有一些代码可以让你随心所欲。我使用了pretty print库（pprint），这样输出看起来很不错。（如果矩阵中的数字是个位数，那么它看起来最漂亮，但是如果有多个数字，对齐就会变得有些混乱。）

它是如何工作的？

因为你只需要改变主对角线上的数字和上下对角线上的数字，我们只需要一个循环。matrix[i][i]总是在主对角线上，因为它向下是i行，而i列。matrix[i][i-1]始终是下邻对角线，因为它向下是i行，而{}列。matrix[i-1][i]总是上相邻的对角线，因为它向下是i-1行，而{}行。在

#!/usr/bin/python
import pprint

matrix = [
    [  0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,],
    [  0,   8,   0,   0,   0,   8,   0,   0,   0,   8,   0,   0,   0,],
    [  0,   0,  13,   0,   0,   0,  13,   0,   0,   0,  13,   0,   0,],
    [  0,   0,   0,  17,   0,   0,   0,  17,   0,   0,   0,  17,   0,],
    [  0,   0,   0,   0,  23,   0,   0,   0,  23,   0,   0,   0,  23,],
    [  0,   8,   0,   0,   0,  31,   0,   0,   0,  31,   0,   0,   0,],
    [  0,   0,  13,   0,   0,   0,  36,   0,   0,   0,  36,   0,   0,],
    [  0,   0,   0,  17,   0,   0,   0,  40,   0,   0,   0,  40,   0,],
    [  0,   0,   0,   0,  23,   0,   0,   0,  46,   0,   0,   0,  46,],
    [  0,   8,   0,   0,   0,  31,   0,   0,   0,  54,   4,   0,   0,],
    [  0,   0,  13,   0,   0,   0,  36,   0,   0,   4,  59,   9,   0,],
    [  0,   0,   0,  17,   0,   0,   0,  40,   0,   0,   9,  63,  13,],
    [  0,   0,   0,   0,  23,   0,   0,   0,  46,   0,   0,  13,  69,]]

print "Original Matrix"
pprint.pprint(matrix)
print

for i in range(len(matrix)):
    matrix[i][i] = 0
    if (i > 0) and (i < (len(matrix))):
        matrix[i][i-1] = 0
        matrix[i-1][i] = 0

print "New Matrix"
pprint.pprint(matrix)

输出：

^{pr2}$