Biopython 简明教程

Biopython - Sequence Alignments

Sequence alignment 是按特定顺序布置两个或多个序列(DNA、RNA 或肽序列)的过程,以识别它们之间的相似区域。

识别相似区域使我们能够推断大量信息,例如物种之间有哪些性状得到了保守,不同物种在遗传上有多接近,物种如何进化等。Biopython 为序列比对提供了广泛的支持。

让我们在本章中学习 Biopython 提供的一些重要功能 -

Parsing Sequence Alignment

Biopython 提供了一个模块 Bio.AlignIO 来读写序列比对。在生物信息学中,有很多可用的格式可用于指定序列比对数据,类似于前面学到的序列数据。Bio.AlignIO 提供类似于 Bio.SeqIO 的 API,不同之处在于 Bio.SeqIO 用作序列数据,而 Bio.AlignIO 用作序列比对数据。

在开始学习之前,我们从互联网下载一个序列比对示例文件。

若要下载示例文件,请执行以下步骤 -

Step 1 - 打开您最喜欢的浏览器并转到 http://pfam.xfam.org/family/browse 网站。它会按字母顺序列出所有 Pfam 家族。

Step 2 - 选择任何种子值较少的家族。它包含最少的数据,使我们能够轻松使用比对。在此,我们选择了/点击了 PF18225,然后它打开转到 http://pfam.xfam.org/family/PF18225 并显示有关它的完整详细信息,包括序列比对。

Step 3 - 转到比对部分并下载斯德哥尔摩格式的序列比对文件 (PF18225_seed.txt)。

让我们尝试使用 Bio.AlignIO 读取下载的序列比对文件,如下所示 -

Import Bio.AlignIO module

>>> from Bio import AlignIO

使用read方法进行读取比对。read方法用于读取给定文件中可用的单个比对数据。如果给定文件包含多个比对,我们可以使用parse方法。parse方法返回可迭代比对对象,类似于Bio.SeqIO模块中的parse方法。

>>> alignment = AlignIO.read(open("PF18225_seed.txt"), "stockholm")

Print the alignment object.

>>> 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
>>>

Multiple Alignments

通常情况下,大多数序列比对文件包含单个比对数据,只要使用 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 Sequence Alignment

Pairwise sequence alignment 一次只比较两个序列,并提供尽可能好的序列比对。 Pairwise 易于理解,并且可以从得到的序列比对中推断出例外。

Biopython提供了一个特殊模块 Bio.pairwise2 ,以使用成对方法识别比对序列。Biopython应用最佳算法找到比对序列,并且它与其他软件相当。

让我们通过使用成对模块编写一个示例来查找两个简单假设序列的序列比对。这将帮助我们理解序列比对的概念以及如何使用Biopython对其进行编程。

Step 1

使用下面给出的命令导入模块 pairwise2

>>> from Bio import pairwise2

Step 2

创建两个序列seq1和seq2:

>>> from Bio.Seq import Seq
>>> seq1 = Seq("ACCGGT")
>>> seq2 = Seq("ACGT")

Step 3

调用方法pairwise2.align.globalxx以及seq1和seq2使用下面这一行代码找到比对:

>>> alignments = pairwise2.align.globalxx(seq1, seq2)

在这里, globalxx 方法执行实际工作并在给定序列中找到所有可能的最佳比对。事实上,Bio.pairwise2提供了一系列方法,它们遵循下面的约定来在不同场景中找到比对。

<sequence alignment type>XY

在这里,序列比对类型是指比对类型,它可以是全局的或局部的。全局类型通过考虑整个序列来找到序列比对。局部类型通过查看给定序列的子集找到序列比对。这样做很乏味,但是可以更好地了解给定序列之间的相似性。

  1. X表示匹配得分。可能的值为x(精确匹配)、m(基于相同字符的得分)、d(带有字符和匹配得分的用户提供的字典)和最后c(提供自定义评分算法的用户定义函数)。

  2. Y表示缺口罚分。可能的值为x(没有缺口罚分)、s(这两个序列有相同的罚分)、d(为每个序列设置不同的罚分)和最后c(提供自定义缺口罚分的用户定义函数)。

因此,localds也是一种有效的方法,它使用局部比对技术、用户提供的匹配字典和用户为两个序列提供的缺口罚分找到序列比对。

>>> test_alignments = pairwise2.align.localds(seq1, seq2, blosum62, -10, -1)

在这里,blosum62是指pairwise2模块中可用的字典以提供匹配得分。-10是指缺口打开罚分,-1是指缺口扩展罚分。

Step 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)

Step 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
>>>

Support for Sequence Alignment Tools

Biopython 通过 Bio.Align.Applications 模块提供大量序列比对工具的接口。下面列出了一些工具 -

  1. ClustalW

  2. MUSCLE

  3. EMBOSS needle and water

让我们编写一个简单的 Biopython 示例,通过最流行的对齐工具 ClustalW 创建序列对齐。

Step 1 − 从 http://www.clustal.org/download/current/ 下载 Clustalw 程序并将其安装。此外,使用“clustal”安装路径更新系统 PATH。

Step 2 − 从 Bio.Align.Applications 模块导入 ClustalwCommanLine。

>>> from Bio.Align.Applications import ClustalwCommandline

Step 3 − 使用 Biopython 包中可用的输入文件 opuntia.fasta,通过调用 ClustalwCommanLine 设置 cmd。 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

Step 4 − 调用 cmd() 将运行 clustalw 命令,并给出结果对齐文件 opuntia.aln 的输出。

>>> stdout, stderr = cmd()

Step 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
>>>