Ich versuche mit BLASTp Homologe einer Reihe von Proteinen zu finden. Ich arbeite mit benutzerdefinierten Datenbanken.
Ich verwende einen Wert von 0,00001 als Schwellenwert.
Ich möchte Abfragen filtern, die Treffer mit >90 % Identitäten haben. Da die BLASTp-Ausgabe auf HSPs basiert, kann ich nicht nach %identities/query filtern, sondern nur nach HSP.
Ich würde gerne wissen, wie das geht und ob ich eine vernünftige Strategie verfolge.
Hier ist ein Ausrichtungsbeispiel: qcovs=100 aber qcovhsp niedriger.
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
Code verwendet
Datenbank machen
makeblastdb -in $Hparasuisfastadatabase -out H_parasuis_strains_gb_ALL.fna_databaseBLAST -dbtype prot -parse_seqids
Führen Sie BLAST aus
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
Danke, Bernhard
Zuallererst, wenn Sie 90%
Identität wollen, können Sie diesen Treffer verwerfen. Keiner der HSP überschreitet diese Schwelle. Da Sie mit Proteinen arbeiten, gibt es darüber hinaus keine Spleißprobleme, und Sie sollten in der Lage sein, ein einziges HSP zu erhalten, das die meisten Abfrage- und Subjektsequenzen umfasst. Vorausgesetzt natürlich, Sie haben ein echtes Homolog.
In Ihrer Ausgabe sehe ich viele kleine, sich überschneidende HSPs, von denen die meisten eine geringe Identität haben. Ich kann nicht sicher sein, ohne die Sequenz zu sehen, aber es ist eine sichere Wette, dass das, was Sie dort haben, Bereiche mit geringer Komplexität/sich wiederholende Bereiche sind, und das ist der Grund, warum Sie so viele separate HSP haben. Das nur halbwegs anständige beginnt an Position 630 der Abfragesequenz und ist nur 595 Reste lang, weniger als die Hälfte Ihres Abfrageproteins. Entweder haben Sie eine sehr divergierende N-terminale Region oder Ihr HSP ist nur eine konservierte Domäne. Auch hier müsste ich die tatsächliche Sequenzausrichtung sehen, um sicher zu sein, aber das sieht nicht wie ein echtes Homolog aus (vorausgesetzt, Ihre Arten sind sich ziemlich ähnlich, was sie sein müssen, wenn Sie eine Identitätsschwelle von 90% verwenden).
Also, immer vorausgesetzt, dass Ihre Spezies nahe genug ist, um anständige Homologe zu erwarten, würde ich einfach die kürzeren HSPs ignorieren und mich mit denen befassen, die mehr als sagen wir 80 % der Länge meiner Abfrage bei> = 90 % Identität darstellen. Kürzere Treffer werden meistens konservierte Domänen oder Regionen mit sich wiederholender/niedriger Komplexität sein. Die von Ihnen gewählten Schwellenwerte hängen von der Art ab, die Sie untersuchen.
Wenn Ihre Arten nicht so nah sind, verwenden Sie BLASTP überhaupt nicht. Stattdessen können Sie etwas wie verwenden hmmer
. Sammeln Sie für jedes Ihrer Suchproteine eine Reihe von Homologen aus verschiedenen Arten, erstellen Sie daraus eine Matrix und verwenden Sie diese Matrix, um Ihre Datenbank zu durchsuchen. Sie könnten auch Selenoprofiles verwenden, das einen ähnlichen Ansatz verwendet.
Biotech
Biotech
Terdon
WYSIWYG
Biotech