2009-11-05 20 views
3

J'ai un fichier de sortie BLAST au format XML. Il s'agit de 22 séquences de requêtes avec 50 occurrences de chaque séquence. Et je veux extraire tous les hits 50x22. C'est le code que j'ai actuellement, mais il extrait seulement les 50 hits de la première requête.BioPython: extraction des ID de séquence à partir d'un fichier de sortie Blast

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

Quelqu'un a des suggestions pour extraire tous les résultats? Je suppose que je dois utiliser autre chose que des alignements. Espérons que c'était clair. Merci!

Jon

+0

Pourquoi le double affichage? http://stackoverflow.com/questions/1684194/python-saving-output-of-a-for-loop-to-file il y a une heure environ? – mjv

+0

Étant donné que la plupart des gens ici n'utilisent probablement pas BioPython, vous pouvez obtenir plus de réponses si vous leur fournissez des liens utiles. –

+0

mjv: La publication précédente concernait la façon d'enregistrer la sortie. C'est à peu près le même code, mais maintenant je veux le changer. Gnibbler: que voulez-vous dire des liens utiles? Comme des liens pour aider à la réponse? J'ai vérifié beaucoup de liens. Comme docs biopython, mais le problème est que j'ai du mal à lire de tels documents – Jon

Répondre

3

Cela devrait permettre d'obtenir tous les enregistrements. La nouveauté par rapport à l'original est le

for blast_record in blast_records 

qui est un langage de python pour parcourir les éléments d'un objet « liste-like », tels que les blast_records (vérification du CBIXML module documentation a montré que parse() retourne en effet un iterator)

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() 
2

J'ai utilisé ce code pour extraire tous les résultats

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" 

ou pour moins de détails

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" 

J'ai utilisé ce site

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