Biopython 序列比对
序列比对 是以特定顺序排列两个或多个序列(DNA、RNA 或蛋白质序列)以识别它们之间的相似区域的过程。
识别相似区域使我们能够推断出很多信息，例如物种之间保守的特征、不同物种在遗传上的接近程度、物种如何进化等。Biopython 为序列比对提供了广泛的支持。
让我们在本章中学习 Biopython 提供的一些重要特性：
解析序列对齐
Biopython 提供了一个模块 Bio.AlignIO 来读取和写入序列比对。在生物信息学中，有许多格式可用于指定类似于早期学习的序列数据的序列比对数据。 Bio.AlignIO 提供类似于 Bio.SeqIO 的 API，不同之处在于 Bio.SeqIO 处理序列数据而 Bio.AlignIO 处理序列比对数据。
在开始学习之前，让我们从网上下载一个样本序列比对文件。
要下载示例文件，请按照以下步骤操作：
步骤 1 ： 打开你喜欢的浏览器，进入 http://pfam.xfam.org/family/browse 网站。它将按字母顺序显示所有 Pfam 系列。
步骤 2 ： 选择任何一个具有较少种子价值的家庭。它包含最少的数据，使我们能够轻松地进行对齐。在这里，我们选择/单击了 PF18225 并打开转到 http://pfam.xfam.org/family/PF18225 并显示有关它的完整详细信息，包括序列比对。
步骤 3 ： 进入比对部分，下载Stockholm格式的序列比对文件(PF18225_seed.txt)。
让我们尝试使用 Bio.AlignIO 读取下载的序列比对文件，如下所示：
导入 Bio.AlignIO 模块
>>> from Bio import AlignIO
使用读取方法读取对齐。 read 方法用于读取给定文件中可用的单个对齐数据。如果给定的文件包含许多对齐，我们可以使用 parse 方法。 parse 方法返回可迭代的对齐对象，类似于 Bio.SeqIO 模块中的 parse 方法。
>>> alignment = AlignIO.read(open("PF18225_seed.txt"), "stockholm")
打印对齐对象。
>>> print(alignment) SingleLetterAlphabet() alignment with 6 rows and 65 columns MQNTPAERLPAIIEKAKSKHDINVWLLDRQGRDLLEQRVPAKVA...EGP B7RZ31_9GAMM/59-123 AKQRGIAGLEEWLHRLDHSEAIPIFLIDEAGKDLLEREVPADIT...KKP A0A0C3NPG9_9PROT/58-119 ARRHGQEYFQQWLERQPKKVKEQVFAVDQFGRELLGRPLPEDMA...KKP A0A143HL37_9GAMM/57-121 TRRHGPESFRFWLERQPVEARDRIYAIDRSGAEILDRPIPRGMA...NKP A0A0X3UC67_9GAMM/57-121 AINRNTQQLTQDLRAMPNWSLRFVYIVDRNNQDLLKRPLPPGIM...NRK B3PFT7_CELJU/62-126 AVNATEREFTERIRTLPHWARRNVFVLDSQGFEIFDRELPSPVA...NRT K4KEM7_SIMAS/61-125 >>>
我们还可以检查对齐中可用的序列(SeqRecord)以及以下：
>>> for align in alignment: ... print(align.seq) ... MQNTPAERLPAIIEKAKSKHDINVWLLDRQGRDLLEQRVPAKVATVANQLRGRKRRAFARHREGP AKQRGIAGLEEWLHRLDHSEAIPIFLIDEAGKDLLEREVPADITA---RLDRRREHGEHGVRKKP ARRHGQEYFQQWLERQPKKVKEQVFAVDQFGRELLGRPLPEDMAPMLIALNYRNRESHAQVDKKP TRRHGPESFRFWLERQPVEARDRIYAIDRSGAEILDRPIPRGMAPLFKVLSFRNREDQGLVNNKP AINRNTQQLTQDLRAMPNWSLRFVYIVDRNNQDLLKRPLPPGIMVLAPRLTAKHPYDKVQDRNRK AVNATEREFTERIRTLPHWARRNVFVLDSQGFEIFDRELPSPVADLMRKLDLDRPFKKLERKNRT >>>
多重对齐
一般来说，大部分序列比对文件都包含单个比对数据，使用起来就足够了 read 方法来解析它。在多序列比对概念中，比较两个或多个序列以获得它们之间的最佳子序列匹配，并在单个文件中产生多序列比对。
如果输入的序列比对格式包含多个序列比对，那么我们需要使用 parse 方法而不是 read 方法如下：
>>> from Bio import AlignIO
>>> alignments = AlignIO.parse(open("PF18225_seed.txt"), "stockholm")
>>> print(alignments)
<generator object parse at 0x000001CD1C7E0360>
>>> for alignment in alignments:
... print(alignment)
...
SingleLetterAlphabet() alignment with 6 rows and 65 columns
MQNTPAERLPAIIEKAKSKHDINVWLLDRQGRDLLEQRVPAKVA...EGP B7RZ31_9GAMM/59-123
AKQRGIAGLEEWLHRLDHSEAIPIFLIDEAGKDLLEREVPADIT...KKP A0A0C3NPG9_9PROT/58-119
ARRHGQEYFQQWLERQPKKVKEQVFAVDQFGRELLGRPLPEDMA...KKP A0A143HL37_9GAMM/57-121
TRRHGPESFRFWLERQPVEARDRIYAIDRSGAEILDRPIPRGMA...NKP A0A0X3UC67_9GAMM/57-121
AINRNTQQLTQDLRAMPNWSLRFVYIVDRNNQDLLKRPLPPGIM...NRK B3PFT7_CELJU/62-126
AVNATEREFTERIRTLPHWARRNVFVLDSQGFEIFDRELPSPVA...NRT K4KEM7_SIMAS/61-125
>>>
在这里，parse 方法返回可迭代的对齐对象，并且可以对其进行迭代以获得实际的对齐。
成对序列比对
成对序列比对 一次仅比较两个序列并提供最佳的序列比对。 Pairwise 很容易理解，并且可以从结果序列比对中推断出来。
Biopython 提供了一个特殊的模块， Bio.pairwise2 使用成对方法识别比对序列。 Biopython 使用最好的算法来找到比对序列，它与其他软件相当。
让我们写一个例子来使用pairwise模块找到两个简单和假设序列的序列比对。这将帮助我们理解序列比对的概念以及如何使用 Biopython 对其进行编程。
步骤 1
导入模块 成对2 使用下面给出的命令：
>>> from Bio import pairwise2
步骤 2
创建两个序列，seq1和seq2：
>>> from Bio.Seq import Seq
>>> seq1 = Seq("ACCGGT")
>>> seq2 = Seq("ACGT")
步骤 3
使用下面的代码行调用pairwise2.align.globalxx方法以及seq1和seq2来查找对齐：
>>> alignments = pairwise2.align.globalxx(seq1, seq2)
Here, globalxx 方法执行实际工作并在给定序列中找到所有可能的最佳比对。实际上，Bio.pairwise2 提供了相当多的方法，它们遵循以下约定来查找不同场景中的比对。
<sequence alignment type>XY
这里，序列比对类型是指比对类型，它可能是 global or 当地的。全球的 类型是通过考虑整个序列来查找序列比对。本地类型也通过查看给定序列的子集来查找序列比对。这将是乏味的，但可以更好地了解给定序列之间的相似性。
-
X 指匹配分数。可能的值是 x(完全匹配)、m(基于相同字符的分数)、d(用户提供的带有字符和匹配分数的字典)，最后是 c(用户定义的函数以提供自定义评分算法)。
-
Y 是指空位罚分。可能的值是 x(无空位罚分)、s(两个序列的罚分相同)、d(每个序列的罚分不同)和最后 c(提供自定义空位罚分的用户定义函数)
因此，localds 也是一种有效的方法，它使用局部比对技术、用户提供的匹配字典和用户提供的两个序列的空位罚分来查找序列比对。
>>> test_alignments = pairwise2.align.localds(seq1, seq2, blosum62, -10, -1)
在这里，blosum62 指的是 pairwise2 模块中提供匹配分数的字典。 -10 是指空位开放罚分，-1 是指空位扩展罚分。
步骤 4
循环遍历可迭代的对齐对象并获取每个单独的对齐对象并打印它。
>>> for alignment in alignments:
... print(alignment)
...
('ACCGGT', 'A-C-GT', 4.0, 0, 6)
('ACCGGT', 'AC--GT', 4.0, 0, 6)
('ACCGGT', 'A-CG-T', 4.0, 0, 6)
('ACCGGT', 'AC-G-T', 4.0, 0, 6)
步骤 5
Bio.pairwise2 模块提供了格式化方法 format_alignment 以更好地可视化结果：
>>> from Bio.pairwise2 import format_alignment >>> alignments = pairwise2.align.globalxx(seq1, seq2) >>> for alignment in alignments: ... print(format_alignment(*alignment)) ... ACCGGT | | || A-C-GT Score=4 ACCGGT || || AC--GT Score=4 ACCGGT | || | A-CG-T Score=4 ACCGGT || | | AC-G-T Score=4 >>>
Biopython 还提供了另一个模块来进行序列比对，Align。该模块提供了一组不同的 API 来简单地设置算法、模式、比赛分数、空位罚分等参数，简单看一下 Align 对象如下：
>>> from Bio import Align >>> aligner = Align.PairwiseAligner() >>> print(aligner) Pairwise sequence aligner with parameters match score: 1.000000 mismatch score: 0.000000 target open gap score: 0.000000 target extend gap score: 0.000000 target left open gap score: 0.000000 target left extend gap score: 0.000000 target right open gap score: 0.000000 target right extend gap score: 0.000000 query open gap score: 0.000000 query extend gap score: 0.000000 query left open gap score: 0.000000 query left extend gap score: 0.000000 query right open gap score: 0.000000 query right extend gap score: 0.000000 mode: global >>>
支持序列比对工具
Biopython 通过 Bio.Align.Applications 模块提供了很多序列比对工具的接口。下面列出了一些工具：
- ClustalW
- MUSCLE
- EMBOSS针和水
让我们在 Biopython 中编写一个简单的示例，通过最流行的比对工具 ClustalW 创建序列比对。
步骤 1 ： 下载 Clustalw 程序 http://www.clustal.org/download/current/ 并安装它。此外，使用“clustal”安装路径更新系统 PATH。
步骤 2 ： 从模块 Bio.Align.Applications 中导入 ClustalwCommanLine。
>>> from Bio.Align.Applications import ClustalwCommandline
步骤 3 ： 通过调用 ClustalwCommanLine 来设置 cmd，输入文件 opuntia.fasta 在 Biopython 包中可用。 https://raw.githubusercontent.com/biopython/biopython/master/Doc/examples/opuntia.fasta
>>> cmd = ClustalwCommandline("clustalw2",
infile="/path/to/biopython/sample/opuntia.fasta")
>>> print(cmd)
clustalw2 -infile=fasta/opuntia.fasta
步骤 4 ：调用 cmd() 将运行 clustalw 命令并给出结果的输出 对齐文件，opuntia.aln。
>>> stdout, stderr = cmd()
步骤 5 ：读取并打印对齐文件如下：
>>> from Bio import AlignIO
>>> align = AlignIO.read("/path/to/biopython/sample/opuntia.aln", "clustal")
>>> print(align)
SingleLetterAlphabet() alignment with 7 rows and 906 columns
TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA
gi|6273285|gb|AF191659.1|AF191
TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA
gi|6273284|gb|AF191658.1|AF191
TATACATTAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA
gi|6273287|gb|AF191661.1|AF191
TATACATAAAAGAAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA
gi|6273286|gb|AF191660.1|AF191
TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA
gi|6273290|gb|AF191664.1|AF191
TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA
gi|6273289|gb|AF191663.1|AF191
TATACATTAAAGGAGGGGGATGCGGATAAATGGAAAGGCGAAAG...AGA
gi|6273291|gb|AF191665.1|AF191
>>>