在以往通用的hmmsearch中，经过一轮hmmsearch（round1）和hmmbuild之后再进行hmmsearch（round2），这样产生的结果round2比round1数量多了一倍，而且假阳性居多，并且e值最显著的序列也有可能是假阳性的错误序列（intreproscan可检查），因此在进行hmmsearch的时候，需要严格的过滤，尤其是在round2过程中，但是考虑到e值最显著的序列也有可能是假阳性结果，那么round2再怎么过滤也无济于事，只能说明是query出了问题。所以关键步骤在于round2的时候使用的query，既然常规操作产生了假阳性，那么在round2的时候就应该舍弃掉冗余序列进行hmmsearch，接下来按照此思路进行实践，看看最后结果如何：
通用的hmmsearch流程代码（翻车代码，作为参考）：
hmmsearch --tblout $i.pk --cut_tc PF00~.hmm.gz $i.pep.fa > $i
grep -v "#" $i.pk | cut -d " " -f1 > $i.ID.txt
perl extract_seq_by_ID.cgi $i.ID.txt $i.pep.fa > $i.ugt.fa
clustalw2 -INFILE=line.fa -TYPE=PROTEIN -OUTFILE=align.aln
python covert.py
hmmbuild --amino new.hmm align.sto
hmmsearch --tblout $i.e.pk -E 0.01 new.hmm $i.pep.fa > $i
grep -v "#" $i.e.pk | cut -d " " -f1 > $i.e.ID.txt
perl extract_seq_by_ID.cgi $i.e.ID.txt $i.pep.fa > $i.e.fa
产生的结果$i.e.fa文件，将序列随机挑选去interproscan验证，发现假阳性居多，因此检测版本重点关注在构建hmm文件的query序列上，原版本将round1的整个序列进行构建，验证版本使用domina的序列，也就是round1结果的部分序列，最后查看结果，并验证假阳性。
修改过后的代码（OK版本）：
hmmsearch --domtblout $i_round1_hmm.out --cut_tc 03hmmFile/PF00~.hmm $i.pep.fa > $i_round1
###主要区别在这各部分，将domina提出，而不是整个序列用于下游分析
grep -v '#' $i_round1_hmm.out | awk 'BEGIN{OFS="\t"} $7 < 1e-10 {print $1,$18,$19}' | sort | uniq > $i_domin.bed
seqtk subseq $i.pep.fa $i_domin.bed > $i_domin.fa
###变动部分结束，后面流程一样
cat domin1.fa domin2.fa domin3.fa > combine_domain.fa
###如果有多个数据的话可以合并domina
clustalw2 -INFILE=$i_domain.fa -TYPE=PROTEIN -OUTFILE=domain_align.aln
hmmbuild --amino domain_align.hmm domain_align.aln
hmmsearch --domtblout $i_round2_hmm.out -E 0.01 domain_align.hmm $i.pep.fa > $i_round2
最后拿到的结果文件中可以将目标$i_round2_hmm.out中的第一列target_name拿去与已知数据库或者interproscan或者NCBI进行验证，看是否是目标基因。
————————————————————————————
代码只展示关键步骤，不做流程化处理，搬运时注意替换路径、文件名等~
