📌  相关文章
📜  Biopython-BLAST概述

📅  最后修改于: 2020-11-09 05:02:05             🧑  作者: Mango

BLAST代表基本局部路线搜索工具。它找到生物学序列之间的相似区域。 Biopython提供了Bio.Blast模块来处理NCBI BLAST操作。您可以在本地连接或Internet连接上运行BLAST。

让我们在下一节中简要了解这两个连接-

通过互联网运行

Biopython提供了Bio.Blast.NCBIWWW模块来调用BLAST的在线版本。为此，我们需要导入以下模块-

>>> from Bio.Blast import NCBIWWW

NCBIWW模块提供qblast函数来查询BLAST在线版本https://blast.ncbi.nlm.nih.gov/Blast.cgi 。 qblast支持在线版本支持的所有参数。

要获得有关此模块的任何帮助，请使用以下命令并了解功能-

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

通常，qblast函数的参数基本上类似于您可以在BLAST网页上设置的不同参数。这使qblast函数易于理解，并减少了使用它的学习曲线。

连接和搜索

要了解连接和搜索BLAST在线版本的过程，让我们通过Biopython对在线BLAST服务器进行简单的序列搜索(在我们的本地序列文件中可用)。

步骤1-在Biopython目录中创建一个名为blast_example.fasta的文件，并提供以下序列信息作为输入

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

>sequence B ggtaagtcctctagtacaaacacccccaatattgtgatataattaaaattatattca
tattctgttgccagaaaaaacacttttaggctatattagagccatcttctttgaagcgttgtc

步骤2-导入NCBIWWW模块。

>>> from Bio.Blast import NCBIWWW

步骤3-使用Python IO模块打开序列文件blast_example.fasta 。

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

步骤4-现在，调用传递序列数据作为主要参数的qblast函数。另一个参数表示数据库(nt)和内部程序(blastn)。

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

blast_results保存我们的搜索结果。可以将其保存到文件中以供以后使用，也可以对其进行解析以获取详细信息。我们将在接下来的部分中学习如何做。

第5步-可以使用Seq对象完成相同的功能，而不是使用整个fasta文件，如下所示-

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

现在，调用传递给Seq对象record.seq作为主要参数的qblast函数。

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

BLAST将自动为您的序列分配一个标识符。

第6步-result_handle对象将具有完整的结果，可以保存到文件中以备后用。

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

我们将在后面的部分中了解如何解析结果文件。

运行独立的BLAST

本节说明如何在本地系统中运行BLAST。如果您在本地系统中运行BLAST，它可能会更快，并且还允许您创建自己的数据库以针对序列进行搜索。

连接BLAST

通常，不建议在本地运行BLAST，这是因为其体积较大，运行软件需要付出额外的努力以及所涉及的成本。在线BLAST足以满足基本和高级目的。当然，有时您可能需要在本地安装它。

考虑到您正在进行频繁的在线搜索，这可能需要大量的时间并且网络流量很大，如果您有专有序列数据或IP相关问题，则建议在本地安装。

为此，我们需要遵循以下步骤-

步骤1-使用给定的链接下载并安装最新的blast二进制文件-ftp://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/

步骤2-使用以下链接下载并解压缩最新的必要数据库-ftp://ftp.ncbi.nlm.nih.gov/blast/db/

BLAST软件在其站点中提供了许多数据库。让我们从blast数据库站点下载alu.n.gz文件，并将其解压缩到alu文件夹中。该文件为FASTA格式。要在blast应用程序中使用此文件，我们需要首先将文件从FASTA格式转换为blast数据库格式。 BLAST提供了makeblastdb应用程序来执行此转换。

使用以下代码片段-

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

运行上面的代码将解析输入文件alu.n并将BLAST数据库创建为多个文件alun.nsq，alun.nsi等。现在，我们可以查询该数据库以查找序列。

我们已经安装在我们的本地服务器的BLAST，也有样本BLAST数据库，阿伦反对它的查询。

步骤3-让我们创建一个样本序列文件来查询数据库。创建文件search.fsa并将以下数据放入其中。

>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

序列数据是从alu.n文件中收集的；因此，它与我们的数据库匹配。

步骤4 -BLAST软件提供了许多应用程序来搜索数据库，我们使用blastn。 blastn应用程序至少需要三个参数db，query和out。 db是指要搜索的数据库； query是要匹配的序列， out是要存储结果的文件。现在，运行以下命令来执行此简单查询-

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

运行上面的命令将搜索并在result.xml文件中提供输出，如下所示(部分数据)-

 

 
   blastn 
   BLASTN 2.7.1+ 
   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.
    
   
   alun 
   Query_1 
   gnl|alu|Z15030_HSAL001056 (Alu-J)
   292 
    
       
         10 
         1 
         -2 
         0 
         0 
         L;m; 
       
    
    
       
         1Query_1 
         gnl|alu|Z15030_HSAL001056 (Alu-J) 
         292 
          
         
            1 
            gnl|alu|Z15030_HSAL001056 
            (Alu-J) 
            Z15030_HSAL001056 
            292
             
               
                 1 
                  540.342 
                  292
                  4.55414e-156 
                  1
                  292 
                  1 
                  292 
                  1
                  1
                  292
                  292 
                  0 
                  292
                  
                  
                     AGGCTGGCACTGTGGCTCATGCTGAAATCCCAGCACGGCGGAGGACGGCGGAAGATTGCTTGAGCCTAGGAGTTTG
                     CGACCAGCCTGGGTGACATAGGGAGATGCCTGTCTCTACGCAAAAGAAAAAAAAAATAGCTCTGCTGGTGGTGCATG
                     CCTATAGTCTCAGCTATCAGGAGGCTGGGACAGGAGGATCACTTGGGCCCGGGAGTTGAGGCTGTGGTGAGCC
                     ACGATCACACCACTGCACTCCAGCCTGGGTGACAGAGCAAGACCCTGTCTCAAAACAAACAAATAA
                   

                  
                     AGGCTGGCACTGTGGCTCATGCTGAAATCCCAGCACGGCGGAGGACGGCGGAAGATTGCTTGAGCCTAGGA
                     GTTTGCGACCAGCCTGGGTGACATAGGGAGATGCCTGTCTCTACGCAAAAGAAAAAAAAAATAGCTCTGCT
                     GGTGGTGCATGCCTATAGTCTCAGCTATCAGGAGGCTGGGACAGGAGGATCACTTGGGCCCGGGAGTTGAGG
                     CTGTGGTGAGCCACGATCACACCACTGCACTCCAGCCTGGGTGACAGAGCAAGACCCTGTCTCAAAACAAAC
                     AAATAA
                  

                  
                     |||||||||||||||||||||||||||||||||||||||||||||||||||||
                     |||||||||||||||||||||||||||||||||||||||||||||||||||||
                     |||||||||||||||||||||||||||||||||||||||||||||||||||||
                     |||||||||||||||||||||||||||||||||||||||||||||||||||||
                     |||||||||||||||||||||||||||||||||||||||||||||||||||||
                     |||||||||||||||||||||||||||
                  
                
            
         
         ......................... 
         ......................... 
         ......................... 
          
          
             
               327 
               80506 
                
               21528364 
               0.46 
               1.28 
               0.85

上面的命令可以使用以下代码在Python内部运行-

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

在这里，第一个是爆炸输出的句柄，第二个是爆炸命令生成的可能的错误输出。

由于我们已将输出文件提供为命令行参数(out =“ results.xml”)并将输出格式设置为XML(outfmt = 5)，因此输出文件将保存在当前工作目录中。

解析BLAST结果

通常，使用NCBIXML模块将BLAST输出解析为XML格式。为此，我们需要导入以下模块-

>>> from Bio.Blast import NCBIXML

现在，直接使用Python open方法打开文件，并使用NCBIXML parse方法，如下所示-

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

这将产生如下输出：

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)