Biopython - Überblick über BLAST

BLAST steht für Basic Local Alignment Search Tool. Es findet Ähnlichkeitsbereiche zwischen biologischen Sequenzen. Biopython bietet das Bio.Blast-Modul für den NCBI BLAST-Betrieb. Sie können BLAST entweder in einer lokalen Verbindung oder über eine Internetverbindung ausführen.

Lassen Sie uns diese beiden Zusammenhänge im folgenden Abschnitt kurz verstehen -

Laufen über das Internet

Biopython bietet das Bio.Blast.NCBIWWW-Modul zum Aufrufen der Online-Version von BLAST. Dazu müssen wir das folgende Modul importieren:

>>> from Bio.Blast import NCBIWWW

Das NCBIWW-Modul bietet eine qblast-Funktion zum Abfragen der BLAST-Onlineversion. https://blast.ncbi.nlm.nih.gov/Blast.cgi. qblast unterstützt alle von der Online-Version unterstützten Parameter.

Um Hilfe zu diesem Modul zu erhalten, verwenden Sie den folgenden Befehl und verstehen Sie die Funktionen -

>>> help(NCBIWWW.qblast) 
Help on function qblast in module Bio.Blast.NCBIWWW: 
qblast(
   program, database, sequence, 
   url_base = 'https://blast.ncbi.nlm.nih.gov/Blast.cgi', 
   auto_format = None, 
   composition_based_statistics = None, 
   db_genetic_code =  None, 
   endpoints = None, 
   entrez_query = '(none)', 
   expect = 10.0, 
   filter = None, 
   gapcosts = None, 
   genetic_code = None, 
   hitlist_size = 50, 
   i_thresh = None, 
   layout = None, 
   lcase_mask = None, 
   matrix_name = None, 
   nucl_penalty = None, 
   nucl_reward = None, 
   other_advanced = None, 
   perc_ident = None, 
   phi_pattern = None, 
   query_file = None, 
   query_believe_defline = None, 
   query_from = None, 
   query_to = None, 
   searchsp_eff = None, 
   service = None, 
   threshold = None, 
   ungapped_alignment = None, 
   word_size = None, 
   alignments = 500, 
   alignment_view = None, 
   descriptions = 500, 
   entrez_links_new_window = None, 
   expect_low = None, 
   expect_high = None, 
   format_entrez_query = None, 
   format_object = None, 
   format_type = 'XML', 
   ncbi_gi = None, 
   results_file = None, 
   show_overview = None, 
   megablast = None, 
   template_type = None, 
   template_length = None
) 
   
   BLAST search using NCBI's QBLAST server or a cloud service provider. 
   
   Supports all parameters of the qblast API for Put and Get. 
   
   Please note that BLAST on the cloud supports the NCBI-BLAST Common 
   URL API (http://ncbi.github.io/blast-cloud/dev/api.html). 
   To use this feature, please set url_base to 'http://host.my.cloud.service.provider.com/cgi-bin/blast.cgi' and 
   format_object = 'Alignment'. For more details, please see 8. Biopython – Overview of BLAST
   
https://blast.ncbi.nlm.nih.gov/Blast.cgi?PAGE_TYPE = BlastDocs&DOC_TYPE = CloudBlast 
   
   Some useful parameters: 
   
   - program blastn, blastp, blastx, tblastn, or tblastx (lower case) 
   - database Which database to search against (e.g. "nr"). 
   - sequence The sequence to search. 
   - ncbi_gi TRUE/FALSE whether to give 'gi' identifier. 
   - descriptions Number of descriptions to show. Def 500. 
   - alignments Number of alignments to show. Def 500. 
   - expect An expect value cutoff. Def 10.0. 
   - matrix_name Specify an alt. matrix (PAM30, PAM70, BLOSUM80, BLOSUM45). 
   - filter "none" turns off filtering. Default no filtering 
   - format_type "HTML", "Text", "ASN.1", or "XML". Def. "XML". 
   - entrez_query Entrez query to limit Blast search 
   - hitlist_size Number of hits to return. Default 50 
   - megablast TRUE/FALSE whether to use MEga BLAST algorithm (blastn only) 
   - service plain, psi, phi, rpsblast, megablast (lower case) 
   
   This function does no checking of the validity of the parameters 
   and passes the values to the server as is. More help is available at: 
   https://ncbi.github.io/blast-cloud/dev/api.html

Normalerweise sind die Argumente der qblast-Funktion grundsätzlich analog zu verschiedenen Parametern, die Sie auf der BLAST-Webseite festlegen können. Dies macht die qblast-Funktion leicht verständlich und verkürzt die Lernkurve, um sie zu verwenden.

Verbinden und Suchen

Um den Prozess des Verbindens und Suchens der BLAST-Online-Version zu verstehen, führen wir eine einfache Sequenzsuche (verfügbar in unserer lokalen Sequenzdatei) gegen den Online-BLAST-Server über Biopython durch.

Step 1 - Erstellen Sie eine Datei mit dem Namen blast_example.fasta im Biopython-Verzeichnis und geben Sie die folgenden Sequenzinformationen als Eingabe ein

Example of a single sequence in FASTA/Pearson format: 
>sequence A ggtaagtcctctagtacaaacacccccaatattgtgatataattaaaattatattcatat
tctgttgccagaaaaaacacttttaggctatattagagccatcttctttgaagcgttgtc 

>sequence B ggtaagtcctctagtacaaacacccccaatattgtgatataattaaaattatattca
tattctgttgccagaaaaaacacttttaggctatattagagccatcttctttgaagcgttgtc

Step 2 - Importieren Sie das NCBIWWW-Modul.

>>> from Bio.Blast import NCBIWWW

Step 3 - Öffnen Sie die Sequenzdatei, blast_example.fasta mit Python IO-Modul.

>>> sequence_data = open("blast_example.fasta").read() 
>>> sequence_data 
'Example of a single sequence in FASTA/Pearson format:\n\n\n> sequence 
A\nggtaagtcctctagtacaaacacccccaatattgtgatataattaaaatt 
atattcatat\ntctgttgccagaaaaaacacttttaggctatattagagccatcttctttg aagcgttgtc\n\n'

Step 4- Rufen Sie nun die qblast-Funktion auf, die Sequenzdaten als Hauptparameter übergibt. Der andere Parameter repräsentiert die Datenbank (nt) und das interne Programm (blastn).

>>> result_handle = NCBIWWW.qblast("blastn", "nt", sequence_data) 
>>> result_handle 
<_io.StringIO object at 0x000001EC9FAA4558>

blast_resultshält das Ergebnis unserer Suche. Es kann zur späteren Verwendung in einer Datei gespeichert und analysiert werden, um die Details zu erhalten. Wir werden im kommenden Abschnitt lernen, wie es geht.

Step 5 - Die gleiche Funktionalität kann auch mit dem Seq-Objekt ausgeführt werden, anstatt die gesamte Fasta-Datei wie unten gezeigt zu verwenden. -

>>> from Bio import SeqIO 
>>> seq_record = next(SeqIO.parse(open('blast_example.fasta'),'fasta')) 
>>> seq_record.id 
'sequence' 
>>> seq_record.seq 
Seq('ggtaagtcctctagtacaaacacccccaatattgtgatataattaaaattatat...gtc', 
SingleLetterAlphabet())

Rufen Sie nun die qblast-Funktion auf, die das Seq-Objekt record.seq als Hauptparameter übergibt.

>>> result_handle = NCBIWWW.qblast("blastn", "nt", seq_record.seq) 
>>> print(result_handle) 
<_io.StringIO object at 0x000001EC9FAA4558>

BLAST weist Ihrer Sequenz automatisch eine Kennung zu.

Step 6 - Das result_handle-Objekt enthält das gesamte Ergebnis und kann zur späteren Verwendung in einer Datei gespeichert werden.

>>> with open('results.xml', 'w') as save_file: 
>>>   blast_results = result_handle.read() 
>>>   save_file.write(blast_results)

Wir werden im späteren Abschnitt sehen, wie die Ergebnisdatei analysiert wird.

Standalone BLAST ausführen

In diesem Abschnitt wird erläutert, wie BLAST im lokalen System ausgeführt wird. Wenn Sie BLAST im lokalen System ausführen, ist es möglicherweise schneller und ermöglicht Ihnen auch, eine eigene Datenbank für die Suche nach Sequenzen zu erstellen.

BLAST anschließen

Im Allgemeinen wird die lokale Ausführung von BLAST aufgrund der Größe, des zusätzlichen Aufwands für die Ausführung der Software und der damit verbundenen Kosten nicht empfohlen. Online BLAST ist ausreichend für grundlegende und fortgeschrittene Zwecke. Natürlich müssen Sie es manchmal lokal installieren.

Bedenken Sie, dass Sie häufig online suchen, was viel Zeit und ein hohes Netzwerkvolumen erfordern kann. Wenn Sie proprietäre Sequenzdaten oder IP-bezogene Probleme haben, wird empfohlen, diese lokal zu installieren.

Dazu müssen wir die folgenden Schritte ausführen:

Step 1- Laden Sie die neueste Blast-Binärdatei über den angegebenen Link herunter und installieren Sie sie - ftp://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/

Step 2- Laden Sie die neueste und notwendige Datenbank über den folgenden Link herunter und entpacken Sie sie - ftp://ftp.ncbi.nlm.nih.gov/blast/db/

Die BLAST-Software bietet viele Datenbanken auf ihrer Website. Laden Sie die Datei alu.n.gz von der Explosionsdatenbank herunter und entpacken Sie sie in den Ordner alu. Diese Datei ist im FASTA-Format. Um diese Datei in unserer Blast-Anwendung zu verwenden, müssen wir zuerst die Datei aus dem FASTA-Format in das Blast-Datenbankformat konvertieren. BLAST bietet eine makeblastdb-Anwendung für diese Konvertierung.

Verwenden Sie das folgende Code-Snippet -

cd /path/to/alu 
makeblastdb -in alu.n -parse_seqids -dbtype nucl -out alun

Durch Ausführen des obigen Codes wird die Eingabedatei alu.n analysiert und die BLAST-Datenbank als mehrere Dateien alun.nsq, alun.nsi usw. erstellt. Jetzt können wir diese Datenbank abfragen, um die Sequenz zu finden.

Wir haben den BLAST auf unserem lokalen Server installiert und haben auch eine BLAST-Beispieldatenbank. alun dagegen abfragen.

Step 3- Lassen Sie uns eine Beispielsequenzdatei erstellen, um die Datenbank abzufragen. Erstellen Sie eine Datei search.fsa und fügen Sie die folgenden Daten ein.

>gnl|alu|Z15030_HSAL001056 (Alu-J) 
AGGCTGGCACTGTGGCTCATGCTGAAATCCCAGCACGGCGGAGGACGGCGGAAGATTGCT 
TGAGCCTAGGAGTTTGCGACCAGCCTGGGTGACATAGGGAGATGCCTGTCTCTACGCAAA 
AGAAAAAAAAAATAGCTCTGCTGGTGGTGCATGCCTATAGTCTCAGCTATCAGGAGGCTG 
GGACAGGAGGATCACTTGGGCCCGGGAGTTGAGGCTGTGGTGAGCCACGATCACACCACT 
GCACTCCAGCCTGGGTGACAGAGCAAGACCCTGTCTCAAAACAAACAAATAA 
>gnl|alu|D00596_HSAL003180 (Alu-Sx) 
AGCCAGGTGTGGTGGCTCACGCCTGTAATCCCACCGCTTTGGGAGGCTGAGTCAGATCAC 
CTGAGGTTAGGAATTTGGGACCAGCCTGGCCAACATGGCGACACCCCAGTCTCTACTAAT 
AACACAAAAAATTAGCCAGGTGTGCTGGTGCATGTCTGTAATCCCAGCTACTCAGGAGGC 
TGAGGCATGAGAATTGCTCACGAGGCGGAGGTTGTAGTGAGCTGAGATCGTGGCACTGTA
CTCCAGCCTGGCGACAGAGGGAGAACCCATGTCAAAAACAAAAAAAGACACCACCAAAGG 
TCAAAGCATA 
>gnl|alu|X55502_HSAL000745 (Alu-J) 
TGCCTTCCCCATCTGTAATTCTGGCACTTGGGGAGTCCAAGGCAGGATGATCACTTATGC 
CCAAGGAATTTGAGTACCAAGCCTGGGCAATATAACAAGGCCCTGTTTCTACAAAAACTT 
TAAACAATTAGCCAGGTGTGGTGGTGCGTGCCTGTGTCCAGCTACTCAGGAAGCTGAGGC 
AAGAGCTTGAGGCTACAGTGAGCTGTGTTCCACCATGGTGCTCCAGCCTGGGTGACAGGG 
CAAGACCCTGTCAAAAGAAAGGAAGAAAGAACGGAAGGAAAGAAGGAAAGAAACAAGGAG 
AG

Die Sequenzdaten werden aus der Datei alu.n gesammelt. Daher stimmt es mit unserer Datenbank überein.

Step 4 - Die BLAST-Software bietet viele Anwendungen zum Durchsuchen der Datenbank und wir verwenden Blastn. blastn application requires minimum of three arguments, db, query and out. db verweist auf die Datenbank gegen zu suchen; query ist die Reihenfolge zu entsprechen und outist die Datei zum Speichern der Ergebnisse. Führen Sie nun den folgenden Befehl aus, um diese einfache Abfrage auszuführen:

blastn -db alun -query search.fsa -out results.xml -outfmt 5

Wenn Sie den obigen Befehl ausführen, wird gesucht und in der Ausgabe ausgegeben results.xml Datei wie unten angegeben (teilweise Daten) -

<?xml version = "1.0"?> 
<!DOCTYPE BlastOutput PUBLIC "-//NCBI//NCBI BlastOutput/EN" 
   "http://www.ncbi.nlm.nih.gov/dtd/NCBI_BlastOutput.dtd">
<BlastOutput> 
   <BlastOutput_program>blastn</BlastOutput_program> 
   <BlastOutput_version>BLASTN 2.7.1+</BlastOutput_version> 
   <BlastOutput_reference>Zheng Zhang, Scott Schwartz, Lukas Wagner, and Webb 
      Miller (2000), "A greedy algorithm for aligning DNA sequences", J 
      Comput Biol 2000; 7(1-2):203-14.
   </BlastOutput_reference> 
   
   <BlastOutput_db>alun</BlastOutput_db> 
   <BlastOutput_query-ID>Query_1</BlastOutput_query-ID> 
   <BlastOutput_query-def>gnl|alu|Z15030_HSAL001056 (Alu-J)</BlastOutput_query-def>
   <BlastOutput_query-len>292</BlastOutput_query-len> 
   <BlastOutput_param> 
      <Parameters> 
         <Parameters_expect>10</Parameters_expect> 
         <Parameters_sc-match>1</Parameters_sc-match> 
         <Parameters_sc-mismatch>-2</Parameters_sc-mismatch> 
         <Parameters_gap-open>0</Parameters_gap-open> 
         <Parameters_gap-extend>0</Parameters_gap-extend> 
         <Parameters_filter>L;m;</Parameters_filter> 
      </Parameters> 
   </BlastOutput_param> 
   <BlastOutput_iterations> 
      <Iteration> 
         <Iteration_iter-num>1</Iteration_iter-num><Iteration_query-ID>Query_1</Iteration_query-ID> 
         <Iteration_query-def>gnl|alu|Z15030_HSAL001056 (Alu-J)</Iteration_query-def> 
         <Iteration_query-len>292</Iteration_query-len> 
         <Iteration_hits> 
         <Hit>
            <Hit_num>1</Hit_num> 
            <Hit_id>gnl|alu|Z15030_HSAL001056</Hit_id> 
            <Hit_def>(Alu-J)</Hit_def> 
            <Hit_accession>Z15030_HSAL001056</Hit_accession> 
            <Hit_len>292</Hit_len>
            <Hit_hsps> 
               <Hsp>
                 <Hsp_num>1</Hsp_num> 
                  <Hsp_bit-score>540.342</Hsp_bit-score> 
                  <Hsp_score>292</Hsp_score>
                  <Hsp_evalue>4.55414e-156</Hsp_evalue> 
                  <Hsp_query-from>1</Hsp_query-from>
                  <Hsp_query-to>292</Hsp_query-to> 
                  <Hsp_hit-from>1</Hsp_hit-from> 
                  <Hsp_hit-to>292</Hsp_hit-to> 
                  <Hsp_query-frame>1</Hsp_query-frame>
                  <Hsp_hit-frame>1</Hsp_hit-frame>
                  <Hsp_identity>292</Hsp_identity>
                  <Hsp_positive>292</Hsp_positive> 
                  <Hsp_gaps>0</Hsp_gaps> 
                  <Hsp_align-len>292</Hsp_align-len>
                  
                  <Hsp_qseq>
                     AGGCTGGCACTGTGGCTCATGCTGAAATCCCAGCACGGCGGAGGACGGCGGAAGATTGCTTGAGCCTAGGAGTTTG
                     CGACCAGCCTGGGTGACATAGGGAGATGCCTGTCTCTACGCAAAAGAAAAAAAAAATAGCTCTGCTGGTGGTGCATG
                     CCTATAGTCTCAGCTATCAGGAGGCTGGGACAGGAGGATCACTTGGGCCCGGGAGTTGAGGCTGTGGTGAGCC
                     ACGATCACACCACTGCACTCCAGCCTGGGTGACAGAGCAAGACCCTGTCTCAAAACAAACAAATAA
                  </Hsp_qseq> 

                  <Hsp_hseq>
                     AGGCTGGCACTGTGGCTCATGCTGAAATCCCAGCACGGCGGAGGACGGCGGAAGATTGCTTGAGCCTAGGA
                     GTTTGCGACCAGCCTGGGTGACATAGGGAGATGCCTGTCTCTACGCAAAAGAAAAAAAAAATAGCTCTGCT
                     GGTGGTGCATGCCTATAGTCTCAGCTATCAGGAGGCTGGGACAGGAGGATCACTTGGGCCCGGGAGTTGAGG
                     CTGTGGTGAGCCACGATCACACCACTGCACTCCAGCCTGGGTGACAGAGCAAGACCCTGTCTCAAAACAAAC
                     AAATAA
                  </Hsp_hseq>

                  <Hsp_midline>
                     |||||||||||||||||||||||||||||||||||||||||||||||||||||
                     |||||||||||||||||||||||||||||||||||||||||||||||||||||
                     |||||||||||||||||||||||||||||||||||||||||||||||||||||
                     |||||||||||||||||||||||||||||||||||||||||||||||||||||
                     |||||||||||||||||||||||||||||||||||||||||||||||||||||
                     |||||||||||||||||||||||||||
                  </Hsp_midline>
               </Hsp> 
            </Hit_hsps>
         </Hit>
         ......................... 
         ......................... 
         ......................... 
         </Iteration_hits> 
         <Iteration_stat> 
            <Statistics> 
               <Statistics_db-num>327</Statistics_db-num> 
               <Statistics_db-len>80506</Statistics_db-len> 
               <Statistics_hsp-lenv16</Statistics_hsp-len> 
               <Statistics_eff-space>21528364</Statistics_eff-space> 
               <Statistics_kappa>0.46</Statistics_kappa> 
               <Statistics_lambda>1.28</Statistics_lambda> 
               <Statistics_entropy>0.85</Statistics_entropy>
            </Statistics>
         </Iteration_stat>
      </Iteration> 
   </BlastOutput_iterations>
</BlastOutput>

Der obige Befehl kann im Python mit dem folgenden Code ausgeführt werden:

>>> from Bio.Blast.Applications import NcbiblastnCommandline 
>>> blastn_cline = NcbiblastnCommandline(query = "search.fasta", db = "alun", 
outfmt = 5, out = "results.xml") 
>>> stdout, stderr = blastn_cline()

Hier ist der erste ein Handle für die Explosionsausgabe und der zweite ist die mögliche Fehlerausgabe, die durch den Explosionsbefehl erzeugt wird.

Da wir die Ausgabedatei als Befehlszeilenargument (out = "results.xml") angegeben und das Ausgabeformat als XML (outfmt = 5) festgelegt haben, wird die Ausgabedatei im aktuellen Arbeitsverzeichnis gespeichert.

Analyse des BLAST-Ergebnisses

Im Allgemeinen wird die BLAST-Ausgabe mithilfe des NCBIXML-Moduls als XML-Format analysiert. Dazu müssen wir das folgende Modul importieren:

>>> from Bio.Blast import NCBIXML

Jetzt, open the file directly using python open method und use NCBIXML parse method wie unten angegeben -

>>> E_VALUE_THRESH = 1e-20 
>>> for record in NCBIXML.parse(open("results.xml")): 
>>>     if record.alignments: 
>>>        print("\n") 
>>>        print("query: %s" % record.query[:100]) 
>>>        for align in record.alignments: 
>>>           for hsp in align.hsps: 
>>>              if hsp.expect < E_VALUE_THRESH: 
>>>                 print("match: %s " % align.title[:100])

Dies erzeugt eine Ausgabe wie folgt:

query: gnl|alu|Z15030_HSAL001056 (Alu-J) 
match: gnl|alu|Z15030_HSAL001056 (Alu-J) 
match: gnl|alu|L12964_HSAL003860 (Alu-J) 
match: gnl|alu|L13042_HSAL003863 (Alu-FLA?) 
match: gnl|alu|M86249_HSAL001462 (Alu-FLA?) 
match: gnl|alu|M29484_HSAL002265 (Alu-J) 

query: gnl|alu|D00596_HSAL003180 (Alu-Sx) 
match: gnl|alu|D00596_HSAL003180 (Alu-Sx) 
match: gnl|alu|J03071_HSAL001860 (Alu-J) 
match: gnl|alu|X72409_HSAL005025 (Alu-Sx) 

query: gnl|alu|X55502_HSAL000745 (Alu-J) 
match: gnl|alu|X55502_HSAL000745 (Alu-J)