Questão:
Por que um hit BLAST muito forte se perde quando eu altero o parâmetro num_alignments, num_descriptions ou max_target_seqs?
voiDnyx
2017-11-14 21:21:32 UTC
view on stackexchange narkive permalink

Isenção de responsabilidade: esta é uma pergunta auto-respondida para fins de documentação e adaptei-a do seguinte github gist. Especialmente dos usuários terrycojones e peterjc , bem como sujaikumar que levantaram o problema.


Tenho um estranho situação. Eu tenho dois bancos de dados BLAST ( 20170821 e 20170824 ) que contêm uma determinada sequência de assunto LC074724 . Quando eu os questiono com uma determinada consulta, Q , em um caso ( 20170824 ) há uma correspondência muito forte e no outro ( 20170821 ) nenhuma correspondência .

Por que não foi encontrada nenhuma correspondência de qualidade inferior embora os bancos de dados sejam quase os mesmos?

  BLAST versão $ blastn -versionblastn: 2.6.0+ Pacote: blast 2.6.0, build 11 de maio de 2017 22:22:40  

Ambos os bancos de dados contêm um assunto idêntico sequência LC074724:

  $ blastdbcmd -entry LC074724 -db 20170824 | md5fcfbc167c7f81dfd75aad1bf3db6220b $ blastdbcmd -entry LC074724 -db 20170821 | md5fcfbc167c7f81dfd75aad1bf3db6220b  

A única diferença para os bancos de dados é 20170821 é apenas 20170824 com 643 sequências extras. O resto é completamente idêntico.

Um responda:
voiDnyx
2017-11-14 21:36:38 UTC
view on stackexchange narkive permalink

Estou tentando descobrir o porquê. Esta será uma leitura mais longa, Tl; dr no final :

Fazendo uma correspondência

Há uma boa correspondência de Q contra LC074724 em 20170824:

  $ blastn -outfmt '6 qaccver saccver bitscore' -db 20170824 -query Q.fasta \ -task blastn -max_hsps 1 | grep LC074724Q LC074724.1 2581  

Observe que o bitcore da correspondência LC074724 é 2581.

Mas não há nenhuma correspondência de LC074724 em 20170821:

  $ blastn -outfmt '6 qaccver saccver bitscore' -db 20170821 -query Q.fasta \ -task blastn -max_hsps 1 | grep LC074724  

A opção -outfmt '6 qaccver saccver bitscore' diz ao BLAST para imprimir o número de acesso da consulta e versão, o acesso do assunto e versão, e o bitscore da partida.

Quais são os melhores resultados para 20170821

Então, quais são os bitscores das partidas Q em 20170821?

  $ blastn - outfmt '6 bitscore' -db 20170821 -query Q.fasta -task blastn -max_hsps 1 | sort -nr | head3301261625902590256325542547254425402540  

O bitcore de Q contra LC074724 em 20170824 (2581) o teria colocado como na 5ª posição entre as 10 melhores correspondências de Q , mas LC074724 não é correspondido de forma alguma em 20170821 ! O menor escore de bits correspondente em 20170821 é 2354.

Qual é a diferença entre os bancos de dados?

Obtenha os ids dos assuntos em cada banco de dados

  $ blastdbcmd -entry all -db 20170824 | egrep '^ >' | cut -c2- | sort > 20170824.ids $ blastdbcmd -entry all -db 20170821 | egrep '^ >' | cut -c2- | sort > 20170821.ids $ wc -l 2017082? .ids 9754 20170821.ids 9111 20170824.ids 18865 total  

20170821 tem 643 sequências que não estão em 20170824:

  $ comm -23 20170821.ids 20170824.ids | wc -l 643  

20170824 não tem sequências que não estejam em 20170821:

  $ comm -13 20170821.ids 20170824.ids | wc -l 0  

e os dois têm 9111 sequências em comum:

  $ comm -12 20170821.ids 20170824.ids | wc -l 9111  

Esses números correspondem aos números de wc -l acima.

Portanto, 20170821 é apenas 20170824 com 643 sequências extras.

Talvez LC074724 corresponda, mas com uma pontuação de bits baixa?

O BLAST mostra apenas os 500 principais assuntos correspondidos por padrão. Vamos pedir para mostrar mais:

  $ blastn -outfmt '6 qaccver saccver bitscore' -db 20170821 -query Q.fasta \ -task blastn -max_hsps 1 -max_target_seqs 10000 > 20170821-10000- Match.txt  

E ....

  $ head 20170821-10000-Match.txtQ Q 3301Q X234A 2616Q CS388973.1 2590Q DM059402.1 2590Q LC074724.1 2581Q KJ843188.1 2576Q AB697496.1 2576Q AB697499.1 2576Q AB697503.1 2576Q AB697508.1 2576  

Pronto! E está na 5ª posição, como deveria estar (com base nos bitscores discutidos acima).

Apenas para confirmar, retirando o -max_target_seqs opção não obtém correspondência:

  $ blastn -outfmt '6 qaccver saccver bitscore' -db 20170821 -query Q.fasta \ -task blastn -max_hsps 1 | grep LC074724 $  

Três opções BLAST com nomes confusos

Sempre achei os nomes e descrições das três opções de linha de comando do BLAST a seguir muito confusos.

Na seção Opções de formatação de blastn -help:

-num_descriptions = 0> Número de sequências de banco de dados para mostrar descrições de uma linha para Não aplicável para outfmt> 4
Padrão = `500 '* Incompatível com: max_target_seqs

-num_alignments = 0> Número de sequências de banco de dados para mostrar alinhamentos para Default = `250 '* Incompatível com: max_target_seqs

E da seção Restringir pesquisa ou resultados:

-max_target_seqs = 1> Número máximo de sequências alinhadas a manter Não aplicável para outfmt < = 4 Padrão = `500 '* Incompatível com: num_descriptions, num_alignments

De nomes de seção, parece que apenas a última opção ( -max_target_seqs ) teria algum efeito na pesquisa, e que as duas primeiras são apenas sobre o que é exibido .

Mas ... usar apenas -num_descriptions ou -num_alignments com um grande valor também obtém uma correspondência:

  $ blastn -db 20170821 -query Q.fasta -task blastn -max_hsps 1 -num_descriptions 10000 | grep LC074724LC074724.1 DNA do vírus da hepatite B, genoma completo, isolado: p621 2581 0.0>LC074724.1 DNA do vírus da hepatite B, genoma completo, isolado: p621 $ blastn -db 20170821 -query Q.fasta -task blastn -max_hs00al | grep LC074724LC074724.1 DNA do vírus da hepatite B, genoma completo, isolado: p621 2581 0.0>LC074724.1 DNA do vírus da hepatite B, genoma completo, isolado: p621  

Então, o BLAST está aparentemente encontrando o LC074724 corresponde mesmo sem a opção -max_target_seqs alta, mas não está exibindo a menos que haja um -num_alignments ou -num_descriptions code alto > valor. ( Tenha em mente que este é o quinto melhor acerto! Não em algum lugar> 500 )

Apenas para confirmar, quando nenhuma das 3 opções é fornecida e nenhuma opção -outfmt também é , a correspondência não foi encontrada:

  $ blastn -db 20170821 -query Q.fasta -task blastn -max_hsps 1 | grep LC074724 $  

Tl; dr:

Em resumo, as opções

  -evalue -max_target_seqs -num_alignments -num_descriptions

são NÃO apenas filtros aplicados aos resultados, mas influenciam diretamente os resultados encontrados .

Isso significa que você pode terminar com uma sequência bacteriana como BEST hit ao usar um valor x , e sem sequências bacterianas de todo em seu top20 ao usar o valor y.

Para citar o serviço de usuário NCBI:

Não consideramos isso um bug, mas concordo que devemos documentar melhor essa possibilidade. Isso pode acontecer porque os limites, incluindo as sequências máximas de destino, são aplicados em uma fase inicial sem intervalo do algoritmo, bem como posteriormente.

E

O cortes são aplicados durante uma fase de extensão sem lacuna e alguns alinhamentos podem melhorar na pontuação durante a extensão com lacuna. A abordagem mais segura é relaxar (aumentar) o limite e lidar com um conjunto de resultados maior.

Esse foi um comportamento totalmente inesperado para mim e, até onde sei, não está documentado em nenhum lugar do BLAST ajuda ou semelhante.



Estas perguntas e respostas foram traduzidas automaticamente do idioma inglês.O conteúdo original está disponível em stackexchange, que agradecemos pela licença cc by-sa 3.0 sob a qual é distribuído.
Loading...