Biopython - Tổng quan về BLAST

BLAST là viết tắt của Basic Local Alignment Search Tool. Nó tìm thấy các vùng tương đồng giữa các trình tự sinh học. Biopython cung cấp mô-đun Bio.Blast để đối phó với hoạt động NCBI BLAST. Bạn có thể chạy BLAST trong kết nối cục bộ hoặc qua kết nối Internet.

Hãy để chúng tôi hiểu ngắn gọn về hai kết nối này trong phần sau:

Chạy qua Internet

Biopython cung cấp mô-đun Bio.Blast.NCBIWWW để gọi phiên bản trực tuyến của BLAST. Để thực hiện việc này, chúng ta cần nhập mô-đun sau:

>>> from Bio.Blast import NCBIWWW

Mô-đun NCBIWW cung cấp chức năng qblast để truy vấn phiên bản trực tuyến BLAST, https://blast.ncbi.nlm.nih.gov/Blast.cgi. qblast hỗ trợ tất cả các tham số được hỗ trợ bởi phiên bản trực tuyến.

Để nhận bất kỳ trợ giúp nào về mô-đun này, hãy sử dụng lệnh dưới đây và hiểu các tính năng -

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

Thông thường, các đối số của hàm qblast về cơ bản tương tự với các tham số khác nhau mà bạn có thể đặt trên trang web BLAST. Điều này làm cho hàm qblast dễ hiểu cũng như giảm bớt thời gian học để sử dụng nó.

Kết nối và Tìm kiếm

Để hiểu quá trình kết nối và tìm kiếm phiên bản BLAST trực tuyến, hãy để chúng tôi thực hiện tìm kiếm theo trình tự đơn giản (có sẵn trong tệp trình tự cục bộ của chúng tôi) trên máy chủ BLAST trực tuyến thông qua Biopython.

Step 1 - Tạo một tệp có tên blast_example.fasta trong thư mục Biopython và cung cấp thông tin trình tự bên dưới làm đầu vào

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

>sequence B ggtaagtcctctagtacaaacacccccaatattgtgatataattaaaattatattca
tattctgttgccagaaaaaacacttttaggctatattagagccatcttctttgaagcgttgtc

Step 2 - Nhập mô-đun NCBIWWW.

>>> from Bio.Blast import NCBIWWW

Step 3 - Mở tệp trình tự, blast_example.fasta sử dụng mô-đun IO python.

>>> 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- Bây giờ, hãy gọi hàm qblast truyền dữ liệu chuỗi làm tham số chính. Tham số khác đại diện cho cơ sở dữ liệu (nt) và chương trình bên trong (blastn).

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

blast_resultsgiữ kết quả tìm kiếm của chúng tôi. Nó có thể được lưu vào một tệp để sử dụng sau này và cũng có thể được phân tích cú pháp để lấy thông tin chi tiết. Chúng ta sẽ tìm hiểu cách thực hiện trong phần tới.

Step 5 - Chức năng tương tự cũng có thể được thực hiện bằng cách sử dụng đối tượng Seq thay vì sử dụng toàn bộ tệp fasta như hình dưới đây -

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

Bây giờ, hãy gọi hàm qblast truyền đối tượng Seq, record.seq làm tham số chính.

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

BLAST sẽ tự động chỉ định một số nhận dạng cho chuỗi của bạn.

Step 6 - Đối tượng result_handle sẽ có toàn bộ kết quả và có thể được lưu thành tệp để sử dụng sau này.

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

Chúng ta sẽ xem cách phân tích cú pháp tệp kết quả trong phần sau.

Chạy BLAST độc lập

Phần này giải thích về cách chạy BLAST trong hệ thống cục bộ. Nếu bạn chạy BLAST trong hệ thống cục bộ, nó có thể nhanh hơn và cũng cho phép bạn tạo cơ sở dữ liệu của riêng mình để tìm kiếm theo trình tự.

Kết nối BLAST

Nói chung, chạy BLAST cục bộ không được khuyến khích do kích thước lớn, cần thêm nỗ lực để chạy phần mềm và chi phí liên quan. BLAST trực tuyến là đủ cho các mục đích cơ bản và nâng cao. Tất nhiên, đôi khi bạn có thể được yêu cầu cài đặt nó cục bộ.

Hãy xem xét bạn đang thực hiện các tìm kiếm trực tuyến thường xuyên có thể đòi hỏi nhiều thời gian và khối lượng mạng cao và nếu bạn có dữ liệu trình tự độc quyền hoặc các vấn đề liên quan đến IP, thì nên cài đặt cục bộ.

Để thực hiện việc này, chúng ta cần làm theo các bước sau:

Step 1- Tải xuống và cài đặt bản nhị phân blast mới nhất bằng liên kết đã cho - ftp://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/

Step 2- Tải xuống và giải nén cơ sở dữ liệu cần thiết và mới nhất bằng liên kết dưới đây - ftp://ftp.ncbi.nlm.nih.gov/blast/db/

Phần mềm BLAST cung cấp rất nhiều cơ sở dữ liệu trong trang web của họ. Hãy để chúng tôi tải xuống tệp alu.n.gz từ trang cơ sở dữ liệu blast và giải nén nó vào thư mục alu. Tệp này ở định dạng FASTA. Để sử dụng tệp này trong ứng dụng blast, trước tiên chúng ta cần chuyển đổi tệp từ định dạng FASTA sang định dạng cơ sở dữ liệu blast. BLAST cung cấp ứng dụng makeblastdb để thực hiện việc chuyển đổi này.

Sử dụng đoạn mã dưới đây -

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

Chạy đoạn mã trên sẽ phân tích cú pháp tệp đầu vào, alu.n và tạo cơ sở dữ liệu BLAST dưới dạng nhiều tệp alun.nsq, alun.nsi, v.v. Bây giờ, chúng ta có thể truy vấn cơ sở dữ liệu này để tìm chuỗi.

Chúng tôi đã cài đặt BLAST trong máy chủ cục bộ của mình và cũng có cơ sở dữ liệu BLAST mẫu, alun để truy vấn chống lại nó.

Step 3- Chúng ta hãy tạo một tệp trình tự mẫu để truy vấn cơ sở dữ liệu. Tạo một tệp search.fsa và đưa dữ liệu bên dưới vào đó.

>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

Dữ liệu trình tự được thu thập từ tệp alu.n; do đó, nó phù hợp với cơ sở dữ liệu của chúng tôi.

Step 4 - Phần mềm BLAST cung cấp nhiều ứng dụng để tìm kiếm cơ sở dữ liệu và chúng tôi sử dụng blastn. blastn application requires minimum of three arguments, db, query and out. db đề cập đến cơ sở dữ liệu chống lại để tìm kiếm; query là trình tự để khớp và outlà tệp để lưu trữ kết quả. Bây giờ, hãy chạy lệnh dưới đây để thực hiện truy vấn đơn giản này -

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

Chạy lệnh trên sẽ tìm kiếm và đưa ra kết quả trong results.xml tệp như được cung cấp bên dưới (một phần dữ liệu) -

<?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>

Lệnh trên có thể được chạy bên trong python bằng đoạn mã dưới đây:

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

Ở đây, cái đầu tiên là một xử lý cho đầu ra blast và cái thứ hai là đầu ra lỗi có thể được tạo ra bởi lệnh blast.

Vì chúng tôi đã cung cấp tệp đầu ra dưới dạng đối số dòng lệnh (out = “results.xml”) và đặt định dạng đầu ra là XML (outfmt = 5), tệp đầu ra sẽ được lưu trong thư mục làm việc hiện tại.

Phân tích cú pháp kết quả BLAST

Nói chung, đầu ra BLAST được phân tích cú pháp dưới dạng định dạng XML bằng cách sử dụng mô-đun NCBIXML. Để thực hiện việc này, chúng ta cần nhập mô-đun sau:

>>> from Bio.Blast import NCBIXML

Hiện nay, open the file directly using python open method và use NCBIXML parse method như dưới đây -

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

Điều này sẽ tạo ra một đầu ra như sau:

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)