python - Needleman-Wunsch实现在cogent和skbio中给出了不同的对齐方式

skbio中的实现与pycogent中的实现相比，给出了一个奇怪的结果。

from cogent.align.algorithm import nw_align as nw_align_cogent
from skbio.alignment import global_pairwise_align_nucleotide as nw_align_scikit

seq_1 = 'ATCGATCGATCG'
seq_2 = 'ATCGATATCGATCG'

print "Sequences: "
print "     %s" % seq_1
print "     %s" % seq_2
print

alignment = nw_align_scikit(seq_1, seq_2)
al_1, al_2 = [alignment.get_seq(_id).__str__() for _id in alignment.ids()]

print "    nw alignment using scikit:"
print "        %s" % al_1
print "        %s" % al_2
print

al_1, al_2 = nw_align_cogent(seq_1, seq_2)

print "    nw alignment using cogent:"
print "        %s" % al_1
print "        %s" % al_2
print

输出如下：

nw alignment using scikit:
    ------ATCGATCGATCG
    ATCGATATCGATCG----

nw alignment using cogent:
    ATCGAT--CGATCG
    ATCGATATCGATCG

最佳答案

这是一个很好的问题，并且与scikit bio和PyCogent中的对齐方式的不同有关。默认情况下，在scikit bio中，终端间隙不会受到惩罚，因为这可能导致一些非常奇怪的排列。这个issue was briefly discussed here和图示here（参见笔记本的最后一个单元格）。
如果您希望获得与PyCogent中的解决方案更相似的解决方案，可以按如下方式将penalize_terminal_gaps=True传递到global_pairwise_align_nucleotide：

alignment = nw_align_scikit(seq_1, seq_2, penalize_terminal_gaps=True)
al_1, al_2 = [alignment.get_seq(_id).__str__() for _id in alignment.ids()]

print "    nw alignment using scikit:"
print "        %s" % al_1
print "        %s" % al_2

输出：

nw alignment using scikit:
        ATCG--ATCGATCG
        ATCGATATCGATCG

您会注意到对齐仍然与从PyCogent获得的对齐不同，但这是一个较小的实现差异：两个结果对齐具有相同的分数（差异在于--是与ATrepeat中的第一个AT对齐，还是与ATATrepeat中的第二个penalize_terminal_gaps=True对齐），并且两个实现在如何选择他们打破了那条领带。
如果你回到你发布的对齐方式（scikit bio的默认设置），你会注意到返回的对齐方式非常好——事实上，如果不惩罚终端间隙（根据定义，因为最佳的评分对齐方式就是它返回的内容），它就是最佳的评分对齐方式。不过，你说得对，这很奇怪，因为scikit bio在这种特定情况下返回的对齐可能不是最与生物学相关的对齐。如果你知道你的序列都是从同一个位置开始，在同一个位置结束，你就可以通过。然而，您的例子是一个很好的例子，在大多数情况下，对于真实的序列，我认为与生物学最相关的对齐将返回默认参数。

关于python - Needleman-Wunsch实现在cogent和skbio中给出了不同的对齐方式，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/25332841/