序列比对：Needleman-Wunsch 算法

 

序列比对：Needleman-Wunsch 算法

在比对两个序列时，不仅要考虑完全匹配的字符，还要考虑一个序列中的空格或间隙（或者，相反地，要考虑另一个序列中的插入部分）和不匹配，在序列比对中，需要找到最优的比对（最优比对大致是指要将匹配的数量最大化，将空格和不匹配的数量最小化）。如果要更正式些，您可以确定一个分数，为匹配的字符添加分数、为空格和不匹配的字符减去分数。

全局和局部序列比对:

全局序列比对尝试找到两个完整的序列 S1 和 S2 之间的最佳比对。如S1=GCCCTAGCG S2=GCGCAATG 如果设定每个匹配字符为1分，每个空格为-2分，每个不匹配为-1分，则下面的比对就是全局最优比对：S1'=GCCCTAGCG S2'=GCGC_AATG,连字符“_”代表空格。在 S2' 中有五个匹配字符，一个空格（或者反过来说，在 S1' 中有一个插入项），有三个不匹配字符。这样得到的分数是 (5×1) + (1×-2) + (3×-1) = 0，这是能够实现的最佳结果。

 

Needleman-Wunsch 算法：

该算法用来计算全局比对。它的思路与 LCS 算法相似。这个算法也使用二维表格，一个序列沿顶部展开，一个序列沿左侧展开。而且也能通过以下三个途径到达每个单元格：1.来自上面的单元格，代表将左侧的字符与空格比对。2.来自左侧的单元格，代表将上面的字符与空格比对。3.来自左上侧的单元格，代表与左侧和上面的字符比对（可能匹配也可能不匹配）。该单元格的值来自于一下3个中的最大值：（1）上方的值-2 （2）左边的值-2 （3）如果该单元格所在的行于所在的列对应的字符相等，则为左上值加1，否则为左上值-1。

 

 

 

#coding=utf_8_decode def Needleman_Wunsch(str1,str2): #字符串为空则返回 if str1=='' or str2=='': return '' #字符串长度 m=len(str1) n=len(str2) #初始化 lcs=[[i*(-2)] for i in range(0,m+1)] lcs[0]=[j*(-2) for j in range(0,n+1)] # for i in range(m): for j in range(n): lcs[i+1].append( max( lcs[i][j]+(1 if str1[i] == str2[j] else -1), lcs[i][j+1]-2, lcs[i+1][j]-2, ) ) #for i in range(m+1): #print lcs[i] i=m-1 j=n-1 common_substr1 = u'' common_substr2 = u'' common_substr1 = u"%s%s" % (str1[i], common_substr1) common_substr2 = u"%s%s" % (str2[j], common_substr2) #回溯 while True: if i == 0 and j == 0: break if str1[i] == str2[j]: if lcs[i-1][j-1]+1>lcs[i-1][j]-2 and lcs[i-1][j-1]+1>lcs[i][j-1]-2: i = i - 1 j = j -1 common_substr1 = u"%s%s" % (str1[i], common_substr1) common_substr2 = u"%s%s" % (str2[j], common_substr2) else: if lcs[i][j+1] > lcs[i+1][j]: i = i-1 common_substr1 = u"%s%s" % (str1[i], common_substr1) common_substr2 = u"%s%s" % ('-', common_substr2) else: j = j-1 common_substr1 = u"%s%s" % ('-', common_substr1) common_substr2 = u"%s%s" % (str2[j], common_substr2) else: if lcs[i-1][j-1]+1>lcs[i-1][j]-2 and lcs[i-1][j-1]+1>lcs[i][j-1]-2: i = i - 1 j = j -1 common_substr1 = u"%s%s" % (str1[i], common_substr1) common_substr2 = u"%s%s" % (str2[j], common_substr2) else: if lcs[i][j+1] > lcs[i+1][j]: i = i-1 common_substr1 = u"%s%s" % (str1[i], common_substr1) common_substr2 = u"%s%s" % ('-', common_substr2) else: j = j-1 common_substr1 = u"%s%s" % ('-', common_substr1) common_substr2 = u"%s%s" % (str2[j], common_substr2) print common_substr1 print common_substr2 Needleman_Wunsch("AGCACACA","ACACTA")