<jupyter_start><jupyter_text>Capítulo 7 - BLAST Todo mundo adora o BLAST, certo? Quer dizer, como pode ser tão fácil comparar sua sequências com milhares de outras sequências conhecidas no mundo? Mas claro, esse capítulo não é para dizer o quão legal o BLAST é (acredito que já sabemos disso). Esse capítulo é sobre o problema do BLAST, pode ser muito difícil lidar com um grande volume de dados gerados por longas execuções e automatizar essas execuções do BLAST em geral.Felizmente, o pessoal do Biopython sabe disso muito bem e por isso desenvolveram ferramentas capazes de lidar com o BLAST e deixar a sua vida como bioinformata muito mais simples. Esse capítulo detalha como usar essas ferramentas e as coisas úteis que você pode realizar com as mesmas.Lidar com BLAST pode ser dividido em duas etapas, ambas podem ser feitas com o Biopython. Primeiramente, executar o BLAST para a(s) sua(s) sequência(s) query e obter muitas saídas. Segundo, analisar a saída do BLAST em Python para análise porterior.A sua introdução ao BLAST provavelmente foi pelo serviço online do NCBI (a minha também foi assim). Na verdade, existem muitos caminhos que você pode executar o BLAST, que podem ser categorizados em diversos caminhos. A distinção mais importante é executar o BLAST localmente (no seu computador) e executar o BLAST remotamente (em outro computador, geralmente nos servidores do NCBI). Vamos começar esse capítulo acessando os servidores do NCBI atavés de um script em Python.*NOTA*: O Capítulo 8 (https://github.com/GTL98/biopython_br/tree/main/Cap%208%20-%20BLAST%20e%20outras%20ferramentas%20de%20pesquisa%20de%20sequ%C3%AAncias) descreve o módulo experimental `Bio.SearchIO`. O pessoal do Biopython quer substituir o módulo `Bio.Blast` mais antigo, pois ele fornece uma estrtura mais geral que trata também de outras ferramentas de busca de sequências relacionadas. Contudo, enquanto é declarado estável, para a produção dos códigos continue usando o módulo `Bio.Blast` para trabalhar com o NCBI BLAST.  7.1 Executando o BLAST pela Internet Usaremos a função **qblast** do módulo `Bio.Blast.NCBIWWW` para usar a versão online do BLAST. Essa função possui três argumentos não opcionais:- O programa do BLAST que será usado deve ser passado como uma string em minúsculo. As opções e descrições dos progrmadas disponíveos podem ser encontrados em: https://blast.ncbi.nlm.nih.gov/Blast.cgi. Atualmente, **qblast** funciona com *blastn*, *blastx*, *tblast* e *tblasx*.- O segundo argumento especifíca o banco de dados que deve ser pesquisado. Novamente, as informações estão dispoíveis em: https://ftp.ncbi.nlm.nih.gov/pub/factsheets/HowTo_BLASTGuide.pdf.- O último argumento é uma string que contém a sequência query. Pode ser a própria sequência em si, a sequência em formato FASTA ou o número de identificação no NCBI.A função **qblast** também aceita outros argumentos de opção que são basicamente análogos aos diferentes parâmetros que podemos definir na página da internet do BLAST. Vamos ver apenas alguns deles:- O argumento **url_base** configura a URL da internet em que o BLAST deve exectar. Por padrão, é conectado ao NCBI, mas pode o NCBI BLAST pode ser executado na nuvem. Em caso de dúvida, veja a documentação da **qblast**.- A função **qblast** pode retornar um resultado BLAST em diversos formatos que podem ser escolhidos pelo argumento **format_type**: *HTML*, *Text*, *ASN.1* ou *XML*. O valor padrão é *XML*, pois esse é o formato esperado pelo analisador, descrito melhor na **Seção 7.3**.- O argumento **expect** configura a expectativa ou o limite de e-value.Para mais informações, pesquise a documentação do NCBI ou use a função **help**:<jupyter_code>from Bio.Blast import NCBIWWW
help(NCBIWWW)<jupyter_output>Help on module Bio.Blast.NCBIWWW in Bio.Blast:

NAME
    Bio.Blast.NCBIWWW - Code to invoke the NCBI BLAST server over the internet.

DESCRIPTION
    This module provides code to work with the WWW version of BLAST
    provided by the NCBI. https://blast.ncbi.nlm.nih.gov/

FUNCTIONS
    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, short_query=None, alignments=500, alignment_view=None, descriptions=500, entrez_links_new_window=None, expect_[...]<jupyter_text>Observe que as configurações padrão o website NCBI BLASt não exatamente iguais com os padrões do QBLAST. Se você obter resultados diferentes, você precisa checar os parâmetros (e.g., o valor limite de e-value e o valor dos gaps).Por exemplo, se você tem uma sequência de nucleotídeos e quer procurar em um banco de dados de nucleotídeo (*nt*) usando o BLASTN e somente com o número de identificação (BE037100.1) ou o número GI (8332116), você pode fazer assim:<jupyter_code># Farei somente essa célula porque o processo é demorado demais para todas as células e no final é o mesmo resultado.
# Somente essa célula fará o que precisamos, mas deixarei as outras como exemplo caso precisem.
from Bio.Blast import NCBIWWW
resultado = NCBIWWW.qblast('blastn', 'nt', 'BE037100.1')<jupyter_output><empty_output><jupyter_text>De modo alternativo, se você tiver a sequência em um arquivo FASTA em seu computados, você apenas precisa abrí-lo e ler o registro como uma string para usar como o argumento de query na função **qblast**:<jupyter_code>from Bio.Blast import NCBIWWW
fasta_string = open('m_cold.fasta').read()
resultado = NCBIWWW.qblast('blastn', 'nt', fasta_string)<jupyter_output><empty_output><jupyter_text>Você pode também ler o arquivo FASTA como um objeto **SeqRecord** e, em seguida, fornecer a sequência desse objeto:<jupyter_code>from Bio import SeqIO
from Bio.Blast import NCBIWWW
registro = SeqIO.read('m_cold.fasta', 'fasta')
resultado = NCBIWWW.qblast('blastn', 'nt', registro.seq)<jupyter_output><empty_output><jupyter_text>Fornecendo apenas a sequência indica que BLAST atribuirá um identificador para a sua sequência automaticamente. Você pode prefirir usar o formato dos objetos **SeqRecord** para fazer uma string FASTA, que incluirá o idetificador existente:<jupyter_code>from Bio import SeqIO
from Bio.Blast import NCBIWWW
registro = SeqIO.read('m_cold.fasta', 'fasta')
resultado = NCBIWWW.qblast('blastn', 'nt', registro.format('fasta'))<jupyter_output><empty_output><jupyter_text>Essa abordagem possui mais sentido se você tiver sequência(s) em formatos de arquivos não FASTA que você pode extrair usando o módulo `Bio.SeqIO` (veja mais no Capítulo 5: https://github.com/GTL98/biopython_br/tree/main/Cap%205%20-%20Sequ%C3%AAncias%20Input_Output).Independente dos argumentos que você use na função **qblast**, você deve ter como resultado em um objeto de identificação (por padrão no formato XML). O próximo passo seria a análise da saída XML em objetos Python que representam os resultados de pesquisa (veja a **Seção 7.3**), mas você pode querer salvar uma cópia da saída em um arquivo. Isso é muito útil ao depurar o código que extrai informações dos resultados do BLAST, porque executar novamente uma pesquisa online é consideravelmente lenta.Precisamos ter um pouco de contato, pois podemos usar **resultado.read()** para ler a saída BLAST apenas uma vez. Chamar **resultado.read()** de novo retorna uma string vazia:<jupyter_code>with open('meu_blast.xml', 'w') as identificador_saida:
    identificador_saida.write(resultado.read())
resultado.close()<jupyter_output><empty_output><jupyter_text>Depois de feito isso, os resultados estão salvos no arquivo **meu_blast.xml** (que está no mesmo repositório deste arquivo) e o identificador original teve todos os seus dados extraídos, portanto podemos fechá-lo. Contudo, a função **parse** do BLAST parser (veja a **Seção 7.3**) usa um objeto semelhante a um identificador de arquivo, então podemos abrir o arquivo salvo para a entrada:<jupyter_code>resultado = open('meu_blast.xml')<jupyter_output><empty_output><jupyter_text>Agora que temos os resultados do BLAST em um identificador novamente, estamos prontos para fazer algo com eles, então podemos ir direto para a seção de análise (**Seção 7.3**). Você pode pular para essa seção agora se quiser.  7.2 Executando o BLAST localmente  7.2.1 Introdução Executar o BLAST localmente possui ao menos duas vantagens:1. BLAST local pode ser mais rápido do que o BLAST pela internet e;2. BLAST local permite que você crie seu próprio banco de dados para pesquisar as sequências.Lidar com seqências privadas ou ainda não publicadas pode ser outro motivo para você usar o BLAST local. Você pode não ter permisão para redistribuir as sequências, então submetê-las ao NCBI como um query do BLAST pode não ser uma opção.Infelizmente, há também as suas desvantagens; instalar o que precisa e configurar corretamente exige algum esforço:1. BLAST local necessita de ferramentas de linhas de comando instaladas e;2. BLAST local requer que base de dados BLAST (grandes) sejam configurados (e potencialmente sejam mantidos atualizados).Para confundir ainda mais as coisas, há diversos pacores BLAST disponíveis e existe outras ferramentas que produzem arquivos de saída igual ao BLAST, como o BLAT.  7.2.2 NCBI BLAST+ autônomo O 'novo' NCBI BLAST+ (https://blast.ncbi.nlm.nih.gov/Blast.cgi?CMD=Web&PAGE_TYPE=BlastDocs&DOC_TYPE=Download) foi lançado em 2009. Ele veio no lugar do antigo NCBI 'legacy' BLAST.Essa seção mostrará brevemente como usar essas ferramentas com o Python. Se você já tiver lido ou testato os exemplos das ferramentas de alinhamento da **Seção 6.5** do Capítulo 6 (https://github.com/GTL98/biopython_br/tree/main/Cap%206%20-%20Objetos%20de%20alinhamento%20de%20sequ%C3%AAncia%20m%C3%BAltipla), tudo o que veremos deverá parecer bem simples. Primeiro, vamos construir uma string de linha de comando (do mesmo modo como você faria em uma linha de comando no prompt se estivesse executando o BLAST autônomo na mão). Podemos executar esse comando através do Python.Por exemplo, usando o arquivo FASTA com uma sequência de nucleotídeos, você pode querer executar um BLASTX (tradução) contra um banco de dados de um banco de dados de proteína não redundantes (NR). Assumindo que você tenha baixado a sequência e instalado o banco de dados NR, você deverá executar:- blastx -query opuntia.fasta -db nr -out opuntia.xml -evalue 0.001 -outfmt 5Isso provavelmente executará o BLASTX contra o banco de dados NR, usando o valor de corte e-value em 0.001 e produzindo um arquivo XML como saída (que podemos então analisar). Dependendo do computador, esse processo pode levar de 5 a 10 minutos para ficar pronto, por isso que é importante salvar todas as informações geradas em um arquivo.De dentro do Biopytonh, podemos usar o NCBI BLASTX com o módulo `Bio.Blast.Applications` para contruir a string de linha de comando e executá-la. O arquivo FASTA usado está no mesmo repositório desde arquivo:<jupyter_code>from Bio.Blast.Applications import NcbiblastxCommandline
help(NcbiblastxCommandline)
blastx_cline = NcbiblastxCommandline(query='opuntia.fasta', db='nr', evalue=0.001, outfmt=5, out='opuntia.xml')
blastx_cline
print(blastx_cline)<jupyter_output>blastx -out opuntia.xml -outfmt 5 -query opuntia.fasta -db nr -evalue 0.001<jupyter_text>Neste exemplo não deve ter nenhuma saída do BLASTX para os terminais, então **stdout** e **stderr** devem ser vazios. Mas mesmo assim, o arquivo **opuntia.xml** foi criado (nesse caso, vazio).Como você deve-se lembrar dos exemplos dos capítulos anteriores, o arquivo **opuntia.fasta** contém sete sequências de DNA, então a saída BLAST XML deve conter múltiplos resultados. Portanto, use o módulo `Bio.Blast.NCBIXML.parse()` para analisá-lo conforme descrito na **Seção 7.3**.*AVISO*: O pessoal do Biopython resolveu retirar essas linhas de comando futuramente. A documentação será atualizada para que seja possível a construção de uma linha de comando diretamente e invocá-la diretamente com módulo `subprocess`.  7.2.3 Outras versões do BLAST O NCBI BLAST+ (escrito em C++) foi lançado em 2009 como um substituto do original NCBI 'legacy' BLAST (escrito em C) que não possui mais atualizações. Houve muitas mudanças: a versão antiga tinha uma única ferramenta de linha de comando **blastlall** que cobria diferentes tipos de BLAST (que com o BLAST+ essa opção é separada em comandos específicos) e todas as opções das linhas de comando foram renomeadas. Os *wrappers* do Biopython para o NCBI 'legacy' BLAST foram descontinuados e retirados a partir da versão 1.73 do Biopython.Você também pode encontrar o Washington University BLAST (WU-BLAST: http://blast.wustl.edu/) e o seu sucessor Advanced Biocomputing BLAST (AB-BLAST, lançado em 2009; serviço pago: https://blast.advbiocomp.com/). Esses pacotes incluem as ferramentas de linha de comando **wu-blastall** e **ab-blastall**, que imitou o **blastall** do NCBI 'legacy' BLAST. O Biopython não possui atualmente *wrappers* para essas ferramentas, mas está apto em analisar qualquer formato de saída que seja compatível com o NCBI.  7.3 Analisando saída do BLAST Como mencionado acima, o BLAST pode gerar diversos formatos de arquivo de saída, como XML, HTML e texto. Originalmente, o Biopython possuia analisadores para saídas em texto e HTML, já que esses eram os únicos formatos de saída oferecidos na época. Infelizmente, a saída do BLAST nesses formatos continuou mudando, prejudicando totalmente os analisadores do Biopython. O analisador HTML BLAST foi removido, enquanto o analisador BLAST de texto está disponível apenas pelo módulo `Bio.SearchIO`. Usar essas saídas é por sua conta e risco, pode funcionar ou não dependendo de qual versão do BLAST você está usando.Acompanhar as mudanças no BLAST tornou-se uma tarefa quase impossível, especialmente com os usuários que usavam diferentes versões do BLAST, agora é recomendado que seja usado os analisadores com saída XML, que podem ser gerados pelas recentes versões do BLAST. A saída XML não é só mais estável que a saída de texto e HTML, mas também é muito mais fácil analisar automaticamente esse tipo de arquivo, tornando o Biopython muito mais estável.Você pode obter a saída BLAST no formato XML de várias maeiras. Para o analisador, pouco importa o caminho que a saída foi gerada, desde que esteja em XML está tudo certo:- Você pode usar o Biopython para executar o BLAST pela internet, como visto na **Seção 7.1**;- Você pode usar o Biopython para executar o BLAST localmente, como visto na **Seção 7.2**;- Você pode fazer uma pesquisa no BLAST pelo site do NCBI manualmente e salvar os resultados. Você deve escolher como formato de arquivo o XML e salvar a página final do BLAST que você quer em um arquivo e;- Você pode executar o BLAST localmente sem o Biopython e salvar a saída em um arquivo. De novo, você deve escolher o formato de arquivo XML para receber os resultados.O ponto importante aqui é que você não precisa sos scripts Biopython para buscar os dados para analisá-los. Fazendo as coisas de uma dessas maneiras, você precisa de um identificador para os resultados. Em Python, um identiicador é apenas uma forma de descrever a entrada para qualquer fonte de informação possa ser recuperada usando **read** ou **readlines** (veja a **Seção 24.1** do Capítulo 24: ***link do capítulo 24***).Se você seguiu o código acima para interagir com o BLAST por meio de um script, então você já tem o identificador **resultado**. Caso você tenha salvo o arquivo XML (como feito nesse arquivo), você pode abri-lo para leitura:<jupyter_code>resultado = open('my_blast.xml')<jupyter_output><empty_output><jupyter_text>Agora que temos o identificador, estamos prontos para analisar a saída. O código para a análise é um tanto quanto pequeno pequeno. Se você espera um único resultado do BLAST, isto é, se você usou uma única sequência query:<jupyter_code>from Bio.Blast import NCBIXML
registro_blast = NCBIXML.read(resultado)<jupyter_output><empty_output><jupyter_text>Ou se você tem muitos resultados, isto é, múltiplas sequências query:<jupyter_code>from Bio.Blast import NCBIXML
registros_blast = NCBIXML.parse(resultado)<jupyter_output><empty_output><jupyter_text>Assim como no `Bio.SeqIO` (Capítulo 5: https://github.com/GTL98/biopython_br/tree/main/Cap%205%20-%20Sequ%C3%AAncias%20Input_Output) e no `Bio.AlignIO` (Capítulo 6: https://github.com/GTL98/biopython_br/tree/main/Cap%206%20-%20Objetos%20de%20alinhamento%20de%20sequ%C3%AAncia%20m%C3%BAltipla), temos um par de funções de entrada, **read** e **parse**; onde **read** é usado quando você tem somente um objeto e **parse** quando você possui mais de um objeto para análise. Mas enquanto antes tinhamos objetos **SeqRecord** e **MultipleSeqAlignment**, agora temos objetos de registro BLAST.Se capaz de lidar com a situação onde o arquivo BLAST possui muita informação, contendo milhares de resultados, **NCBIXML.parse()** retorna um iterador. De um jeito simples: um iterador permite que você percorra a saída do BLAST, recuperando os registros do BLAST um por um para cada resultado de pesquisa do BLAST:<jupyter_code>resultado = open('my_blast.xml')
from Bio.Blast import NCBIXML
registros_blast = NCBIXML.parse(resultado)
registro_blast = next(registros_blast)
registro_blast = next(registros_blast)<jupyter_output><empty_output><jupyter_text>Ou você pode usar um loop **for**:<jupyter_code>for registro_blast in registros_blast:
    # Faça alguma coisa
    pass<jupyter_output><empty_output><jupyter_text>Observe que você pode percorrer os registros BLAST apenas uma vez. Geralmente, de cada registro BLAST você salvaria as informações de seu interesse. Se você quiser salvar e retornar aos registros BLAST, você pode converter o iterador em uma lista:<jupyter_code>registros_blast = list(registros_blast)<jupyter_output><empty_output><jupyter_text>Agora você pode acessar cada registro BLAST na lista com a indexação usual. Se o seu arquivo BLAST for consideravelmente grande, você pode pode ter problemas de memória tentanto salvar todos os registros em uma lista.Normalmente, você executará uma pesquisa do BLAST por vez. Logo, tudo o que você precisa é selecionar o primeiro (e somente ele) registro BLAST do **registros_blast**:<jupyter_code>from Bio.Blast import NCBIXML
resultado = open('my_blast.xml')
registros_blast = NCBIXML.parse(resultado)
registro_blast = next(registros_blast)<jupyter_output><empty_output><jupyter_text>Ou de um modo mais elegante:<jupyter_code>from Bio.Blast import NCBIXML
resultado = open('my_blast.xml')
registro_blast = NCBIXML.read(resultado)<jupyter_output><empty_output><jupyter_text>Mas agora eu imagino que você possa estar se perguntando o que há em um registro BLAST.  7.4 A classe de registro BLAST Um registro BLAST contém tudo o que você possa querer extrair de uma saída BLAST. Por agora, vamos ver uma exemplo de como pegar as informações de um relatório BLAST, mas se você quiser fazer algo em particular que não esteja descrito aqui, veja as informações da classe detalhadamente no código ou na documentação gerada automaticamente (as docstrings possuem muitas informações sobre o que está armazenado e o que cada pedaço da informação faz).Para continuar com nosso exemplo, vamos mostrar algumas informações resumidas sobre todas as informações presentes no relatório BLAST maiores que um determinado limite. O código abaixo faz isso:<jupyter_code>resultado = open('my_blast.xml')
from Bio.Blast import NCBIXML
registro_blast = NCBIXML.read(resultado)
LIMITE_E_VALUE = 0.04
for alinhamento in registro_blast.alignments:
    for hsp in alinhamento.hsps:
        if hsp.expect < LIMITE_E_VALUE:
            print('****Alinhamento****')
            print(f'Sequência: {alinhamento.title}')
            print(f'Tamanho: {alinhamento.length}')
            print(f'e-value: {hsp.expect}')
            print(hsp.query[0:75] + '...')
            print(hsp.match[0:75] + '...')
            print(hsp.sbjct[0:75] + '...')<jupyter_output>****Alinhamento****
Sequência: gi|262205317|ref|NR_030195.1| Homo sapiens microRNA 520b (MIR520B), microRNA
Tamanho: 61
e-value: 4.91307e-23
CCCTCTACAGGGAAGCGCTTTCTGTTGTCTGAAAGAAAAGAAAGTGCTTCCTTTTAGAGGG...
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||...
CCCTCTACAGGGAAGCGCTTTCTGTTGTCTGAAAGAAAAGAAAGTGCTTCCTTTTAGAGGG...
****Alinhamento****
Sequência: gi|301171311|ref|NR_035856.1| Pan troglodytes microRNA mir-520b (MIR520B), microRNA
Tamanho: 60
e-value: 1.71483e-22
CCTCTACAGGGAAGCGCTTTCTGTTGTCTGAAAGAAAAGAAAGTGCTTCCTTTTAGAGGG...
||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||...
CCTCTACAGGGAAGCGCTTTCTGTTGTCTGAAAGAAAAGAAAGTGCTTCCTTTTAGAGGG...
****Alinhamento****
Sequência: gi|270133242|ref|NR_032573.1| Macaca mulatta microRNA mir-519a (MIR519A), microRNA
Tamanho: 85
e-value: 2.54503e-20
CCCTCTACAGGGAAGCGCTTTCTGTTGTCTGAAAGAAAAGAAAGTGCTTCCTTTTAGAGGG...
||||||| ||||||||||||||||| |||||||||||||||||||||||||||||||||||...
CCCTCTAGAGGGAAGCGCTTTCTGTGGTCTGAAAGAAAAGAAAGTGC[...]
{'fraction_non_alphanumeric': 0.07747206310451686, 'fraction_numerical': 0.014930979434688704, 'mean_word_length': 5.505497861942578, 'pattern_counts': {'":': 0, '<': 67, '<?xml version=': 0, '>': 66, 'https://': 11, 'lorem ipsum': 0, 'www.': 0, 'xml': 12}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '17', 'dir': 'jupyter-structured-clean-dedup', 'id': '414077', 'license': 'MIT', 'n_tokens_mistral': 7682, 'n_tokens_neox': 7101, 'n_words': 2694, 'path': 'Cap 7- BLAST/Cap 7 - BLAST.ipynb', 'repo_name': 'GTL98/biopython_br', 'size': '21298'}