Transcription Fatcor (TF), TCDB注释步骤(含批量)

11.Transcription Fatcor

新建文件夹

mkdir 11.TF && cd 11.TF

需要interpro5.tsv

perl /media/aa/DATA2/bin/interpro2tf_for_Fungi.pl /media/aa/DATA/SZQ2/bj/functional_annotation/05.InterProScan/pep70/$i.interpro.tsv --out_prefix TF

批量操作

批量 pep70

for i in `cat /media/aa/DATA/SZQ2/bj/functional_annotation/70list.txt`
do
    echo "perl /media/aa/DATA2/bin/interpro2tf_for_Fungi.pl /media/aa/DATA/SZQ2/bj/functional_annotation/05.InterProScan/pep70/$i.interpro.tsv --out_prefix $i.TF"
done > command.interpro2tf_for_Fungi.list
ParaFly -c command.interpro2tf_for_Fungi.list -CPU 48

批量 pepmy

for i in `cat /media/aa/DATA/SZQ2/bj/functional_annotation/pepmylist.txt`
do
    echo "perl /media/aa/DATA2/bin/interpro2tf_for_Fungi.pl /media/aa/DATA/SZQ2/bj/functional_annotation/05.InterProScan/pepmy/$i.interpro.tsv --out_prefix $i.TF"
done > command.interpro2tf_for_Fungi.list
ParaFly -c command.interpro2tf_for_Fungi.list -CPU 48

14.TCDB 转运蛋白

本地BLASTP比对注释

对于TCDB而言,是可以免费下载的,下载的链接:

TCDB » Downloads

打开TCDB FastASequences,为蛋白序列文件,下载至本地可以使用blastp比对。

在“>”之后的序列标识符中,1001796365代表该蛋白质序列在GeneBank数据库中的编号,4.F.1.1.5 是TC Number, 代表该蛋白质所属的转运蛋白家族,CDP-alcohol phosphatidyltransferase是对该转运蛋白家族功能的具体描述,[Marinobacterexcellens] 是该蛋白序列的来源物种。

进入conda环境

conda activate pfam_scan

1)diamond blastp

diamond blastp --db /media/aa/DATA/SZQ2/TCDB/tcdb --query /media/aa/DATA/SZQ2/bj/functional_annotation/pep70/$i.fasta --out $i.tcdb10.xml --outfmt 5 --sensitive --max-target-seqs 20 --evalue 1e-10 --id 20 --tmpdir /dev/shm --index-chunks 1

批量操作 

批量 pep70

for i in `cat /media/aa/DATA/SZQ2/bj/functional_annotation/70list.txt`
do
    echo "diamond blastp --db /media/aa/DATA/SZQ2/TCDB/tcdb --query /media/aa/DATA/SZQ2/bj/functional_annotation/pep70/$i.fasta --out $i.tcdb10.xml --outfmt 5 --sensitive --max-target-seqs 20 --evalue 1e-10 --id 20 --tmpdir /dev/shm --index-chunks 1"
done > command.tcdb.list
ParaFly -c command.tcdb.list -CPU 48

批量 pepmy

for i in `cat /media/aa/DATA/SZQ2/bj/functional_annotation/pepmylist.txt`
do
    echo "diamond blastp --db /media/aa/DATA/SZQ2/TCDB/tcdb --query /media/aa/DATA/SZQ2/bj/functional_annotation/pepmy/$i.fasta --out $i.tcdb10.xml --outfmt 5 --sensitive --max-target-seqs 20 --evalue 1e-10 --id 20 --tmpdir /dev/shm --index-chunks 1"
done > command.tcdb.list
ParaFly -c command.tcdb.list -CPU 48

2)parsing_blast_result.pl

# 新建文件夹
mkdir tcdb10.tab && cd tcdb10.tab
/media/aa/DATA2/bin/parsing_blast_result.pl --evalue 1e-10 --HSP-num 1 --out-hit-confidence --suject-annotation ../$i.tcdb10.xml > $i.tcdb10.tab

批量操作 

批量 pep70

for i in `cat /media/aa/DATA/SZQ2/bj/functional_annotation/70list.txt`
do
    echo "/media/aa/DATA2/bin/parsing_blast_result.pl --evalue 1e-10 --HSP-num 1 --out-hit-confidence --suject-annotation ../$i.tcdb10.xml > $i.tcdb10.tab"
done > command.tcdb.list
ParaFly -c command.tcdb.list -CPU 48

批量 pepmy

for i in `cat /media/aa/DATA/SZQ2/bj/functional_annotation/pepmylist.txt`
do
    echo "/media/aa/DATA2/bin/parsing_blast_result.pl --evalue 1e-10 --HSP-num 1 --out-hit-confidence --suject-annotation ../$i.tcdb10.xml > $i.tcdb10.tab"
done > command.tcdb.list
ParaFly -c command.tcdb.list -CPU 48

Blast分析后的注释结果中一条序列含有多种结果,后续的分析需要提取得到基因最可能的注释结果。

如“/media/aa/DATA2/bin/gene_annotation_from_SwissProt.pl的操作,如下:

3)比对结果中筛选每个query的最佳subject

进入conda环境

conda activate jcvi
python -m jcvi.formats.blast best -n 1 $i.tcdb10.tab

批量操作 

批量 pep70

for i in `cat /media/aa/DATA/SZQ2/bj/functional_annotation/70list.txt`
do
    echo "python -m jcvi.formats.blast best -n 1 $i.tcdb10.tab"
done > command.jcvi.list
ParaFly -c command.jcvi.list -CPU 48

批量 pepmy

for i in `cat /media/aa/DATA/SZQ2/bj/functional_annotation/pepmylist.txt`
do
    echo "python -m jcvi.formats.blast best -n 1 $i.tcdb10.tab"
done > command.jcvi.list
ParaFly -c command.jcvi.list -CPU 48

4)整理best文件

mkdir best && cd best
cp ../*.tcdb10.tab.best ./

Usage: /home/chenlianfu/chenlianfu_scripts/parsing_blast_result.pl [options] blast.out > blast.tab 对BLAST的xml或tab格式的结果进行解析和过滤,得到更准确的BLAST结果。结果为表格形式(BLAST outfmt6),结果按query序列的ID排序,每个query序列的比对结果按得分排序。 --type default: xml 设置输入BLAST结果文件的类型。可以设置为xml或tab两种类型。 若是tab格式,则BLAST结果中没有query与subject的序列长度信息,默认设置下无法使用--subject-coverage和--query-coverage参数的覆盖率阈值对结果进行过滤。在设置--db-subject输入数据库FASTA文件后可以使用--subject-coverage参数进行过滤;在设置--db-query输入query序列FASTA文件后可以使用--query-coverage参数进行过滤。 若是xml格式,结果文件中包query和subject长度信息,从而不需要使用--db-subject和--db-query参数输入FASTA序列文件。 --no-header 添加该参数则不输出表头。 --max-hit-num default: 20 设置允许的最大hit数量。 --evalue default: 1e-5 设置HSP的evalue阈值。 --identity default: 0.05 设置HSP的identity阈值。 --CIP default: 0.2 设置cumulative identity percentage阈值(这里依然使用了比值,单位不是%,所以其值要设置不大于1,默认值0.2表示20%阈值),对Hit进行过滤。CIP = 所有HSPs的一致位点之和 / 所有HSPs的比对长度之和。 --subject-coverage default: 0.2 设置所有HSPs对subject序列总体的覆盖率阈值。该参数阈值在文献中也被称为CALP(cumulative alignment length percentage),即 sum of all HSPs / subject length。 --db-subject 输入数据库的FASTA文件,以获取subject序列长度信息。 --query-coverage default: 0.2 设置所有HSPs对query序列总体的覆盖率阈值。该参数阈值在文献中也被称为CALP(cumulative alignment length percentage),即 sum of all HSPs / query length。 --db-query 输入query序列的FASTA文件,以获取query序列长度信息。 --percentage-of-top-bitscore default: 100 使用bitscore得分对hit进行过滤,设置输出hits的bitscore得分和最高得分相差不超过最高得分的百分数。hit若有多个HSPs,则取最高的HSP得分作为hit的得分;若数据库非常大,则推荐将设置该参数值设置为10,则能极大减少比对结果,保留最准确的结果;若数据库比较小,则推荐设置该参数值为50,或使用默认值;使用该参数来减少比对结果,优于仅使用最优比对结果。 --HSP-num default: max 若一个hit有多个HSPs,该参数设置输出得分指定数目个最高的HSPs。默认输出所有的HSPs。 --out-hit-confidence 添加该参数,则在表格结果第13、14和15列分别输出Hit的CIP、CALP_query、CALP_subject值。 --suject-annotation 若--type参数的值是xml,添加该参数可以生效,则额外增加最后一列suject annotation注释结果。
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值