Categorias
Introdução à Programação para Bioinformática com Perl

BLAST

Capítulo 11

Este conteúdo faz parte do livro “Introdução à Programação para Bioinformática com Perl“. 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 (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.
blastpblastp 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.
tblastnCompara 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.
tblastxCompara 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.

BioPerl trata o BLAST como um Wrapper (uma ferramenta que pode ser acessada através de outra ferramenta) e apresenta o módulo Bio::SearchIO que permite que os resultados de BLAST possam ser analisados.

Instalando a Suíte BLAST+

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, BioPerl 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.

Para Windows faça o download do arquivo de extensão “exe”. Para MacOS faça o download do arquivo de extensão “dmg”. Para Linux você pode fazer o download do arquivo compactado extensão “tar.gz”. Ubuntu permite ainda a instalação via linha de comando:

sudo apt-get install ncbi-blast+

Executando BLAST via Perl

BLAST pode ser executado via linha de comando. Para isso deve-se determinar qual aplicativo deverá ser executado, blastn, blastp, blastx, tblastn ou tblastx, seguido do nome de um arquivo com a sequência que será pesquisada e de um arquivo com as sequências alvo ou do nome base de dados.

blastn -query sequencia_buscada.fasta -subject sequencias_alvo.fasta

BLAST ainda permite que parâmetros sejam utilizados para incrementar buscas. Veja alguns parâmetros abaixo:

dbNome do banco de dados BLAST. Input: string; valor padrão: none (nenhum).
queryNome de um arquivo para busca. Input: string; valor padrão: requer o nome de um arquivo.
outNome do arquivo onde os resultados de BLAST serão salvos. Input: string; valor padrão: stdout.
evalueE-value. Input: real; valor padrão: 10.0.
subjectArquivo com sequências onde será realizada a busca. Input: string; valor padrão: none.
num_threadsNúmero de threads (CPUs) usadas em buscas pelo BLAST. Input: inteiro; valor padrão: 1.
word_sizeTamanho 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;
matrixDefine a matriz de pontuação usada no alinhamento (usado em alinhamento de proteínas). Input: string; valor padrão: “BLOSSUM 62”.
outfmtOpções de visualização de alinhamentos. Pode variar de 0 a 13. 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.
Fonte: https://www.ncbi.nlm.nih.gov/books/NBK279684/

BioPerl oferece o módulo Bio::Tools::Run::StandAloneBlast para execução de buscas com BLAST. Entretanto, esse módulo se tornou depreciado, uma vez que NCBI passou a não oferecer suporte ao programa blastall, que foi substituído pela suíte BLAST+. Assim, BioPerl passou a oferecer o módulo Bio::Tools::Run::StandAloneBlastPlus.

Devido à diversidade de módulos do BioPerl e de versões do BLAST, além do fato de que BLAST pode ser executado por linha de comando ou pelo Website <http://blast.ncbi.nlm.nih.gov/Blast.cgi> decidimos ensinar um método global para executar BLAST via Perl e focar na análise de resultados. Observe abaixo como realizar uma comparação entre dois arquivos. Neste exemplo faremos a comparação do “arquivo4.fasta” (criado anteriormente) com ele mesmo (para que possamos ter certeza que haverá pelo menos uma sequência igual).

my $blast = "blastn";
my $query = "arquivo4.fasta";
my $subject   = "arquivo4.fasta";
my $opcoes = "";

my $comando = "$blast -query $query -subject $subject $op-coes";

open (my $fh,"$comando |") or die("Não foi possível execu-tar o comando: $comando: $!\n");

while(my $l = <$fh>){
	print $l;
}

Observe que primeiro definimos qual aplicativo do BLAST iriamos utilizar (blastn). A seguir, definimos a sequência que seria buscada (query) e a sequência alvo (subject). Lembre-se que em uma comparação real, as sequências estariam em arquivos diferentes, mas para nosso teste queremos ter certeza que haja ao menos um resultado (hit). Definimos uma variável para receber parâmetros opcionais. Em seguida, criamos uma variável para armazenar o comando que será executado na última linha. O resultado será gravado em um arquivo cujo apelido está armazenado na variável $fh. Veja abaixo um trecho do resultado:

BLASTN 2.2.31+
Query= NC_009934 Acaryochloris marina MBIC11017 plasmid pREB9, complete sequence.
Length=2133
Subject= NC_009934 Acaryochloris marina MBIC11017 plasmid pREB9, complete sequence.
Length=2133
 Score = 3940 bits (2133),  Expect = 0.0
 Identities = 2133/2133 (100%), Gaps = 0/2133 (0%)
 Strand=Plus/Plus
Query  1     GCTAGATCCTGCTAACAGACCAGTGCCTGATTGCTGGATGAGGTT  60
             |||||||||||||||||||||||||||||||||||||||||||||
Sbjct  1     GCTAGATCCTGCTAACAGACCAGTGCCTGATTGCTGGATGAGGTT  60
# [...]

Uma das principais vantagens do BioPerl são os módulos que permitem a análise de resultados de diversos programas, dentre eles BLAST. Perl fornece o módulo Bio::SearchIO para ler e gravar resultados de BLAST. Veja como ficaria o script alterado para podermos analisar os resultados.

use Bio::SearchIO;

my $blast = "blastn";
my $query = "arquivo4.fasta";
my $subject = "arquivo4.fasta";
my $opcoes = "-outfmt 5";

my $comando = "$blast -query $query -subject $subject $op-coes";

open (my $fh,"$comando |") or 
die("Não foi possível executar o comando: $comando: $!\n");

my $in = Bio::SearchIO->new(
-format => 'blastxml', 
-fh => $fh
);

while(my $r = $in->next_result){

    print "Numero de resultados: ".$r->num_hits."\n";
    print "Nome da query: ".$r->query_name."\n";
    print "Descricao: ".$r->query_description."\n";
    print "Tamanho da query: ".$r->query_length."\n";

    # Para cada resultado
    while(my $h = $r->next_hit ) {
      print "Hit: ".$h->name."\n";

      while(my $hsp = $h->next_hsp ) {

    	     print "Tamanho HSP: ". $hsp->length('total')."\n";
            print "E-value: ", $hsp->evalue."\n";
            print "Bit score: ", $hsp->score."\n";
            print "Regiao de hit (query): ".
            $hsp->query->start." - ".$hsp->query->end."\n";
            print "Regiao de hit (subject): ".
            $hsp->hit->start." - ".$hsp->hit->end."\n";

    	}

    }

    print "------------------------------------------\n";

}

A primeira alteração no script foi a declaração do módulo Bio::SearchIO. Outra modificação foi a inserção do parâmetro “-outfmt 5” na variável $opcoes. Esse parâmetro permite que o BLAST retorne seus resultados em formato XML, o que evita possíveis problemas com incompatibilidade de versões do BLAST.

O método new do módulo Bio::SearchIO recebe como argumentos o formato do arquivo enviado “blastxml” e o nome da variável que armazena o apelido do arquivo. A partir desse ponto, três laços são executados para imprimir os campos desejados:

• Laço 1 (Results): através do método next_results, esse la-ço percorre todos os resultados;
• Laço 2 (Hits): através do método next_hit, esse laço per-corre todos os hits (resultados com acertos). Cada resulta-do pode apresentar vários hits;
• Laço 3 (HSP): através do método next_hsp, esse laço per-corre todos os HSP (High-scoring Segment Pairs ou pares de segmentos de alta pontuação);

No exemplo demonstrado, imprimimos algumas informações como: nome e descrição da query, tamanho da query, total de hits e posições de regiões similares.

Result, Hit e HSP são as três maneiras que Bio::SearchIO fornece para percorrer resultados de BLAST. Além dos métodos apresentados no script anterior, cada um deles possui outros métodos que poderão ser úteis. Segue abaixo a lista de métodos.

Tabela de métodos do objeto result

MétodoRetorna
algorithm Aplicativo utilizado (ex.: blastn).
algorithm_versionVersão do aplicativo.
query_name Nome dado a query. Geralmente consta no cabeçalho do arquivo.
query_accessionCódigo de acesso da query (quando retirada de bancos de dados públicos, como GenBank).
query_lengthTamanho da query.
query_descriptionDescrição da query.
database_name Nome do banco de dados utilizado (quando utilizada o parâmetro “-db”).
database_lettersTotal de caracteres na base de dados.
database_entriesNúmero de entradas de bancos de dados.
available_statisticsEstatísticas utilizadas.
available_parametersParâmetros utilizados.
num_hitsTotal de acertos.
hitsLista todos os hits.
rewindReseta o iterador (next_hit).
Fonte: <http://www.bioperl.org/wiki/HOWTO:SearchIO> (adaptado).

Tabela de métodos do objeto hit

MétodoRetorna
nameTítulo dado ao hit.
lengthTamanho do hit.
accessionCódigo de acesso da query (quando retirada de bancos de dados públicos, como NCBI).
descriptionDescrição do hit.
algorithmAplicativo utilizado (ex.: blastn)
raw_scorePontuação obtida pelo alinhamento.
significanceSignificância do hit.
bitsBits do hit.
hspsLista com hsps.
num_hspsNúmero de hsps.
locusNome do locus.
accession_numberNúmero de acesso.
rewindReseta o iterador (next_hsp)
Fonte: http://www.bioperl.org/wiki/HOWTO:SearchIO (adaptado).

Tabela de métodos do objeto hsp

MétodoRetorna
algorithmAplicativo utilizado (ex.: blastn)
evalueE-value obtido.
expectMesmo que e-value.
frac_identicalFração idêntica.
frac_conservedFração conservada.
gapsNúmero de gaps (regiões em que há uma quebra na identidade).
query_stringString que armazena região da query.
hit_stringString que armazena região do hit.
homology_stringString que armazena região homóloga.
length(‘total’)Tamanho de HSP (incluindo gaps).
length(‘hit’)Tamanho da região do hit (menos gaps).
length(‘query’)Tamanho da região do hit (menos gaps).
hsp_lengthTamanho de HSP (mesmo que total).
num_conservedNúmero de resíduos conservados.
num_identicalNúmero de resíduos identicos.
rankRank do HSP.
scorePontuação.
bits Pontuação em bits.
range(‘query’)Posição de início e fim (query), como array.
range(‘hit’)Posição de início e fim (hit), como array.
percent_identityPercentual de identidade.
start(‘query’)Posição de início do alinhamento na query.
end(‘query’)Posição de fim do alinhamento na query.
start(‘hit’)Posição de início do alinhamento no hit.
end(‘hit’)Posição de fim do alinhamento no hit.
matches(‘hit’)Número de idênticos e conservados, como array.
matches(‘query’)Número de idênticos e conservados, como array.
get_alnObjeto Bio::SimpleAlign.
Fonte: http://www.bioperl.org/wiki/HOWTO:SearchIO (adaptado)

Resultados de BLAST em formato tabular

Até agora, analisamos resultados de BLAST em formato textual (0) e XML (5), entretanto o formato tabular (6) pode auxiliar na visualização em paralelo de diversos resultados. O formato tabular pode ainda ser aberto através de editores de planilhas, como MS Office Excel ou Calc. Pode ainda ser facilmente analisado com expressões regulares em Perl. 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:

qseqidNome da sequência usada na busca (query). Pode ser utilizada mais de uma sequência (arquivo Multi-FASTA).
qlenTamanho da sequência usada na busca (query).
sseqidNome da sequência a qual é realizada a busca (subject). Pode ser utilizada mais de uma sequência (arquivo Multi-FASTA).
slenTamanho da sequência a qual é realizada a busca (subject).
qstartPosição de início do alinhamento na query.
qendPosição de término do alinhamento na query.
sstartPosição de início do alinhamento no subject.
sendPosição de término do alinhamento no subject.
evalueValor mínimo esperado para o e-value.
scoreValor de pontuação.
lengthIndica tamanho do alinhamento.
pident Indica o percentual de correspondências idênticas (nucleotídeos/aminoácidos).
nidentIndica o número de correspondências idênticas (nucleotídeos/aminoácidos).
mismatchIndica o número de mismatches (nucleotídeos/aminoácidos diferentes em um alinhamento).
positiveNúmero de nucleotídeos/aminoácidos idênticos.
gapopenNúmero de gaps abertos.
gaps Número de gaps existentes.
pposPercentual de nucleotídeos/aminoácidos classificados como positivos.
bitscorePontuação bit.
Fonte: https://www.ncbi.nlm.nih.gov/books/NBK279684/

No exemplo a seguir vamos substituir o formato para tabular e faremos a análise utilizando funções nativas do Perl.

my $blast = "blastn";
my $query = "arquivo4.fasta";
my $subject = "arquivo4.fasta";
my $opcoes = "-outfmt 6";

my $comando = "$blast -query $query -subject $subject $op-coes";

open (my $fh,"$comando |") or die("Não foi possível execu-tar o comando: $comando: $!\n");

while(my $linha = <$fh>){

  my @l = split("\t",$linha);
  my $query2 = $l[0];
  my $subject2 = $l[1];
  my $identidade = $l[2];
  print "Query: $query2\nSubject: $subject2\n
  Identidade: $identidade\n";

}
#Query: NC_009934
#Subject: NC_009934
#Identidade: 100.00

Observe que para cada linha aplicamos a função split. Essa função divide uma string na posição indicada pelo primeiro argumento e armazena o resultado em um array. No nosso exemplo, imprimimos as posições 0, 1 e 2 do array, que armazenam nome da query, nome do subject e identidade do alinhamento.

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 importante a criação de bancos de dados para agilizar a busca.

Infelizmente, BioPerl 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 “exemplo4.fasta”em um banco de dados (Linux e MacOS X; para Windows acrescente o endereço completo do programa makeblastdb).

makeblastdb -in exemplo4.fasta -dbtype 'nucl' -out exemplo4_db -title exemplo4_db

Você perceberá que três novos arquivos foram criados: exemplo4_db.nhr, exemplo4_db.nin e exemplo4_db.nsq.

Agora vamos modificar o script criado anteriormente, substituir o parâmetro subject por db e acrescentar o nome do banco de dados criado anteriormente: exemplo4_db.

Agora execute o script. Você notou que o script pode ser executado mais rápido? Acreditamos 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:

Nota do autor
Prefácio

Capítulo 1
Introdução ao Perl

Capítulo 2
Comandos condicionais

Capítulo 3
Strings

Capítulo 4
Arrays

Capítulo 5
Laços de repetição

Capítulo 6
Manipulando arquivos

Capítulo 7
Sub-rotinas

Capítulo 8
“O guia de sobrevivência para expressões regulares em Perl”

Capítulo 9
Introdução ao BioPerl

Capítulo 10
Sequências

Capítulo 11
BLAST

Capítulo 12
Estruturas de proteínas

Capítulo 13
Hierarquia do BioPerl

Epílogo
Referências bibliográficas
Sobre os autores

Por favor, nos cite:

MARIANO, DIEGO CÉSAR BATISTA; de MELO-MINARDI, R. C. . Introdução à Programação para Bioinformática com Perl. 1. ed. North Charleston, SC (EUA): CreateSpace Independent Publishing Platform, 2016. v. 2. 200p .

Por Diego Mariano

Doutor em Bioinformática pela Universidade Federal de Minas Gerais com atuação na área de ciência de dados e aprendizado de máquina aplicados ao aperfeiçoamento de enzimas usadas na produção de biocombustíveis. Mestre em Bioinformática, também pela UFMG, atuando na área de desenvolvimento de sistemas Web para montagem de genomas. Atualmente realiza estágio pós-doutoral no Departamento de Ciência da Computação da UFMG com foco em desenvolvimento de sistemas Web para Bioinformática, análise exploratória e visualização de dados. Tem conhecimentos nas linguagens: PHP, JavaScript, Python, R, Perl, HTML, CSS e SQL.