重叠记分矩阵生物塞顿

Question

我有一个包含DNA序列和序列名称的FASTA文件，我需要制作一个重叠分数矩阵。我在Biopython中找到了模块pairwise2，它似乎做得很好。除了我的序列已经对齐，当我使用pairwise2时，它再次尝试对齐序列，这需要很长的时间，而且每次对齐都会得到相同的重叠分数。所以我的问题是，在不尝试重新对齐序列的情况下，如何获得重叠分数？ 以下是我目前所掌握的情况：

from Bio.Alphabet import IUPAC
from Bio import SeqIO
from Bio import pairwise2

fasta_file = SeqIO.parse('unambiguous.fasta', 'fasta', alphabet=IUPAC.ambiguous_dna)

all_seq = []
for seq_record in fasta_file:
    all_seq += [str(seq_record.seq)]

compare = pairwise2.align.globalms(all_seq[0], all_seq[1], 2, -1, -1, 0)
print(compare)

我只使用了FASTA文件中的第一个和第二个序列作为试用。正如你在脚本中看到的，匹配应该得到2分，不匹配和差距-1。当两个序列在同一位置上有间隔时，0应该是奖励。我知道把0放在第四位不会给我想要的结果，但我还没有解决这个问题的办法。在这一点上，对齐问题似乎更大。 那么，有没有人有过pairwise2或其他python/biopython模块的经验，可以让我得到重叠分数吗？在

Answer 1

据我所知，unambiguous.fasta包含对齐的遗传序列。您可以使用符合您需要的评分功能对其进行评分：

from itertools import starmap, combinations


def score(seq1, seq2):
    def score_(a, b):
        return (0 if a == b == "-" # both are gaps
                else -1 if a != b  # mismatch or gap
                else 2)            # match

    return sum(starmap(score_, zip(seq1, seq2)))

您可能需要修改它以忽略基础不明确的位置，就像人们通常所做的那样。这里有一种比较所有序列的简便方法：

^{pr2}$

一旦执行，scores（注意它是一个惰性迭代器）将生成分数对矩阵的平坦上三角。score应该可以很快地工作，尽管您可能希望在使用Cython或Numba时重新实现，以防有数千个序列（即要计算数百万个比较）。在

在python2.x上，您可能希望将zip替换为izip。在