#安装hmmer
conda install hmmer
#下载pfam数据库,这里文件下载到哪里都可,记住库地址就行
wget ftp://ftp.ebi.ac.uk/pub/databases/Pfam/releases/Pfam33.1/Pfam-A.hmm.gz
gzip -d Pfam-A.hmm.gz
得到 PFAM 数据库的 HMM 文件。 HMM 文件是文本文件,需要将其变成二进制格式,以加快运算速度,同时进行压缩,并建立成索引数据库。
hmmpress Pfam-A.hmm
3.使用hmmscan 进行Pfam注释
Pfam 数据库中每个编号代表一个蛋白质家族。Pfam 分 A 和 B 两个数据库,其中 A 数据库是经过手工校正的高质量数据库, B 数据库虽然质量低些,依然可以用来寻找蛋白质家族的保守位点。Pfam 30.0以后的版本,Pfam-B就去除掉了,所以我们只关心PfamA就好了。
#因为我要查看结构域的大小,方便我去人工删除不满足大小的基因(如只包含结构域片段的基因)
hmmscan -o out.txt --tblout out.tbl --domtblout out.dom --noali -E 1e-5 /pfam/Pfam-A.hmm multifile.fasta
附上参数解读
$ hmmscan [-options]
-h
显示帮助信息
-o FILE
将结果输出到指定的文件中。默认是输出到标准输出。
–tblout FILE
将蛋白质家族的结果以表格形式输出到指定的文件中。默认不输出该文件。
–domtblout FILE
将蛋白结构域的比对结果以表格形式输出到指定的文件中。默认不输出该文件。该表格中包含query序列起始结束位点与目标序列起始结束位点的匹配信息。
–acc
在输出结果中包含 PF 的编号,默认是蛋白质家族的名称。
–noali
在输出结果中不包含比对信息。输出文件的大小则会更小。
-E FLOAT default:10.0
设定 E_value 阈值,推荐设置为 1e-5 。
-T FLOAT
设定 Score 阈值。
–domE FLOAT default:10.0
设定 E_value 阈值。该参数和 -E 参数类似,不过是 domain 比对设定的值。
–cpu
多线程运行的CPU。默认应该是大于1的,表示支持多线程运行。但其实估计一般一个hmmscan程序利用150%个CPU。并且若进行并行化调用hmmscan,当并行数高于4的时候,会报错:Fatal exception (source file esl_threads.c, line 129)。这时,设置–cpu的值为1即可。
pfamscan 的使用_【HMMSCAN】使用pfam数据库对多序列文件进行结构域注释
安装要求
64-bit Linux
Perl 5 (default on most Linux distributions)
Python 3 (InterProScan 5.30-69.0 onwards)
Java JDK/JRE version 11 (InterProScan 5.37-76.0 onwards)
具体如何配置安装环境查看官方文档。
conda install interproscan
测试安装是否成功,interproscan.sh在interproscan安装的文件夹下,如需在其他的文件夹下操作,可以使用文件夹的绝对路径
./interproscan.sh -i test_proteins.fasta -f tsv
fasta 格式的蛋白或核酸序列,序列中不能含有 · - 或 * 等非法字符。
一般常用的参数有这些:
-appl,–applications 用于指定使用Interpro中哪些数据库,默认全部数据库
-b,–output-file-base 用于指定输出文件的路径or文件夹,默认是输入文件的路径
-f,–formats 用于指定输出文件的后缀,蛋白序列默认输出TSV, XML and GFF3
-i,–input 输入文件,一般要为fasta格式,不要带有其他特殊符号
下载COG数据库的whog和fun.txt 文件。whog文件包含COG编号和COG数据库中序列名的对应关系,也包含COG编号和25个大类的对应关系;fun.txt是25个大类的描述性信息。我们根据这2个文件的信息来编写程序对Blast的结果进行处理,得到COG注释。
$wget ftp://ftp.ncbi.nih.gov/pub/COG/COG/whog
$wget ftp://ftp.ncbi.nih.gov/pub/COG/COG/fun.txt