利用Needleman–Wunsch算法进行DNA序列全局比对

本文详细介绍了如何使用Needleman-Wunsch算法进行DNA序列的全局比对,并提供了一个具体的Python实现示例。该算法适用于生物信息学中序列比对的基本任务,能够帮助理解不同序列间的相似性和差异。

生物信息学原理作业第二弹:利用Needleman–Wunsch算法进行DNA序列全局比对。

具体原理:https://en.wikipedia.org/wiki/Needleman%E2%80%93Wunsch_algorithm

利用Needleman–Wunsch算法进行DNA序列全局比对

转载请保留出处!

贴上python代码:

 1 # -*- coding: utf-8 -*-
 2 """
 3 Created on Sat Nov 25 18:20:01 2017
 4 
 5 @author: zxzhu
 6 后需修改:
 7 1.加命令行参数
 8 2.给出多种比对结果
 9 """
10 
11 import numpy as np
12 import pandas as pd
13 sequence1 = 'AACGTACTCA'
14 sequence2 = 'TCGTACTCA'
15 s1 = ''
16 s2 = ''
17 gap = -4
18 score_matrix = pd.read_excel('score.xlsx')        #score matrix
19 best_matrix = np.empty(shape= (len(sequence2)+1,len(sequence1)+1),dtype = int)
20 
21 def get_match_score(s1,s2):
22     score = score_matrix[s1][s2]
23     return score
24 
25 for i in range(len(sequence2)+1):
26     for j in range(len(sequence1)+1):
27         if i == 0:
28             best_matrix[i][j] = gap * j
29         
30         elif j == 0:
31             best_matrix[i][j] = gap *i
32         else:
33             match = get_match_score(sequence2[i-1],sequence1[j-1])
34             gap1_score = best_matrix[i-1][j]+gap
35             gap2_score = best_matrix[i][j-1]+gap
36             match_score = best_matrix[i-1][j-1]+match
37             best_matrix[i][j] = max(gap1_score,gap2_score,match_score)
38 print(best_matrix)
39 i,j = len(sequence2),len(sequence1)
40 while(i>0 or j>0):
41     match = get_match_score(sequence2[i-1],sequence1[j-1])
42     if i>0 and j>0 and best_matrix[i][j] == best_matrix[i-1][j-1]+match:
43         s1 += sequence1[j-1]
44         s2 += sequence2[i-1]
45         i-=1;j-=1
46     elif i>0 and best_matrix[i,j] == best_matrix[i-1,j]+gap:
47         s1+='-'
48         s2+=sequence2[i-1]
49         i-=1
50     else:
51         s1+=sequence1[j-1]
52         s2+='-'
53         j-=1
54 print(s1[::-1]+'\n'+s2[::-1])

后面会加入命令行。

多种结果这里只取了一种,这个问题有待解决。

如果有其他的方法我会及时添加。

转载于:https://www.cnblogs.com/zxzhu/p/7903706.html

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值