Este conteúdo faz parte do livro “Introdução à Programação para Bioinformática com Biopython“. Você pode adquirir a versão impressa desse livro aqui ou a versão para Kindle aqui. Para nos citar, consulte este link.
BLAST (Basic Local Alignment Search Tool ou Ferramenta de Busca de Alinhamento Local Básico) é uma ferramenta de comparação entre sequências primárias (nucleotídeos ou aminoácidos) que utiliza um algoritmo de programação dinâmica para alinhamento local. BLAST apresenta uma metodologia de grande eficiência para comparação par-a-par entre sequências, quanto comparação de uma sequência com uma grande base de dados.
A versão mais recente do BLAST foi introduzida em 2009 pelo NCBI. A suíte de aplicativos stand-alone BLAST+ foi construída com a linguagem C++, o que permitiu uma maior performance e uma melhoria nos recursos das aplicações. Assim, a suíte apresenta cinco aplicativos: blastn, blastp, blastx, tblastn e tblastx.
blastn | Compara uma sequência de nucleotídeos contra uma base de dados de sequência de nucleotídeos. |
blastp | Compara uma sequência de aminoácidos contra uma base de dados de sequência de proteínas. |
blastx | Compara uma sequência de nucleotídeos traduzida, em todos os quadros de leitura, contra uma base de dados de sequência de proteínas. |
tblastn | Compara uma sequência de proteína contra uma base de dados de sequência de nucleotídeos traduzida dinamicamente em todos os quadros de leitura. |
tblastx | Compara os seis quadros de leitura traduzidos de uma sequência de nucleotídeos contra os seis quadros de leitura traduzidos de uma base de dados de sequência de nucleotídeos. |
Em Bioinformática, BLAST é uma ferramenta tão difundida que poderíamos afirmar que BLAST é a ferramenta favorita de 11 em cada 10 bioinformatas (BLAST S2). Brincadeiras a parte, Biopython apresenta o pacote Bio.Blast que permite o uso do BLAST tanto com bases de dados locais quanto online.
Executando BLAST pela internet
BLAST pode ser executado online pelo Website oficial no NCBI <http://blast.ncbi.nlm.nih.gov>, entretanto o módulo Bio.Blast possui o método qblast que permite a execução de buscas através de scripts em Python.
A sintaxe do método é: “NCBIWWW.qblast(parâmetros)”. Sendo os parâmetros:
program | blastn, blastp, blastx, tblastn ou tblastx. |
database | Base de dados a qual será feita a busca (ex.: nr). |
sequence | Sequência que será buscada. |
descriptions | Número de descrições que serão exibidas (valor padrão: 500). |
alignments | Número de alinhamentos que serão exibidos (valor padrão: 500). |
expect | Valor de cutoff esperado (valor padrão: 10.0). |
matrix_name | Especificar uma matriz de pontuação (PAM30, PAM70, BLOSUM80, BLOSUM45). |
filter | Valor “nome” desativa filtragem. |
format_type | “HTML”, “Text”, “ASN.1” ou “XML” (valor padrão: “XML”). |
hitlist_size | Número de resultados exibidos (valor padrão: 50). |
megablast | TRUE/FALSE: Ativa o algoritmo MEGABLAST (exclusivo para blastn). |
service | plain, psi, phi, rpsblast, megablast. |
Para testarmos as buscas na base de dados nt do NCBI com BLAST através do método qblast( ), primeiro crie um arquivo chamado exemplo.fasta.
>beta-glucosidase_parcial_mRNA
CGTGTGCTTGCGCTACTGCTTGGTCTCACTGTGTGAA
CAACTGCTGCCGCCGCCTGCCGGCGGTGTGCTC
Agora vamos criar um script chamado blast.py (script c11_s1.py no GitHub).
from Bio.Blast import NCBIWWW
from Bio import SeqIO
arquivo = SeqIO.read("exemplo.fasta", format="fasta")
print("Iniciando busca no NCBIWWW..." )
resultado = NCBIWWW.qblast("blastn", "nt", \
arquivo.seq, format_type="Text")
print(resultado.read())
print("Busca concluida." )
Se você obteve ao final da execução a mensagem “Busca concluida”, sua consulta foi realizada com sucesso. Se o script travou com a mensagem “Iniciando busca no NCBIWWW…”, seu script provavelmente obteve algum problema ao conectar-se a banco de dados do NCBI. Recomendamos que aguarde alguns minutos. Se ainda não obtiver nenhum resultado tente interromper a execução do script (comando CONTROL + C) e executá-lo novamente (se estiver utilizando o Sublime Text pode ser necessário finalizar o aplicativo e abri-lo novamente). Lembre-se que é necessária conexão com a internet. Agora vamos entender melhor o que foi realizado nesse script.
No início do código, importamos o módulo NCBIWWW do pacote Bio.Blast e também o pacote SeqIO. O arquivo FASTA foi lido e aplicado à variável arquivo. O método qblast( ) recebeu como parâmetros o programa blastn, a base de dados nt, a sequência (seq) armazenada na variável arquivo e o tipo de formato a qual seria exibido o resultado: “Text”. O resultado da consulta foi armazenado na variável resultado e foi impresso através do comando “print resultado.read( )”. Observe que inserimos uma mensagem antes da execução do comando qblast( ), indicando o início da busca, e outra após a execução, avisando sobre o término. A impressão de mensagens na tela é importante para detectarmos em qual ponto do código há maior tempo de execução. Por se tratar de um comando que faz uma consulta remota à base de dados do NCBI, não podemos definir com precisão qual o tempo total de execução. Em nossos testes, obtivemos consultas realizadas entre 9 s e 267 s, entretanto esse número irá variar de acordo com a conexão com a internet do usuário.
O formato Text nos dá uma visualização simples dos resultados, entretanto podemos obter uma visualização mais eficiente utilizando o formato XML. Entretanto, a análise de arquivos no formato XML é um pouco mais complexa.
Analisando resultados de BLAST no formato XML
Vamos alterar o script anterior para obter como resultado um arquivo no formato XML. Por fim, vamos ler esse arquivo XML e imprimi-lo formatado na tela.
from Bio.Blast import NCBIWWW
from Bio.Blast import NCBIXML
from Bio import SeqIO
arquivo = SeqIO.read("exemplo.fasta", format="fasta")
print("Iniciando busca no NCBIWWW..." )
resultado = NCBIWWW.qblast("blastn", "nt", \
arquivo.seq, format_type="XML")
print("Busca concluida. Salvando resultados..." )
saida = open("blast_resultado.xml", "w")
saida.write(resultado.read())
saida.close()
# Lendo arquivos XML
arquivo_xml = open("blast_resultado.xml","r")
dados = NCBIXML.parse(arquivo_xml)
item = next(dados)
i = 1
for a in item.alignments:
for hsp in a.hsps:
print('Alinhamento',i)
print('Sequencia: '+a.title)
print('Tamanho: ',a.length)
print('Score: ',hsp.score)
print('Gaps: ',hsp.gaps)
print(hsp.query )
print(hsp.match)
print(hsp.sbjct)
print("\n")
i+=1
print("Executado com sucesso." )
Para a análise de arquivos em XML foi necessário importar o módulo NCBIXML. Esse módulo possui dois métodos: read para leitura de resultados de BLAST em arquivos XML e parse para análise iterativa de resultados. Em nosso exemplo anterior, após o fechamento do arquivo XML utilizamos a função open( ) para reabrí-lo novamente. Você deve estar pensando que essa não é a melhor maneira de realizar esta tarefa, e tem razão, afinal por que fecharíamos um arquivo e na próxima linha o abriríamos novamente? O programa perde em desempenho, entretanto essa é a forma mais didática de explicar como analisar arquivos XML.
Após a abertura do arquivo XML os comandos “dados = NCBIXML.parse(arquivo_xml)” e “item = next(dados)” permitem a ler o arquivo XML iterativamente e receber o próximo item em uma iteração, respectivamente. É declarada uma variável i que recebe o valor 1. Essa variável será utilizada como contador na impressão dos alinhamentos. Em seguida, dois laços de repetição percorrem os itens da lista de alinhamentos e para cada alinhamento é analisado o conteúdo de hsp (high-scoring segment pairs ou pares de segmentos de alta pontuação), onde estão contidos os resultados de pontuações de alinhamento (score), falhas no alinhamento (gaps), sequências buscadas (query), sequências correspondentes nos bancos de dados (sbjct) e as partes semelhantes entre sequências (match).
Se tudo ocorrer sem erros, o usuário recebe na tela mensagem: “Executado com sucesso.”.
Realizando BLAST local
Executar o BLAST localmente permite que suas buscas sejam realizadas com menor tempo. Ainda permite que você construa bases de dados personalizadas. Para executar o BLAST local, Biopython requer que a suíte de aplicativos BLAST+ esteja instalada. A última versão da suíte BLAST+ pode ser obtida em: <ftp://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/> e tem suporte para sistemas operacionais Windows, Linux e MacOS.
Biopython trata o BLAST como um Wrapper (uma ferramenta que pode ser acessada através de outra ferramenta), e fornece métodos de acesso através do módulo Applications do pacote Blast. Para isso, Biopython permite a construção de linhas de comando que executarão o BLAST, e por fim, essas serão executadas através do Python.
Como exemplo, vamos comparar duas sequências de nucleotídeos. Faça uma cópia do arquivo exemplo.fasta, criado anteriormente (ou obtê-lo no diretório exemplos do GitHub). Pode chamá-lo de exemplo2.fasta. Faça algumas modificações no arquivo, como por exemplo, insira um número dois no cabeçalho e altere o quinto nucleotídeo de “T” para “A”. Agora vamos compará-los com BLAST. Veja:
from Bio.Blast.Applications import *
comando_blastn = NcbiblastnCommandline( \
query="exemplo.fasta", subject="exemplo2.fasta", \
outfmt=0, out="out.txt")
print(comando_blastn)
# Executando
stdout, stderr = comando_blastn()
# Abrindo resultado
blast_result = open("out.txt","r")
linhas = blast_result.read()
print(linhas)
É possível determinar se o BLAST foi executado com sucesso se você obtiver como resultado algo similar a:
/usr/local/ncbi/blast/bin/blastn -out out.txt -outfmt 0 -query exemplo.fasta -subject exemplo2.fasta
BLASTN 2.2.27+
Query= beta-glucosidase_parcial_mRNA
...
Entretanto, se você obtiver como resultado: “blastn: command not found”, provavelmente o comando blastn não está declarado na variável $PATH. Você pode corrigir esse problema inserindo no seu script uma variável com localização do executável blastn e declarando essa variável no método NcbiblastnCommandline( ).
Em geral, BLAST+ é instalado nos diretórios:
- “/usr/bin/blastn” no Linux;
- “/usr/local/ncbi/blast/bin/blastn” no MacOS;
- “C:\Program Files\NCBI\blast-2.2.31+\bin\blastn” no Windows.
Entretanto, o endereço dos diretórios pode variar de acordo com a versão do BLAST ou do Sistema Operacional. Você pode tentar encontrá-lo via terminal utilizando os seguintes comandos:
- “whereis blastn” ou “locate blastn” no Linux;
- ‘mdfind “kMDItemFSName=blastn”’ no MacOS;
- no Windows abra um diretório qualquer, pressione F3 e busque por “blastn”.
Por exemplo, no MacOS o script correto ficaria desta forma:
from Bio.Blast.Applications import *
blastn = "/usr/local/ncbi/blast/bin/blastn"
comando_blastn = NcbiblastnCommandline(cmd=blastn, \
query="exemplo.fasta", subject="exemplo2.fasta", \
outfmt=0, out="out.txt")
print(comando_blastn)
# Executando
stdout, stderr = comando_blastn()
# Abrindo resultado
blast_result = open("out.txt","r")
linhas = blast_result.read()
print(linhas)
Caso deseje executar esse script em um sistema operacional diferente, altere o valor da variável blastn.
Nesse script gostaríamos de destacar a maneira com a qual a linha de comando foi executada. A classe NcbiblastnCommandline( ) inicialmente gera a linha de comando. Assim, a execução só é realizada com a chamada do objeto criado, e quando isso ocorre duas variáveis recebem os resultados dessa chamada: (i) stdout, que recebe a saída da execução, e (ii) stderr, que recebe o log de erros.
Parâmetros do BLAST
Parâmetros do BLAST podem ser aplicados a classe NcbiblastnCommandline( ). Veja alguns parâmetros abaixo:
db | Nome do banco de dados BLAST. Input: string; valor padrão: none (nenhum). |
query | Nome de um arquivo para busca. Input: string; valor padrão: requer o nome de um arquivo. |
out | Nome do arquivo onde os resultados de BLAST serão salvos. Input: string; valor padrão: stdout. |
evalue | E-value. Input: real; valor padrão: 10.0. |
subject | Arquivo com sequências onde será realizada a busca. Input: string; valor padrão: none. |
num_threads | Número de threads (CPUs) usadas em buscas pelo BLAST. Input: inteiro; valor padrão: 1. |
word_size
| Tamanho mínimo do fragmento de sequência usado na busca. Input: valor inteiro; valores padrões: (i) blastn: 11; (ii) blastp: 3; (iii) blastx: 3; (iv) tblastn: 3; (v) tblastx: 3; (vi) rpsblast: 40; (vii) megablast: 28; |
matrix | Define a matriz de pontuação usada no alinhamento (usado em alinhamento de proteínas). Input: string; valor padrão: “BLOSSUM 62”. |
outfmt | Opções de visualização de alinhamentos. Pode variar de 0 a 14. Input: string; valor padrão: “0”.
0 = par a par (textual), 1 = query-anchored mostrando identidades, 2 = query-anchored sem identidades, 3 = flat query-anchored (mostra identidade), 4 = flat query-anchored (sem identidades), 5 = XML, 6 = tabular, 7 = tabular com linhas comentadas, 8 = Textual ASN.1, 9 = Binário ASN.1, 10 = Valores separados por vírgula, 11 = Formato de arquivo BLAST (ASN.1), 12 = Saída JSON Seqalign, 13 = Saída JSON Blast, 14 = Saída XML2 Blast |
Resultados de BLAST em formato tabular
Até agora, analisamos resultados de BLAST em formato textual par a par e em formato XML, entretanto o formato tabular (6) pode auxiliar na visualização em paralelo de diversos resultados, podendo ser aberto através de editores de planilhas, como MS Office Excel ou o Calc, ou facilmente analisado através de scripts. Por padrão o formato tabular exibe os seguintes resultados: ‘qseqid sseqid pident length mismatch gapopen qstart qend sstart send evalue bitscore‘, entretanto outras opções podem ser adicionadas. Vamos entender o que cada um desses elementos representa:
qseqid | Nome da sequência usada na busca (query). Pode ser utilizada mais de uma sequência (arquivo Multi-FASTA). |
qlen | Tamanho da sequência usada na busca (query). |
sseqid | Nome da sequência a qual é realizada a busca (subject). Pode ser utilizada mais de uma sequência (arquivo Multi-FASTA). |
slen | Tamanho da sequência a qual é realizada a busca (subject). |
qstart | Posição de início do alinhamento na query. |
qend | Posição de término do alinhamento na query. |
sstart | Posição de início do alinhamento no subject. |
send | Posição de término do alinhamento no subject. |
evalue | Valor mínimo esperado para o evalue. |
score | Valor de pontuação. |
length | Indica tamanho do alinhamento. |
pident | Indica o percentual de correspondências idênticas (nucleotídeos/aminoácidos). |
nident | Indica o número de correspondências idênticas (nucleotídeos/aminoácidos). |
mismatch | Indica o número de mismatches (nucleotídeos/aminoácidos diferentes em um alinhamento). |
positive | Número de nucleotídeos/aminoácidos idênticos. |
gapopen | Número de gaps abertos. |
gaps | Número de gaps existentes. |
ppos | Percentual de nucleotídeos/aminoácidos classificados como positivos. |
bitscore | Pontuação bit. |
No exemplo abaixo, vamos gravar em um arquivo os atributos: qlen, qstart, qend, slen, sstart, send, positive, gaps, mismatch e pident.
from Bio.Blast.Applications import *
comando_blastn = NcbiblastnCommandline(query= \
"example.fasta", subject="example2.fasta", \
outfmt="'6 qlen qstart qend slen sstart send \
positive gaps mismatch pident'", out="out.txt")
print(comando_blastn)
# Executando
stdout, stderr = comando_blastn()
# Abrindo resultado
blast_result = open("out.txt","r")
linhas = blast_result.read()
print(linhas)
Observe a forma a qual o parâmetro outfmt é inserido na classe NcbiblastnCommandline( ). O valor contido em outfmt é registrado com aspas duplas seguido de aspas simples. Isso ocorre pois Biopython requer aspas duplas no parâmetro, e registra o conteúdo da variável como do tipo inteiro. Ao inserir as aspas simples após as aspas duplas, forçamos Biopython a enviar o argumento como string. Isso só é possível devido ao fato de estamos utilizando um wrapper e do BLAST aceitar como argumentos valores inteiros e strings no campo outfmt.
Se você executou o script com os arquivos de exemplo criados anteriormente, você deve ter obtido um resultado similar a este:
blastn -out out.txt -outfmt '6 qlen qstart qend slen sstart send positive gaps mismatch pident' -query exemplo.fasta -subject exemplo2.fasta
70 1 70 70 1 70 69 0 1 98.57
Na primeira linha temos a impressão do comando criado. A seguir os resultados na respectiva ordem: tamanho da sequência query (70), posição primeiro nucleotídeo do alinhamento na query (1), posição do último nucleotídeo do alinhamento na query (70), tamanho da sequência subject (70), posição primeiro nucleotídeo do alinhamento no subject (1), posição do último nucleotídeo do alinhamento no subject (70), número de nucleotídeos idênticos (69), número de gaps (0), número de mismatches (1) e percentual de identidade (98,57%).
Outras classes do Bio.Blast.Applications
Nos exemplos demonstrados anteriormente, apresentamos alinhamentos de nucleotídeos com a ferramenta blastn, entretanto é possível realizar alinhamentos utilizando as outras ferramentas da suíte BLAST+. Para isso é necessário apenas algumas alterações nos métodos utilizados. Segue abaixo, algumas classes do módulo Applications do pacote Bio.Blast:
Ncbiblastp Commandline( ) | Permite alinhamentos de sequências de proteínas com outras sequências de proteínas (blastp). |
Ncbiblastx Commandline( ) | Permite a comparação de uma sequência de nucleotídeos traduzidas, em todos os quadros de leitura, contra uma base de dados de sequências de proteínas (blastx). |
Ncbitblastn Commandline( ) | Permite a comparação de uma sequência de proteína contra uma base de dados de sequência de nucleotídeos traduzida dinamicamente em todos os quadros de leitura (tblastn). |
Ncbitblastx Commandline( ) | Permite a comparação de seis quadros de leitura traduzidos de uma sequência de nucleotídeos contra os seis quadros de leitura traduzidos de um banco de dados sequência de nucleotídeos (tblastx). |
Ncbipsiblast Commandline( ) | Permite o uso do programa psiblast. Exemplo: NcbipsiblastCommandline(cmd=’psiblast’, help=True) |
Ncbirpsblast Commandline( ) | Permite o uso do programa rpsblast. Exemplo: NcbirpsblastCommandline(cmd=rpsblast, help=True) |
Ncbirpstblastn Commandline( ) | Permite o uso do programa rpstblastn. Exemplo: NcbirpstblastnCommandline(cmd=rpstblastn, help=True) |
Ncbideltablast Commandline( ) | Permite o uso do programa deltablast (para proteínas). Exemplo: NcbideltablastCommandline(cmd=’deltablast’, query=’rosemary.pro’, db=’nr’, evalue=0.001, remote=True) |
Ncbiblastforma tterCommand line( ) | Permite a conversão de arquivos ASN.1 para outros formatos de saída BLAST. Exemplo: NcbiblastformatterCommandline(cmd=’blast_formatter’, out=’exemplo.xml’, outfmt=5, archive=’exemplo.asn’) |
Para mais informações sobre o pacote Bio.Blast, acesse: <http://biopython.org/DIST/docs/tutorial/Tutorial.html#htoc84>.
Construindo seu próprio banco de dados (opcional)
Em alguns casos, comparações entre uma sequência (query) e um grande quantidade de sequências (subject), como por exemplo um arquivo Multi-FASTA com centenas de milhares de sequências, pode demandar uma grande quantidade de tempo. Assim, torna-se viável a criação de bancos de dados que facilitem o trabalho de busca.
Infelizmente, Biopython não fornece métodos para criação de bancos de dados, entretanto é possível utilizar o comando makeblastdb da suíte BLAST+ para realizar essa tarefa.
Como exemplo, vamos converter nosso arquivo exemplo.fasta em um banco de dados (Linux e MacOS X; para windows acrescente o endereço completo do programa makeblastdb).
makeblastdb -in exemplo.fasta -dbtype 'nucl' -out exemplo_db -title exemplo_db
Você perceberá que três novos arquivos foram criados: exemplo_db.nhr, exemplo_db.nin e exemplo_db.nsq.
Agora vamos modificar o script usado anteriormente, substituindo o parâmetro subject por db e acrescentando o nome do banco de dados criado anteriormente: exemplo_db. Não se esqueça de alterar o nome da query, afinal criamos nossa base de dados com o arquivo exemplo.fasta, logo faz mais sentido compararmos o arquivo exemplo2.fasta com o banco de dados criado.
from Bio.Blast.Applications import *
comando_blastn = NcbiblastnCommandline( \
query="example2.fasta", db="exemplo_db", \
outfmt="'6 qlen qstart qend slen sstart send \
positive gaps mismatch pident'", out="out.txt")
print(comando_blastn)
# Executando
stdout, stderr = comando_blastn()
# Abrindo resultado
blast_result = open("out.txt","r")
linhas = blast_result.read()
print(linhas)
Você notou que o script pode ser executado mais rápido? É claro que não! Diferenças só serão perceptíveis com grandes bancos de dados.
Quer aprender mais? Conheça nossos cursos profissionalizantes à partir de R$19,99:
Capítulo 1
Introdução ao Python
Capítulo 2
Comandos condicionais
Capítulo 3
Laços de repetição
Capítulo 4
Trabalhando com strings
Capítulo 5
Listas
Capítulo 6
Manipulando arquivos
Capítulo 7
Funções
Capítulo 8
Princípios da orientação a objetos
Capítulo 9
Introdução ao Biopython
Capítulo 10
Sequências
Capítulo 11
BLAST
Capítulo 12
PDB
Capítulo 13
Visualização de dados em Python
Capítulo 14
Outras coisas interessantes que se pode fazer com Biopython
Capítulo 15
Hierarquia do Biopython