BioPython:从Blast输出fi中提取序列id

2024-04-28 20:21:19 发布

您现在位置:Python中文网/ 问答频道 /正文

我有一个XML格式的BLAST输出文件。它是22个查询序列,每个序列报告50次点击。我想提取所有的50x22点击。这是我目前拥有的代码,但它只从第一个查询中提取50个匹配项。

from Bio.Blast import NCBIXM
blast_records = NCBIXML.parse(result_handle)
blast_record = blast_records.next()

save_file = open("/Users/jonbra/Desktop/my_fasta_seq.fasta", 'w')

for alignment in blast_record.alignments:
    for hsp in alignment.hsps:
            save_file.write('>%s\n' % (alignment.title,))
save_file.close()

有人有什么建议可以提取所有的点击率吗?我想我得用些别的东西而不是排列。 希望这是清楚的。谢谢!

乔恩


Tags: 文件inforsave格式报告序列xml
2条回答

我用这个代码提取了所有的结果

from Bio.Blast import NCBIXML
for record in NCBIXML.parse(open("rpoD.xml")) :
    print "QUERY: %s" % record.query
    for align in record.alignments :
        print " MATCH: %s..." % align.title[:60]
        for hsp in align.hsps :
            print " HSP, e=%f, from position %i to %i" \
                % (hsp.expect, hsp.query_start, hsp.query_end)
            if hsp.align_length < 60 :
                 print "  Query: %s" % hsp.query
                 print "  Match: %s" % hsp.match
                 print "  Sbjct: %s" % hsp.sbjct
            else :
                 print "  Query: %s..." % hsp.query[:57]
                 print "  Match: %s..." % hsp.match[:57]
                 print "  Sbjct: %s..." % hsp.sbjct[:57]


print "Done"

或是为了少一些细节

from Bio.Blast import NCBIXML
for record in NCBIXML.parse(open("NC_003197.xml")) :
    #We want to ignore any queries with no search results:
    if record.alignments :
        print "QUERY: %s..." % record.query[:60]
        for align in record.alignments :
            for hsp in align.hsps :
                print " %s HSP, e=%f, from position %i to %i" \
                % (align.hit_id, hsp.expect, hsp.query_start, hsp.query_end)
print "Done"

我用过这个网站

http://www2.warwick.ac.uk/fac/sci/moac/currentstudents/peter_cock/python/rpsblast/

这会得到所有的记录。与原始版本相比,新颖之处在于

for blast_record in blast_records

这是一个python习惯用法,用于遍历“类列表”对象(如blast_记录)中的项(检查CBIXML module documentation显示parse()确实返回一个迭代器)

from Bio.Blast import NCBIXM
blast_records = NCBIXML.parse(result_handle)

save_file = open("/Users/jonbra/Desktop/my_fasta_seq.fasta", 'w')

for blast_record in blast_records:
  for alignment in blast_record.alignments:
      for hsp in alignment.hsps:
            save_file.write('>%s\n' % (alignment.title,))
  #here possibly to output something to file, between each blast_record
save_file.close()

相关问题 更多 >