使用pairwise2的成对比对
Please note that Bio.pairwise2 was deprecated in Release 1.80. 作为替代方案，请考虑使用 Bio.Align.PairwiseAligner （在第章中描述 双序列比对 ).
Bio.pairwise2 contains essentially the same algorithms as water (local) and needle (global) from the EMBOSS suite（见上文）并返回相同的结果。的 pairwise2 模块最近在速度和内存消耗方面进行了一些优化（Biopython版本>1.67），因此对于短序列（全局比对：~2000个残基，局部比对~600个残基），使用起来更快（或同样快） pairwise2 而不是打电话给TMF ' water 或 needle 通过命令行工具。
假设您想在上面相同的两个血红蛋白序列之间进行全局成对比对 (HBA_HUMAN , HBB_HUMAN ）存储在 alpha.faa 和 beta.faa :
>>> from Bio import pairwise2
>>> from Bio import SeqIO
>>> seq1 = SeqIO.read("alpha.faa", "fasta")
>>> seq2 = SeqIO.read("beta.faa", "fasta")
>>> alignments = pairwise2.align.globalxx(seq1.seq, seq2.seq)
如您所见，我们将对齐函数称为 align.globalxx .棘手的部分是函数名称的最后两个字母（这里： xx ），用于解码匹配（和不匹配）和差距的分数和处罚。第一个字母解码比赛比分，例如 x 意味着匹配计算1，而不匹配则没有成本。与 m 可以定义匹配或不匹配的常规值（有关完整详细信息，请参见 Bio.pairwise2 ).第二封信解码了缺口的成本; x 意味着根本没有缺口成本， s 可以对打开和扩大差距进行不同的处罚。所以， globalxx 意味着仅计算两个序列之间的匹配。
我们的可变 alignments 现在包含对给定条件具有相同最佳分数的比对列表（至少一个）。在我们的示例中，这是80个不同的对齐，得分为72 (Bio.pairwise2 将返回多达1000个对齐）。看看这些路线之一：
>>> len(alignments)
80
>>> print(alignments[0])
Alignment(seqA='MV-LSPADKTNV---K-A--A-WGKVGAHAG...YR-', seqB='MVHL-----T--PEEKSAVTALWGKV----...Y-H', score=72.0, start=0, end=217)
每个比对都是一个命名的多元组，由两个比对序列、分数、比对的开始和结束位置组成（在全局比对中，开始始终是0，结束始终是比对的长度）。 Bio.pairwise2 具有如下功能 format_alignment 想要更好的收件箱：
>>> print(pairwise2.format_alignment(*alignments[0]))
MV-LSPADKTNV---K-A--A-WGKVGAHAG---EY-GA-EALE-RMFLSF----PTTK-TY--F...YR-
|| | | | | | |||| | | ||| | | | | |...|
MVHL-----T--PEEKSAVTALWGKV-----NVDE-VG-GEAL-GR--L--LVVYP---WT-QRF...Y-H
Score=72
自Biopython 1.77以来，所需的参数可以随关键字一起提供。最后一个例子现在也可以写为：
>>> alignments = pairwise2.align.globalxx(sequenceA=seq1.seq, sequenceB=seq2.seq)
更好的对齐通常是通过惩罚差距来实现的：打开差距的成本更高，扩大现有差距的成本更低。对于氨基酸序列，匹配分数通常编码在如下矩阵中 PAM 或 BLOSUM .因此，通过使用BLOSUM62矩阵以及空位开放罚分10和空位延伸罚分0.5（使用 globalds ):
>>> from Bio import pairwise2
>>> from Bio import SeqIO
>>> from Bio.Align import substitution_matrices
>>> blosum62 = substitution_matrices.load("BLOSUM62")
>>> seq1 = SeqIO.read("alpha.faa", "fasta")
>>> seq2 = SeqIO.read("beta.faa", "fasta")
>>> alignments = pairwise2.align.globalds(seq1.seq, seq2.seq, blosum62, -10, -0.5)
>>> len(alignments)
2
>>> print(pairwise2.format_alignment(*alignments[0]))
MV-LSPADKTNVKAAWGKVGAHAGEYGAEALERMFLSFPTTKTY...KYR
|| |.|..|..|.|.|||| ......|............|.......||.
MVHLTPEEKSAVTALWGKV-NVDEVGGEALGRLLVVYPWTQRFF...KYH
Score=292.5
这种比对具有与我们之前使用相同序列和相同参数使用CLARSS针获得的相同分数。
局部对齐的调用与该函数类似 align.localXX ，其中XX再次代表匹配和间隔函数的两个字母代码：
>>> from Bio import pairwise2
>>> from Bio.Align import substitution_matrices
>>> blosum62 = substitution_matrices.load("BLOSUM62")
>>> alignments = pairwise2.align.localds("LSPADKTNVKAA", "PEEKSAV", blosum62, -10, -1)
>>> print(pairwise2.format_alignment(*alignments[0]))
3 PADKTNV
|..|..|
1 PEEKSAV
Score=16
在最近的Biopython版本中， format_alignment 将仅打印局部对齐的对齐部分（以及以1为基础的符号形式的开始位置，如上面的示例所示）。如果您也对序列的非对齐部分感兴趣，请使用关键字参数 full_sequences=True :
>>> from Bio import pairwise2
>>> from Bio.Align import substitution_matrices
>>> blosum62 = substitution_matrices.load("BLOSUM62")
>>> alignments = pairwise2.align.localds("LSPADKTNVKAA", "PEEKSAV", blosum62, -10, -1)
>>> print(pairwise2.format_alignment(*alignments[0], full_sequences=True))
LSPADKTNVKAA
|..|..|
--PEEKSAV---
Score=16
请注意，正如Smith & Waterman所定义的那样，局部对齐必须具有正值（>0）。因此， pairwise2 如果没有获得>0的分数，则可能不返回任何对齐。此外， pairwise2 不会报告由于在任一站点上添加零分扩展而导致的对齐。在下一个示例中，丝氨酸/天冬氨酸（S/D）和赖氨酸/天冬酰胺（K/N）对都具有0的匹配分数。如你所见，对齐的部分没有被延伸：
>>> from Bio import pairwise2
>>> from Bio.Align import substitution_matrices
>>> blosum62 = substitution_matrices.load("BLOSUM62")
>>> alignments = pairwise2.align.localds("LSSPADKTNVKKAA", "DDPEEKSAVNN", blosum62, -10, -1)
>>> print(pairwise2.format_alignment(*alignments[0]))
4 PADKTNV
|..|..|
3 PEEKSAV
Score=16
匹配代码不是提供完整的匹配/不匹配矩阵， m 允许轻松定义一般匹配/不匹配值。下一个示例使用5/-4的匹配/不匹配分数和2/0.5的差距罚分（开放/扩展），使用 localms :
>>> alignments = pairwise2.align.localms("AGAACT", "GAC", 5, -4, -2, -0.5)
>>> print(pairwise2.format_alignment(*alignments[0]))
2 GAAC
| ||
1 G-AC
Score=13
的一个有用的关键字参数 Bio.pairwise2.align 函数 score_only .如果设置为 True 它将仅返回最佳对齐的分数，但是时间要短得多。它还允许在内存错误出现之前对更长的序列进行比对。另一个有用的关键字参数是 one_alignment_only=True 这也将导致一些速度增加。
不幸的是， Bio.pairwise2 （目前）不适用于Biopython的多序列比对对象。然而，该模块具有一些有趣的高级功能：您可以定义自己的匹配和间隙函数（有兴趣测试仿射对数间隙成本？），两个序列的缺口罚分和末端缺口罚分可能不同，序列可以作为列表提供（如果您有由多个字符编码的残基，这很有用）等。这些特征很难（如果有的话）用其他比对工具实现。有关更多详细信息，请参阅该模块的API文档 Bio.pairwise2 .