Global Alignment(全局比对)--从算法(Needleman-Wunsch)到python实现
时间: 2024-01-03 13:03:16 浏览: 27
全局比对是一种常用的序列比对方法,它可以比对两个序列的整个长度,并给出它们的相似程度。其中 Needleman-Wunsch 是一种经典的全局比对算法。
Needleman-Wunsch 算法的本质是动态规划,它将全局比对的问题拆分成若干个子问题,并通过计算每个子问题的得分来求解最终的比对结果。具体实现步骤如下:
1. 初始化得分矩阵(score matrix)并计算第一行和第一列的得分;
2. 逐行逐列计算得分矩阵中的每个元素的得分,根据得分矩阵中每个元素的左上、上、左三个元素的得分值,计算当前元素的得分;
3. 回溯得分矩阵,得到最优比对结果。
以下是 Needleman-Wunsch 算法的 Python 实现:
```python
# Needleman-Wunsch 算法实现
import numpy as np
def needleman_wunsch(seq1, seq2, match_score=1, mismatch_score=-1, gap_score=-1):
# 初始化得分矩阵
m, n = len(seq1), len(seq2)
score_matrix = np.zeros((m+1, n+1))
for i in range(1, m+1):
score_matrix[i,0] = score_matrix[i-1,0] + gap_score
for j in range(1, n+1):
score_matrix[0,j] = score_matrix[0,j-1] + gap_score
# 填充得分矩阵
for i in range(1, m+1):
for j in range(1, n+1):
score1 = score_matrix[i-1,j-1] + (match_score if seq1[i-1] == seq2[j-1] else mismatch_score)
score2 = score_matrix[i-1,j] + gap_score
score3 = score_matrix[i,j-1] + gap_score
score_matrix[i,j] = max(score1, score2, score3)
# 回溯得分矩阵,得到最优比对结果
align1, align2 = '', ''
i, j = m, n
while i > 0 and j > 0:
score, score1, score2, score3 = score_matrix[i,j], score_matrix[i-1,j-1], score_matrix[i-1,j], score_matrix[i,j-1]
if score == score1 + (match_score if seq1[i-1] == seq2[j-1] else mismatch_score):
align1 = seq1[i-1] + align1
align2 = seq2[j-1] + align2
i, j = i-1, j-1
elif score == score2 + gap_score:
align1 = seq1[i-1] + align1
align2 = '-' + align2
i -= 1
else:
align1 = '-' + align1
align2 = seq2[j-1] + align2
j -= 1
while i > 0:
align1 = seq1[i-1] + align1
align2 = '-' + align2
i -= 1
while j > 0:
align1 = '-' + align1
align2 = seq2[j-1] + align2
j -= 1
return align1, align2
```
在上面的代码中,我们使用了 `numpy` 库来创建得分矩阵,这使得代码更加简洁和高效。在实现中,我们还可以通过调整 `match_score`、`mismatch_score` 和 `gap_score` 的值来控制比对的结果。