luria的个人博客分享 http://blog.sciencenet.cn/u/luria

博文

使用Biopython解析Blast结果

已有 2010 次阅读 2018-7-19 15:57 |个人分类:Metagenomics|系统分类:科研笔记| blast, python, biopython, blast格式, xml

同一个比对结果:

format 0:

format 6:

format 5:

1种格式适合人看,机器识别起来较困难,第6种(或第7种)格式使用也较多,但是遗憾的是信息不全。第5种格式是标准的XML格式,人看起来不方便,但机器似乎很喜欢识别这种格式,这篇博文重点使用Bio.Blast.NCBIXML解析Blast5种格式(XML格式)。

for record in NCBIXML.parse(open(blast_fmt5.xml))

其中recorddescriptionsalignmentsmultiple_alignment列表对象, 还有BlastHeaderDatabaseReportParameters方法。

列表对象可以通过循环获取到,例如:

for element in record. descriptions

其中elementtitlescorebits enum_alignments对象组成。

recorddesciptions对象中,给出了Subject IDbitscore值,evalue值等信息

 for element in record.alignments

其中elementtitlehit_idhit_deflengthhsps组成。

其中titledescriptions里的title相同,hit_id不是subject id而且自编的一个idhit_defsubject id

其中最重要的是hsps,它可以再用一个循环获得。

for hsp in element.hsps

 

按照以下代码将其中querymatchsbjct打印出来(一行一个)如下:

代码如下:

当不清楚如何的对应关系时,可以查看xml文件,如下:

例如:当不清楚subject序列id号是什么时,可以

grep [Subject ID] blast5.xml

 

但是很坑的是,很难找到query id,因为它并没有在文档中给出,甚至连help(record)也看不到,而且用上面grep [query id] blast5.xml的方式搜索出来的是query-def,但是用query_def会报错,因为根本就没有这个属性!

原来,它藏在record中的属性中,通过record.query提取出来!

还有一个SearchIO更好用,先mark一下,等整理出来再分享



http://blog.sciencenet.cn/blog-2970729-1124895.html

上一篇:Nanopore basecalling之Albacore
下一篇:采用Stacks2分析RADseq数据[无参]

0

该博文允许注册用户评论 请点击登录 评论 (0 个评论)

数据加载中...
扫一扫,分享此博文

Archiver|手机版|科学网 ( 京ICP备14006957 )

GMT+8, 2019-10-17 06:39

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部