教 程 目 录

序列比对是按特定顺序排列两个或多个序列(DNA,RNA或蛋白质序列)的过程,以确定它们之间的相似区域.

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

让我们在本章中了解Biopython提供的一些重要特性 :

解析序列比对

Biopython提供了模块,Bio.AlignIO读取和写入序列比对.在生物信息学中,有许多格式可用于指定与早期学习的序列数据类似的序列比对数据. Bio.AlignIO提供类似于Bio.SeqIO的API,除了Bio.SeqIO处理序列数据,Bio.AlignIO处理序列比对数据.

在开始学习之前,让我们从Internet下载样本序列比对文件.

要下载样本文件,请按照以下步骤 : 去;

步骤1 : 打开您喜欢的浏览器并转到 http://pfam.xfam.org/family/browse 网站.它将按字母顺序显示所有Pfam系列.

第2步 : 选择具有较少种子值的任何一个家庭.它包含最少的数据,使我们能够轻松地进行对齐.在这里,我们选择/点击了PF18225并打开了 http://pfam.xfam.org/family/PF18225 并显示有关它的完整详细信息,包括序列比对.

步骤3 : 转到对齐部分并以斯德哥尔摩格式(PF18225_seed.txt)下载序列比对文件.

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

导入Bio.AlignIO模块>>> from Bio import AlignIO

使用read方法读取对齐. 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

>>>

多个对齐

通常,大多数序列比对文件都包含单个比对数据,它足以使用读取方法来解析它.在多序列比对概念中,比较两个或更多序列以获得它们之间的最佳子序列匹配,并在单个文件中导致多序列比对.

如果输入序列比对格式包含多个序列比对,然后我们需要使用parse方法而不是read方法,如下所示 :>>> from Bio import AlignIO

>>> alignments = AlignIO.parse(open("PF18225_seed.txt"), "stockholm")

>>> print(alignments)

>>> 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方法返回可迭代的对齐对象,可以迭代它以获得实际的对齐.

成对序列对齐

成对序列比对一次只比较两个序列,并提供最佳的序列比对. 成对易于理解,并且可以从结果序列比对中推断出来.

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

让我们编写一个例子,使用成对模块找到两个简单和假设序列的序列比对.这将有助于我们理解序列比对的概念以及如何使用Biopython对其进行编程.

步骤1

导入模块pairwise2使用下面给出的命令 :>>> 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)

此处,globalxx方法执行实际工作并找到最佳效果给定序列中可能的比对.实际上,Bio.pairwise2提供了一系列遵循以下惯例的方法,以在不同情况下找到对齐.XY

这里,序列比对类型是指对齐类型,可以是global或local.全局类型通过考虑整个序列来发现序列比对.本地类型也通过查看给定序列的子集来寻找序列比对.这将是乏味的,但更好地了解给定序列之间的相似性.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步 : 从 http://www.clustal.org/download/current/下载Clustalw程序并安装它.另外,使用"clustal"安装路径更新系统PATH.

步骤2 : 从模块Bio.Align.Applications导入ClustalwCommanLine.>>> from Bio.Align.Applications import ClustalwCommandline

步骤3 : 通过使用输入文件调用ClustalwCommanLine来设置cmd,在Biopython包中提供opuntia.fasta.  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

>>>

相关新手教程:

vpa函数python_Biopython - 序列比对相关推荐

  1. 学习MATLAB的第一天,梳理一些见到的函数。1.matlab中sin、cos、tan三角函数问题。2.abs函数。3.vpa函数。4.disp函数。5.class函数。6.logical函数。

    1.matlab中sin.cos.tan三角函数问题. 在MATLAB中三角函数sin.cos.tan都是以弧度为单位的.例如sin()在括号中输入的数系统默认为输入的是弧度值.若想要输入角度值,可以 ...

  2. MATLAB中的vpa函数简单实用记录——精度控制

    vpa函数有两种语法格式: vpa(x) vpa(x,d) 下面是MATLAB帮助文档上的解释: vpa(x) uses variable-precision floating-point arith ...

  3. matlab指令vpa(j10),matlab中vpa函数

    Matlab中矩阵函数_IT/计算机_专业资料.Matlab中矩阵函数 矩阵转置... Matlab 中 solve 函数主要是用来求解线性方程组的解析解或者精确解.对于得 出的结果是符号变量,可以通 ...

  4. oracle表,视图,存储过程,函数,序列.....查询

    查询存储过程,函数,序列,表,视图的名字 select object_name from user_objects where object_type = 'PROCEDURE' select obj ...

  5. MATLAB解隐函数方程时符号表达式转化为数值的方法-用vpa函数

    今天在解决一个小问题时,遇到解隐函数方程,中间涉及一个解的传递问题,才好好研究了一下 syms这个语句的一些语法规则,最终用vpa这个函数解决了符号表达式到double数值的转化. syms是符号函数 ...

  6. 关于φ与Φ函数与序列中分数个数的讨论

    关于 φ\varphiφ与Φ\PhiΦ函数与序列中分数个数的讨论 序列0/n,1/n⋯(n−1)/n0/n,1/n\cdots (n-1)/n0/n,1/n⋯(n−1)/n中分数个数 考虑序列0n,1 ...

  7. vpa函数python_关于取MATLAB的有效位数问题,以及vpa函数

    vpa函数 设用于设置数据的有效位数, 但是如果用这个方式处理数据后返回的则是 sys变量, 处理很不方便 digits() 设置了以后,后面的数据都是按这种设置的显示 format:设置输出格式 对 ...

  8. Matlab中vpa一直在忙,matlab vpa 函数是什么意思?

    matlab控制运算精度用的是digits和vpa这两个函数 digits用于规定运算精度,比如: digits(20); 这个语句就规定了运算精度是20位有效数字.但并不是规定了就可以使用,因为实际 ...

  9. 【Matlab学习笔记】控制运算精度digits和vpa函数

    matlab控制运算精度用的是digits和vpa这两个函数 digits用于规定运算精度,比如: digits(20); 这个语句就规定了运算精度是20位有效数字.但并不是规定了就可以使用,因为实际 ...

最新文章

  1. 如何实现远程控制你的电脑? 网穿软件
  2. 10版微机监测怎么显示服务器,铁路信号网络版微机监测系统的研究
  3. 升级鸿蒙系统的手机名单,倒计时2天!首批鸿蒙OS适配名单确定,你的手机在列吗?...
  4. shell调用python函数_shell调用python函数
  5. 32 位的有符号整数_leetcode 7 整数反转
  6. 面试求职中你需要了解的Java面向对象
  7. 体检异常率98%?数据分析告诉你如今的90后身体状况到底有多差?
  8. buffer string builder简单说明
  9. 99年毕业设计获优的程序-图书管理程序 续
  10. c++简单的加法函数
  11. 数据结构:线性表理论题目集
  12. 经典图像分割方法总结
  13. python中调用shell命令
  14. 从 VDN 到 QMIX的学习笔记
  15. 深入解析MySQL索引原理
  16. Java实现文件搜索
  17. 2022-11-18 mysql列存储引擎-assert failed on i < m_idx.size() at rc_attr.h:342-问题分析
  18. NIT考试感想与复习unity基础
  19. 微信域名防红防屏蔽技术,微信域名总是被封要怎么解决
  20. 计算机课外活动兴趣小组内容,学校课外兴趣小组活动总结

热门文章

  1. 厚物科技PXIe/PXI一体化测控平台HW-1043d
  2. 做数据分析有前景吗?
  3. 后端程序员转行前端,强烈推荐这6个前端UI框架,第二款小程序UI框架颜值最高!
  4. 不良资产评估方法改进研究
  5. Java上传文件到minio
  6. amd速龙黑苹果内核补丁_Lilu.kext v1.4.2 黑苹果内核扩展补丁
  7. linux运维监控内容,Linux运维工程师要掌握的常用监控指标总结
  8. linux bios设置重置,如何将BIOS重置为默认设置 | MOS86
  9. 移相全桥PWM发波要求
  10. 春田花花幼稚园校歌 (普通话版)铃声 春田花花幼稚园校歌 (普通...