基因家族分析之同源基因的寻找

Blast进行同源基因的寻找

参考博客：

• http://boyun.sh.cn/bio/?p=1849
• https://www.bioinfo-scrounger.com/archives/673/

基于蛋白的比对结果，寻找某一个蛋白家族的同源基因，使用如下的参数

• identity >30%；
• e-value <1e-10；
• score>200
• overlap >60%

首先对感兴趣的基因家族蛋白序列建立索引

makeblastdb -in test.fsa -parse_seqids -dbtype prot -out test_db

然后使用blastp比对到建立好索引的数据库

-outfmt "6 qseqid sseqid pident length mismatch gapopen qstart qend sstart send evalue bitscore qcovs"

blastp -query Bju.chr.modified_id.pep.fa -db MAGL_pep -outfmt "6 qseqid sseqid pident length mismatch gapopen qstart qend sstart send evalue bitscore qcovs" -max_hsps 1  -num_alignments 1 -evalue 1e-10 -num_threads 30  -out Bju.MAGL.v2.txt

如果需要blast比对返回一个最优的比对结果，需要控制-max_target_seqs , -num_alignments 和 -max_hsps 选项：

-max_target_seqs <Integer, >=1>Maximum number of aligned sequences to keepNot applicable for outfmt <= 4* Incompatible with: num_descriptions, num_alignments 
-num_alignments <Integer, >=0>Number of database sequences to show alignments for* Incompatible with: max_target_seqs

---

分割NR子库

NCB blast-2.8版本可支持用NCBI自带代码分割的NR子库的索引作为比对的库，使用比较方便

NR库也要重新下载了ftp://ftp.ncbi.nlm.nih.gov/blast/db/v5/

如果只想比对到单一物种人9606

blastp –db nr –query query.fasta –taxids 9606 –outfmt 6 –out blast.outfm6

比对NR子库哺乳动物的话，需要先建个哺乳动物子库tax_id索引

get_species_taxids.sh -t 40674 > 40674.txids

将序列比对至NR哺乳动物子库

blastp –db nr –query query.fasta –taxidlist 40674.txids –outfmt 6 –out blast.outfm6