Estoy tratando de encontrar homólogos de un conjunto de proteínas usando BLASTp. Estoy trabajando con bases de datos personalizadas.
Estoy usando un valor de 0.00001 como umbral.
Me gustaría filtrar consultas que tengan resultados con más del 90 % de identidades. Debido a que la salida de BLASTp se basa en HSP, no puedo filtrar por %identities/query, solo por HSP.
Me gustaría saber cómo hacer esto y también si estoy siguiendo una estrategia razonable.
Aquí hay un ejemplo de alineación: qcovs=100 pero qcovhsp inferior.
qseqid sseqid pident length mismatch gapopen qstart qend sstart send evalue bitscore qcovs qcovhsp
HPNK_01698 HAPS_0519 81.88 596 75 5 630 1225 615 1177 0.0 889 100 49
HPNK_01698 HAPS_0519 49.17 301 115 8 84 366 201 481 2e-56 214 100 23
HPNK_01698 HAPS_0519 53.64 261 61 6 436 684 616 828 6e-49 191 100 20
HPNK_01698 HAPS_0519 46.61 251 62 3 332 510 584 834 6e-46 181 100 15
HPNK_01698 HAPS_0519 53.27 214 79 4 1 194 1 213 1e-45 180 100 16
HPNK_01698 HAPS_0519 55.96 218 60 8 550 764 643 827 1e-40 164 100 18
HPNK_01698 HAPS_0519 51.56 225 61 7 516 731 642 827 1e-38 157 100 18
HPNK_01698 HAPS_0519 49.57 230 77 6 484 713 643 833 1e-37 154 100 19
HPNK_01698 HAPS_0519 57.89 76 26 1 364 433 760 835 1e-13 76.3 100 6
Código utilizado
hacer base de datos
makeblastdb -in $Hparasuisfastadatabase -out H_parasuis_strains_gb_ALL.fna_databaseBLAST -dbtype prot -parse_seqids
Ejecutar BLAST
blastp -db H_parasuis_strains_gb_ALL.fna_databaseBLAST -query 'out_2.fasta' -out HPNK_selected_vs_H_parasuis_strainss.tblastn -evalue 0.00001 -outfmt "6 qseqid sseqid pident length mismatch gapopen qstart qend sstart send evalue bitscore qcovs qcovhsp" -max_target_seqs 50
gracias, bernardo
En primer lugar, si quieres 90%identidad, puedes descartar este golpe. Ninguno de los HSP supera ese umbral. Además, dado que está trabajando con proteínas, no hay problemas de empalme involucrados y debería poder obtener un único HSP que abarque la mayoría de las secuencias de consulta y asunto. Asumiendo, por supuesto, que tienes un verdadero homólogo.
En su salida, veo muchos HSP pequeños superpuestos, la mayoría de los cuales tienen una identidad baja. No puedo estar seguro sin ver la secuencia, pero es una apuesta segura que lo que tiene allí son regiones repetitivas/de baja complejidad y es por eso que tiene tantos HSP separados. El único medio decente comienza en la posición 630 de la secuencia de consulta y tiene solo 595 residuos, menos de la mitad de su proteína de consulta. O tiene una región N-terminal muy divergente o su HSP es solo un dominio conservado. Nuevamente, necesitaría ver la alineación de la secuencia real para estar seguro, pero eso no parece un verdadero homólogo (suponiendo que sus especies estén razonablemente cerca, lo que debe ser si está usando un umbral de identidad del 90%).
Por lo tanto, siempre suponiendo que su especie esté lo suficientemente cerca como para esperar homólogos decentes, simplemente ignoraría los HSP más cortos y trataría con aquellos que representan más del 80 % de la longitud de mi consulta con una identidad >=90 %. Los hits más cortos serán dominios conservados o regiones repetitivas/de baja complejidad. Los umbrales que elija dependen de la especie que esté estudiando.
Si su especie no es tan cercana, no use BLASTP en absoluto. En su lugar, puede usar algo como hmmer. Recopile un conjunto de homólogos de varias especies para cada una de sus proteínas de consulta, construya una matriz usándolas y use esa matriz para buscar en su base de datos. También podría usar Selenoprofiles que usa un enfoque similar.
biotecnología
biotecnología
terdón
WYSIWYG
biotecnología