python实现动态规划全局比对

python实现动态规划全局比对

• • 代码

代码

import numpy as np

def huisu(m,n):
    num=matrix[m][n]-matrix[m-1][n-1]
    if (num==match and seq_A[n-1] == seq_B[m-1]) or (num ==mis_match and seq_A[n-1] != seq_B[m-1]):
        matched_A.insert(0,seq_A[n-1])
        matched_B.insert(0,seq_B[n-1])
        m-=1
        n-=1
    elif (matrix[m][n]-matrix[m][n-1]) == gap_symbol or m==0:
        matched_A.insert(0,seq_A[n-1])
        matched_B.insert(0,'-')
        n-=1
    else :
        matched_A.insert(0,'-')
        matched_B.insert(0,seq_B[n-1])
        m-=1
    return(m,n)

seq_A=input("SequenceA:")
seq_B=input("SequenceB:")
match=int(input("Match:"))
mis_match=int(input("Mismatch:"))
gap_symbol=int(input("Gap_symbol:"))

length_A=len(seq_A)
length_B=len(seq_B)
matrix=np.zeros((length_B+1,length_A+1))
for i in range(1,length_A+1):
    matrix[0][i]=matrix[0][i-1]+gap_symbol
for i in range(1,length_B+1):
    matrix[i][0]=matrix[i-1][0]+gap_symbol

for m in range(1,length_B+1):
    for n in range(1,length_A+1):
        upd=matrix[m-1][n]+gap_symbol
        leftd=matrix[m][n-1]+gap_symbol
        if seq_B[m-1] == seq_A[n-1] :
            upleftd=matrix[m-1][n-1]+match
        else:
            upleftd=matrix[m-1][n-1]+mis_match
        matrix[m][n]=max(upd,leftd,upleftd)
#print(matrix)
print("\t0\t",end="")
for i in seq_A:
    print(i,end="\t")
print("\n")
for m in range(0,length_B+1):
    if m > 0:
        print(seq_B[m-1],end="\t")
    else:
        print("0",end="\t")
    for n in range(0,length_A+1):
        print(matrix[m][n],end="\t")
    print("\n")

m=length_B
n=length_A
matched_A=[]
matched_B=[]
while(m>0 or n>0):
    (m,n)=huisu(m,n)
print("\n")
for i in matched_A:
    print(i,end=" ")
print("\n")
for i in matched_B:
    print(i,end=" ")
print("\n")

声明：写得不好，别喷

本文地址：https://blog.csdn.net/weixin_45214599/article/details/110481212