Find homologous sequences using iterative BLAST (PSI-BLAST), profile HMMs (HMMER), and reciprocal best hit analysis. Use when identifying orthologs, distant homologs, or protein family members where standard BLAST is not sensitive enough.
Find homologous sequences using iterative BLAST (PSI-BLAST), profile HMMs (HMMER), and reciprocal best hit analysis. Use when identifying orthologs, distant homologs, or protein family members where standard BLAST is not sensitive enough.
tool_type
mixed
primary_tool
BLAST+
Sequence Similarity Searches
Advanced methods for finding homologous sequences beyond standard BLAST.
PSI-BLAST (Position-Specific Iterated BLAST)
Builds a position-specific scoring matrix (PSSM) through iterations to find distant homologs.
Basic PSI-BLAST
psiblast -query protein.fasta -db nr -out results.txt -num_iterations 3
awk 'FNR==NR {a[$1]=$2; next} $2 in a && a[$2]==$1 {print $1"\t"$2}' \
A_vs_B.txt B_vs_A.txt > reciprocal_best_hits.txt
Python RBH Script
deffind_rbh(forward_blast, reverse_blast):
'''Find reciprocal best hits from BLAST results'''
forward = {}
withopen(forward_blast) as f:
for line in f:
parts = line.strip().split('\t')
query, subject = parts[0], parts[1]
if query notin forward:
forward[query] = subject
reverse = {}
withopen(reverse_blast) as f:
for line in f:
parts = line.strip().split('\t')
query, subject = parts[0], parts[1]
if query notin reverse:
reverse[query] = subject
rbh = []
for a, b in forward.items():
if b in reverse and reverse[b] == a:
rbh.append((a, b))
return rbh
rbh_pairs = find_rbh('A_vs_B.txt', 'B_vs_A.txt')
for a, b in rbh_pairs:
print(f'{a}\t{b}')
Delta-BLAST
Uses conserved domain database for more sensitive initial search.
deltablast -query protein.fasta -db nr -rpsdb cdd_delta -out results.txt
PHI-BLAST (Pattern-Hit Initiated)
Search with a pattern plus sequence.
phi_pattern="G-x(2)-[ST]-x-[RK]"
phiblast -query protein.fasta -db nr -pattern "$phi_pattern" -out results.txt
Iterative Search with Biopython
from Bio.Blast import NCBIWWW, NCBIXML
withopen('query.fasta') as f:
query = f.read()
result_handle = NCBIWWW.qblast('psiblast', 'nr', query, expect=0.001, word_size=3)
withopen('psiblast_result.xml', 'w') as out:
out.write(result_handle.read())
result_handle.close()
withopen('psiblast_result.xml') as f:
records = NCBIXML.parse(f)
for record in records:
for alignment in record.alignments:
for hsp in alignment.hsps:
if hsp.expect < 1e-10:
print(f'{alignment.hit_def[:50]}: E={hsp.expect}')
HMMER with Biopython
from Bio import SearchIO
results = SearchIO.parse('hmmsearch_output.txt', 'hmmer3-text')
for query_result in results:
print(f'Query: {query_result.id}')
for hit in query_result:
print(f' Hit: {hit.id}, E-value: {hit.evalue}')
for hsp in hit:
print(f' Domain: {hsp.bitscore} bits')
