来源于:https://github.com/zhangrengang/TEsorter
TE的鉴定可以使用LTRreviewer 或/和 repeatMasker,此处是对鉴定的TE进行分类和可视化。
实质是对TE的进一步分类,因为目前的TE分类比较粗浅Copia和Gypsy,实质上可以分的更加详细。
TEsorter只是分类器,其实是基于dante开发的。分类的依据数据库可以是REXdb或GyDB,REXdb的分类更加详细,所以TEsorter作者推荐使用这个库。
dante本身就是用于基于domain的转座子注释工具。dante github
不同数据库对Copia的分类,来源于REXdb论文
不同数据库对Gypsy的分类
image.png
Copia包括这些小类Ale , AlesiaAlesia, AngelaAngela, BiancaBianca, BrycoBryco, LycoLyco, Gymco I–IV, IkerosIkeros, IvanaIvana, OsserOsser, SIRESIRE, TARTAR, and Tork clades and the
Gypsy包括CRM, ChlamyvirChlamyvir, GaladrielGaladriel, Tcn1Tcn1, ReinaReina, TekayTekay, AthilaAthila, Tat I–III, OgreOgre, RetandRetand, PhygyPhygy, and Selgy clades.
TEsorter把Copia和Gypsy分类为如下:
- capsid protein (GAG),
- aspartic proteinase(AP), 或者PROT
- integrase (INT),
- reverse transcriptase (RT)
- RNase H (RNaseH) 或者RH
提取的时候注意使用的数据库不同关键词略有区别rexdb和gydb的域名有些不同:PROT(rexdb)=AP(gydb),RH(rexdb)=RNaseH(gydb)
安装TEsorter
- python >3
-
hmmscan 3.3x: be compatible with HMMER3/f database format. quickly install by
conda install hmmer
-
blast+: quickly install by
conda install blast
- TEsorter:
git clone https://github.com/zhangrengang/TEsorter
cd TEsorter
python setup.py install
直接使用conda安装也可以conda install -c bioconda tesorter
测试TEsorter
TEsorter TEsorter/test/rice6.9.5.liban
使用参数讲解
- -p 24 cpu数量
- -db rexdb-plant 植物数据库 ,也可以是gydb
- -st nucl 默认是nucl,可选值是:nucl,prot 分别表示DNA nuclear和蛋白质protein。
- -pre prefix 指定输出文件的前缀
从基因组中提取 TE 序列用于 TEsorter
当您只有基因组序列时,以下是从广泛使用的软件的输出中提取 TE 序列的示例。
- 从RepeatMasker输出中提取所有 TE 序列:
#基因组RepeatMasker的下游分析,TE的分类
#指定RepeatMasker的输出文件所在的路径里的基因组.out文件
genome_out="/share/home/repeatsequence/Repeat_all_result/C.s.genome.fa.out"
#基因组位置
genome="/share/home/repeatsequence/Repeat_all_result/C.s.genome.fa"
#指定输出的前缀
prefix="C.s"
pwd="/share/home/TEsort"
cd ${pwd}
if ! [ -e Genome_TEsort ];then
mkdir Genome_TEsort
fi
cd Genome_TEsort
# run RepeatMasker, which will generate a *.out file.
#RepeatMasker [options] genome.fa
# extract sequences
RepeatMasker.py out2seqs $genome_out $genome > ${prefix}.whole_genome_te.fa
#获取分类信息(计算非常快)
TEsorter ${prefix}.whole_genome_te.fa -db rexdb-plant -p 20
#提取domains,此处是提取所有的。 也可以只提取一种,例如GAG 或PROT,此处只提取保守的三种RH,RT,INT
concatenate_domains.py ${prefix}.rexdb-plant.cls.pep INT > ${prefix}.rexdb-plant.cls.pep.INT.aln
concatenate_domains.py ${prefix}.rexdb-plant.cls.pep TPase > ${prefix}.rexdb-plant.cls.pep.TPase.aln
cat ${prefix}.rexdb-plant.cls.pep.INT.aln ${prefix}.rexdb-plant.cls.pep.TPase.aln > ${prefix}.rexdb-plant.cls.pep.INT_TPase.faa
mafft --auto ${prefix}.rexdb-plant.cls.pep.INT_TPase.faa > ${prefix}.rexdb-plant.cls.pep.INT_TPase.aln
#构建进化树
iqtree2 -s ${prefix}.rexdb-plant.cls.pep.INT_TPase.aln -B 1000 -T 20 -m MFP --bnni
#进化树的可视化(可视化可能会失败)
LTR_tree.R ${prefix}.rexdb-plant.cls.pep.INT_TPase.aln.treefile ${prefix}.rexdb-plant.cls.tsv ${prefix}.rexdb-plant.cls.pep.INT_TPase.aln.pdf
- 从LTR_retriever输出中提取所有完整的 LTR-RTs 序列:
run LTR_retriever, which generate two *.pass.list files.
LTR_retriever -genome genome.fa [options]
自动化脚本修改前3行的参数即可运行
#用于LTR-RTs的下游分析,TE的分类
LTR-RTs_genome="/share/home/LTR/C.s.genome.chr.fa" #指定LTRreviewer的输出文件所在的路径里的基因组文件
prefix="C.s" #指定输出的前缀
if ! [ -e TEsort ];then
mkdir TEsort
fi
cd TEsort
LTR_retriever.py get_full_seqs ${LTR-RTs_genome} >${prefix}.ltr.fa
#获取LTR-RTs的分类信息(计算非常快)
TEsorter ${LTR-RTs_genome}.ltr.fa -db rexdb-plant -p 20
#提取domains,此处是提取所有的。 也可以只提取一种,例如GAG 或PROT或者RH
concatenate_domains.py ${prefix}.ltr.fa.rexdb-plant.cls.pep RH RT INT > ${prefix}.ltr.fa.rexdb-plant.cls.pep.full.aln
#构建进化树
iqtree2 -s ${prefix}.ltr.fa.rexdb-plant.cls.pep.full.aln -B 1000 -T 20 -m MFP --bnni
#进化树的可视化(可视化可能会失败)
LTR_tree.R ${prefix}.ltr.fa.rexdb-plant.cls.pep.full.aln.treefile ${prefix}.ltr.fa.rexdb-plant.cls.tsv ${prefix}.ltr.fa.rexdb-plant.cls.pep.full.aln.treefile.pdf
domains的不同种类的保守型不一致,RH RT INT的保守性比较好,所以一般用来构建进化树。GAG和PROT的保守性就没有前3者的稳定。来源于TEsorter和REXdb论文.
支持对 TE 多蛋白序列(示例)或基因蛋白序列进行分类:
TEsorter RepeatPeps.lib -st prot -p 20
从1.4版本可以直接使用基因组来鉴定
TEsorter genome.fasta -genome -p 20
输出结果
rice6.9.5.liban.rexdb.domtbl HMMScan raw output
rice6.9.5.liban.rexdb.dom.faa protein sequences of domain, which can be used for phylogenetic analysis.
rice6.9.5.liban.rexdb.dom.tsv inner domains of TEs/LTR-RTs, which might be used to filter domains based on their scores and coverages.
rice6.9.5.liban.rexdb.dom.gff3 domain annotations in `gff3` format
rice6.9.5.liban.rexdb.cls.tsv TEs/LTR-RTs classifications
Column 1: raw id
Column 2: Order, e.g. LTR
Column 3: Superfamily, e.g. Copia
Column 4: Clade, e.g. SIRE
Column 5: Complete, "yes" means one LTR Copia/Gypsy element with full GAG-POL domains.
Column 6: Strand, + or - or ?
Column 7: Domains, e.g. GAG|SIRE PROT|SIRE INT|SIRE RT|SIRE RH|SIRE; `none` for pass-2 classifications
rice6.9.5.liban.rexdb.cls.lib fasta library for RepeatMasker
rice6.9.5.liban.rexdb.cls.pep the same sequences as `rice6.9.5.liban.rexdb.dom.faa`, but id is changed with classifications.
根据需要提取对应的domain来进行比对后构建进化树
提取的时候注意使用的数据库不同关键词略有区别rexdb和gydb的域名有些不同:PROT(rexdb)=AP(gydb),RH(rexdb)=RNaseH(gydb)
我此处使用的是rexdb,所以提取所有domain的时候使用的是PROT和RH
网友评论