2009-12-20 4 views
6

NCBI xml BLAST 파일에서 첫 번째 히트 만 추출하려고합니다. 다음으로는 첫 번째 HSP 만 얻고 싶습니다. 마지막 단계에서 나는 최고 점수에 기초하여 이것을 얻고 싶다. XML NCBI BLAST 파일에서 첫 번째 적중 요소를 추출하는 방법은 무엇입니까?

<?xml version="1.0"?> 
<!DOCTYPE BlastOutput PUBLIC "-//NCBI//NCBI BlastOutput/EN" "http://www.ncbi.nlm.nih.gov/dtd/NCBI_BlastOutput.dtd"> 
<BlastOutput> 
    <BlastOutput_program>blastx</BlastOutput_program> 
    <BlastOutput_version>blastx 2.2.22 [Sep-27-2009]</BlastOutput_version> 
    <BlastOutput_reference>~Reference: Altschul, Stephen F., Thomas L. Madden, Alejandro A. Schaffer, ~Jinghui Zhang, Zheng Zhang, Webb Miller, and David J. Lipman (1997), ~&quot;Gapped BLAST and PSI-BLAST: a new generation of protein database search~programs&quot;, Nucleic Acids Res. 25:3389-3402.</BlastOutput_reference> 
    <BlastOutput_db>/Applications/blast/db/viral1.protein.faa</BlastOutput_db> 
    <BlastOutput_query-ID>lcl|1_0</BlastOutput_query-ID> 
    <BlastOutput_query-def>DSAD-090629_plate11A01a.g1 CHROMAT_FILE: DSAD-090629_plate11A01a.g1 PHD_FILE: DSAD-090629_plate11A01a.g1.phd.1 CHEM: term DYE: big TIME: Thu Sep 17 15:33:59 2009 TEMPLATE: DSAD-090629_plate11A01a DIRECTION: rev</BlastOutput_query-def> 
    <BlastOutput_query-len>1024</BlastOutput_query-len> 
    <BlastOutput_param> 
    <Parameters> 
     <Parameters_matrix>BLOSUM62</Parameters_matrix> 
     <Parameters_expect>1e-05</Parameters_expect> 
     <Parameters_gap-open>11</Parameters_gap-open> 
     <Parameters_gap-extend>1</Parameters_gap-extend> 
     <Parameters_filter>F</Parameters_filter> 
    </Parameters> 
    </BlastOutput_param> 
    <BlastOutput_iterations> 
    <Iteration> 
     <Iteration_iter-num>1</Iteration_iter-num> 
     <Iteration_query-ID>lcl|1_0</Iteration_query-ID> 
     <Iteration_query-def>DSAD-090629_plate11A01a.g1 CHROMAT_FILE: DSAD-090629_plate11A01a.g1 PHD_FILE: DSAD-090629_plate11A01a.g1.phd.1 CHEM: term DYE: big TIME: Thu Sep 17 15:33:59 2009 TEMPLATE: DSAD-090629_plate11A01a DIRECTION: rev</Iteration_query-def> 
     <Iteration_query-len>1024</Iteration_query-len> 
     <Iteration_stat> 
     <Statistics> 
      <Statistics_db-num>68007</Statistics_db-num> 
      <Statistics_db-len>19518578</Statistics_db-len> 
      <Statistics_hsp-len>0</Statistics_hsp-len> 
      <Statistics_eff-space>0</Statistics_eff-space> 
      <Statistics_kappa>0.041</Statistics_kappa> 
      <Statistics_lambda>0.267</Statistics_lambda> 
      <Statistics_entropy>0.14</Statistics_entropy> 
     </Statistics> 
     </Iteration_stat> 
     <Iteration_message>No hits found</Iteration_message> 
    </Iteration> 
    <Iteration> 
<Iteration> 
     <Iteration_iter-num>6</Iteration_iter-num> 
     <Iteration_query-ID>lcl|6_0</Iteration_query-ID> 
     <Iteration_query-def>DSAD-090629_plate11A05a.g1 CHROMAT_FILE: DSAD-090629_plate11A05a.g1 PHD_FILE: DSAD-090629_plate11A05a.g1.phd.1 CHEM: term DYE: big TIME: Thu Sep 17 15:33:59 2009 TEMPLATE: DSAD-090629_plate11A05a DIRECTION: rev</Iteration_query-def> 
     <Iteration_query-len>1068</Iteration_query-len> 
     <Iteration_hits> 
     <Hit> 
      <Hit_num>1</Hit_num> 
      <Hit_id>gnl|BL_ORD_ID|23609</Hit_id> 
      <Hit_def>gi|38707884|ref|NP_945016.1| Putative ribose-phosphate pyrophosphokinase [Enterobacteria phage Felix 01]</Hit_def> 
      <Hit_accession>23609</Hit_accession> 
      <Hit_len>293</Hit_len> 
      <Hit_hsps> 
      <Hsp> 
       <Hsp_num>1</Hsp_num> 
       <Hsp_bit-score>49.2914</Hsp_bit-score> 
       <Hsp_score>116</Hsp_score> 
       <Hsp_evalue>5.15408e-06</Hsp_evalue> 
       <Hsp_query-from>580</Hsp_query-from> 
       <Hsp_query-to>792</Hsp_query-to> 
       <Hsp_hit-from>202</Hsp_hit-from> 
       <Hsp_hit-to>273</Hsp_hit-to> 
       <Hsp_query-frame>-1</Hsp_query-frame> 
       <Hsp_identity>26</Hsp_identity> 
       <Hsp_positive>45</Hsp_positive> 
       <Hsp_gaps>2</Hsp_gaps> 
       <Hsp_align-len>73</Hsp_align-len> 
       <Hsp_qseq>MHIIGDVE--GRTCILVDDMVDTAGTLCHAAKALKERGAAKVYAYCTHPVLSGRAIENIENSVLDELVVTNTI</Hsp_qseq> 
       <Hsp_hseq>MRILDDVDLTDKTVMILDDICDGGRTFVEAAKHLREAGAKRVELYVTHGIFS-KDVENLLDNGIDHIYTTNSL</Hsp_hseq> 
       <Hsp_midline>M I+ DV+ +T +++DD+ D T AAK L+E GA +V Y TH + S + +EN+ ++ +D + TN++</Hsp_midline> 
      </Hsp> 
      </Hit_hsps> 
     </Hit> 
     <Hit> 
      <Hit_num>2</Hit_num> 
      <Hit_id>gnl|BL_ORD_ID|2466</Hit_id> 
      <Hit_def>gi|51557505|ref|YP_068339.1| large tegument protein [Suid herpesvirus 1]</Hit_def> 
      <Hit_accession>2466</Hit_accession> 
      <Hit_len>3084</Hit_len> 
      <Hit_hsps> 
      <Hsp> 
       <Hsp_num>1</Hsp_num> 
       <Hsp_bit-score>48.9062</Hsp_bit-score> 
       <Hsp_score>115</Hsp_score> 
       <Hsp_evalue>6.70494e-06</Hsp_evalue> 
       <Hsp_query-from>369</Hsp_query-from> 
       <Hsp_query-to>875</Hsp_query-to> 
       <Hsp_hit-from>2312</Hsp_hit-from> 
       <Hsp_hit-to>2465</Hsp_hit-to> 
       <Hsp_query-frame>-2</Hsp_query-frame> 
       <Hsp_identity>52</Hsp_identity> 
       <Hsp_positive>70</Hsp_positive> 
       <Hsp_gaps>4</Hsp_gaps> 
       <Hsp_align-len>173</Hsp_align-len> 
      <Hsp_qseq>APESQEPGASTWRSSTSVVKKGQPSQK*CTSSVTSKAVPASWSTTWSTLPAPCATPPKR*KSAAPPRSTPTAPTRCCPAAPSRTSRIPSWTSWWSPTPSRCPLRRSPARVFASSTSPR-SSPKRSAASATKNRSAP---CSAKRNWPDHTAPPRAGLFALPPEAGRKPQGGLV</Hsp_qseq> 
      <Hsp_hseq>APPAQKPPAQPATAAATTAPKATPQTQPPTRAQTQTAPPPPSAAT-----AAAQVPPQ------PPSSQPAAKPRGAPPAPPAPP--PPSAQTTLPRPAAPPAPPPPS---AQTTLPRPAPPPPSAPAATPTPPAPGPAPSAKKSDGDRIVEPKAG---APPDVRDAKFGGKV</Hsp_hseq> 
      <Hsp_midline>AP +Q+P A ++ + K P + T + T A P + T  A PP+  PP S P A R P AP  P  P P+ P P+ A +T PR + P SA +AT AP SAK++ D P+AG PP+  GG V</Hsp_midline> 
     </Hsp> 
     </Hit_hsps> 
    </Hit> 
    </Iteration_hits> 
    <Iteration_stat> 
    <Statistics> 
     <Statistics_db-num>68007</Statistics_db-num> 
     <Statistics_db-len>19518578</Statistics_db-len> 
     <Statistics_hsp-len>0</Statistics_hsp-len> 
     <Statistics_eff-space>0</Statistics_eff-space> 
     <Statistics_kappa>0.041</Statistics_kappa> 
     <Statistics_lambda>0.267</Statistics_lambda> 
     <Statistics_entropy>0.14</Statistics_entropy> 
    </Statistics> 
    </Iteration_stat> 
</Iteration> 

은 기본적으로 각 쿼리 검색이 반복 요소를 만듭니다 는 XML 파일의 샘플 분명히 여기 일을 확인합니다. 각 반복에는 여러 개의 HSP가있을 수있는 여러 개의 히트가있을 수 있습니다. 첫 번째 히트 만 얻고 싶습니다. 각 반복에서 첫 번째 HSP입니다. BLAST에서 히트가 발견되지 않으면 반복을 무시하고 싶습니다. 나는이 간단한 코드 흥분 : 어떤 도움이 크게 appriciated 될

#!/usr/bin/env python 
from elementtree.ElementTree import parse 
from elementtree import ElementTree as ET 
file = open("/Applications/blast/blanes_viral_nr_results.xml", "r") 
save_file = open("/Applications/blast/Blast_parse_ET.txt", 'w') 
tree = parse(file) 
elem = tree.getroot() 
print elem 
Per_ID =() 

save_file.write('>%s\t%s\t%s\t%s\t%s\t%s\t\n\n\n\n' % ("It_Num\t", "It_ID\t", "Hit_Def\t", "Num\t", "ID\t", "ACC\t")) 
iteration = tree.findall('BlastOutput_iterations/Iteration') 
for iteration in iteration: 
    for hit in iteration.findall('Iteration_hits/Hit'): 
    It_Num = iteration.findtext('Iteration_iter-num') 
    It_ID = iteration.findtext('Iteration_query-def') 
    Hit_Def = hit.findtext('Hit_def') 
    Num = hit.findtext('Hit_num') 
    ID = hit.findtext('Hit_id') 
    DEF = hit.findtext('Hit_def') 
    ACC = hit.findtext('Hit_accession') 
    save_file.write('>%s\t%s\t%s\t%s\t%s\t%s\t' % (It_Num, It_ID[12:26], Hit_Def[1:10], Num, ID, ACC,)) 
    for hsp in hit.findall('Hit_hsps'): 
     HSPN = hsp.findtext('Hsp/Hsp_num') 
     identities = hsp.findtext('Hsp/Hsp_identity') 
     #print 'id: ', identities.rjust(4), 
     length = hsp.findtext('Hsp/Hsp_align-len') 
     #print 'len:', length.rjust(4), 
     Per_ID = int(identities) * 100.0/int(length) 
     #print hsp.findtext('Hsp/Hsp_qseq')[:50] 
     #print hsp.findtext('Hsp/Hsp_midline')[:50] 
     #print hsp.findtext('Hsp/Hsp_hseq')[:50] 
     save_file.write('%s\t%s\t%s\%st\n' % ('***', '%', HSPN, Per_ID)) 
    save_file.write('n\n' %()) 

을!

답변

9

자신 만의 파서를 만드는 것이 재미있을 수 있습니다. 이미 BLAST xml 파일을 파싱 할 수있는 패키지가 있습니다. 원하는 경우 BLAST 인스턴스의 중간 호출을 수행 할 수도 있습니다.

주요 사이트

은 여기에 있습니다 : http://biopython.org/wiki/Biopython

하고 XML BLAST 파서는 여기에 있습니다 : http://biopython.org/DIST/docs/tutorial/Tutorial.html#htoc82

뭔가 같은 :

from Bio.Blast import NCBIXML 
with open('xml/results/file') as handle: 
    all_records = NCBIXML.parse(handle) 
    first_record = all_records.next() 

작동합니다. 필자는 일반적으로 BioPython 파서와 작성자를 좋아하지만 클래스 구조 조직을 좋아하지 않습니다. 그래서 저는 일반적으로 파서를 사용하고 필요한 정보를 구조로 추출합니다. 도움이 YMMV

희망.

+0

덕분에 최고 히트 싱글과 blastoutput의 표로 양식을 볼 수 있습니다 MS Excel에서 XML 출력 파일을 엽니 다. 나는 과거와 비슷한 것을 시도했지만 지금은 제안을 시도했지만 특정 항목을 호출하려고하면 다음 마사지를 계속받습니다. first_record = blast_records [0] TypeError : 'generator'객체 감사합니다, –

+1

나는 대답을 바꿨습니다 ... 나는 잠시 동안 biopython을 사용하지 않았고 발전기를 사용하기 위해 리팩토링 한 것을 잊었습니다. 이 새로운 것을 시도하십시오. – JudoWill

6

나는 Judo가 제안한 것보다 두 번째로 더 자세하게 설명합니다. BioPython 파서를 사용하면 훨씬 더 똑똑하고 똑똑해집니다. 이 더 당신에게 약간을 얻어야한다 :

from Bio.Blast import NCBIXML 
blast = NCBIXML.parse(open('results.xml','rU')) 
for record in blast: 
    if record.alignments: 
     # to print the "best" matches e-score 
     print record.alignments[0].hsps[0].expect 
     # to print the "best" matches bit-score 
     print record.alignments[0].hsps[0].score 
     break 

이 첫 번째 쿼리 (첫 번째와 가장 일치하는 반환) 후 중지됩니다. 동일한 파일 내에서 다른 검색어에 대한 결과가 필요할 것으로 판단되기 때문에 마지막 줄에서 break을 삭제하면됩니다. 난 당신이 최고 인기를 얻으려면 다음 기본 요구 사항을 이해한다면

1

/쿼리 단백질/염기 sequence.Why의 HSP는 서식 NR/NT를 데이터베이스 시스템에 독립 폭발을 설치 해달라고. 입력 옵션

blastall -p {blast programme blastp for protein,blastn for nucleotide} -d {database} -i {input query} -v 1{for top hit} -b 1{alignment of the top hit with query} -m 7{xml blast output} -o example.xml 

은 각 질의 시퀀스에 답장을