「基因组注释」构建重复序列数据库

本文参考自Repeat Library Construction-Advanced,整体思路一致,但是所用软件有所不同。

流程主要分析MITE和LTR,先根据其结构特征进行注释,之后根据同源信息进行注释,最后进行整合。

主要用到如下软件:

涉及如下几个环境变量,可以根据不同项目进行更改

  • REFERENCE: 用于注释的基因组序列文件,FASTA格式
  • SPECIES: 物种名
  • THREADS: 线程数

MITE

使用MITE-Hunter鉴定Miniature inverted TEs (MITEs),使用方法阅读「基因组注释」MITE-Hunter鉴别基因组的MITE序列

perl MITE_Hunter_manager.pl \
  -i $REFERENCE \
  -g $SPECIES \
  -n 5 \
  -P 1 \
  -S 12345678 \
  -c $THREADS &

将这一步输出文件的"Step8.*fa"和"Step8_singlet.fa"进行合并,作为潜在MITE序列,命名为MITE.lib

cat *Step8.*fa *Step8_singlet.fa > MITE.lib

手工检查候选MITE中的TSD和TIR,将模棱两可的TSD和TIR归为未知序列。

LTR retrotransposons

在植物基因组中的所有重复序列中,LTR逆转座子是其中比例最高的一类结构,因此我们需要尽可能得到高可信的LTR信息。

这一步使用了LTR_retriever分析流程,它整合LTRharverstLTR_FINDER的输出结果,然后得到更可信的LTR-RT序列。

#LTRharvest
gt suffixerator \
  -db $REFERENCE \
  -indexname $SPECIES \
  -tis -suf -lcp -des -ssp -sds -dna
gt ltrharvest \
  -index  $SPECIES \
  -similar 85 -vic 10 -seed 20 -seqids yes \
  -minlenltr 100 -maxlenltr 7000 -mintsd 4 -maxtsd 6 \
  -motif TGCA -motifmis 1  > ltr.harvest.scn &
# LTR_FINDER
ltr_finder -D 15000 -d 1000 -L 7000 -l 100 -p 20 -C -M 0.9 $REFERENCE > ltr.finder.scn &
LTR_retriever -genome $REFERENCE  -inharvest ltr.harvest.scn -infinder ltr.finder.scn -threads $THREADS 

重命名$REFERENCE.LTRlib.fa

mv $REFERENCE.LTRlib.fa LTR.lib

,原教程步骤如下:

  1. 使用LTRharvest收集候选的LTR序列
  2. 使用LTRdigest寻找有PPT(poly purine tract)或PBS(primer binding site)的序列
  3. 对候选序列进行过滤,包括局部串联重复(如着丝粒重复),近期基因重复引起的局部基因聚类,两个不同转座因子距离处于邻近位置
  4. 识别嵌套插入的序列
  5. 建立代表性的LTR序列

RepeatModeler 鉴定其余重复序列

这一步用RepeatMasker以两步构建的文库作为额外的数据库对基因组进行重复序列屏蔽,然后用RepeatModeler从头鉴定基因组中的重复序列。

如果你输入基因组过大,那么建议将其拆分成多个小文件,然后每个文件都单独用RepeatMasker屏蔽序列,最后进行合并。

cat $REFERENCE.LTRlib.fa MITE.fa > MITE_LTR.lib
RepeatMasker -lib MITE_LIB.lib -dir . $REFERECE

输出文件是$REFERENCE.masked

将其中以N标记的重复序列(或基因组上的gap)都删掉

tr -d 'nN' < $REFERENCE.masked | seqkit seq > rmmasked.fa

用RepeatModeler鉴定其余的重复序列

BuildDatabase -name rmdb -engine ncbi rmmasked.fa
nohup /opt/biosoft/RepeatModeler-open-1.0.11/RepeatModeler -data rmdb >& um.out  &

输出结果是consensi.fa.classified, 在RM_xxx文件下

进一步可以按照标识符中是否为unknown将consensi.fa.classified进行拆分

seqkit grep -nrp 'Unknown' consensi.fa.classified > repeatmodeler_unknowns.fasta
seqkit grep -vnrp 'Unknown' consensi.fa.classified > repeatmodeler_identities.fasta

在转座酶数据库中搜索repeatmodeler_unknowns序列,如果能够匹配,则归到repeatmodeler_identities

wget http://www.hrt.msu.edu/uploads/535/78637/Tpases020812.gz
gunzip 
makeblastdb -in Tpases020812 -dbtype prot -out Tpases020812

blastx -query repeatmodeler_unknowns.fasta -db Tpases020812 -evalue 1e-10 -num_descriptions 10 -out modelerunknown_blast_results.txt

perl transposon_blast_parse.pl --blastx modelerunknown_blast_results.txt --modelerunknown repeatmodeler_unknowns.fasta

输出结果

  • identified_elements.txt
  • unknown_elements.txt
mv  unknown_elements.txt  ModelerUnknown.lib
cat  identified_elements.txt  repeatmodeler_identities.fasta  > ModelerID.lib

过滤基因片段

到目前位置,我们已经构建了如下重复序列数据库

  • MITE.lib: MITE重复序列数据库
  • LTR.lib: LTR重复序列数据库
  • ModelerID.lib: 已知分类重复数据库
  • ModelerUnknown.lib: 未知分类重复数据库

如果为了进一步提高重复序列的可靠性,可以将上述序列分别和植物蛋白数据库进行比对。以ModelerUnknown.lib为例

http://www.hrt.msu.edu/uploads/535/78637/alluniRefprexp070416.gz
gunzip alluniRefprexp070416.gz
makeblastdb -in alluniRefprexp070416 -dbtype prot -out alluniRefprexp070416
blastx -query ModelerUnknown.lib -db alluniRefprexp070416  -evalue 1e-10 -num_descriptions 10 \
-num_threads 20 -out ModelerUnknown.lib_blast_results.txt

然后用ProtExcluder进行过滤

/opt/biosoft/ProtExcluder1.1/ProtExcluder.pl  -f 50 ModelerUnknown.lib_blast_results.txt ModelerUnknown.lib

输出结果是Modelerunknown.libnoProtFinal

最后,MITE.lib, LTR.lib 和 ModelerID.lib 合并成 KnownRepeats.lib。 KnownRepeats.lib 和 Modelerunknown.lib 合并成 allRepeats.lib。

KnownRepeats.lib准确性较高,但是不一定有新的重复序列家族,allRepeats.lib 全面但不一定准。

参考文献

Campbell, M. S., Law, M., Holt, C., Stein, J. C., Moghe, G. D., Hunagel, D. E., Lei, J., Achawanantakun, R., Jiao, D., Lawrence, C. J., Ware, D., Shiu, S-H., Childs, K. L., Sun, Y., Jiang, N, Yandell, M. 2014. MAKER-P: a tool-kit for the rapid creation, management, and quality control of plant genome annotations. Plant Physiology 164 513-524.

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 159,015评论 4 362
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 67,262评论 1 292
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 108,727评论 0 243
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 43,986评论 0 205
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 52,363评论 3 287
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 40,610评论 1 219
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 31,871评论 2 312
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 30,582评论 0 198
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 34,297评论 1 242
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 30,551评论 2 246
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 32,053评论 1 260
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 28,385评论 2 253
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 33,035评论 3 236
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 26,079评论 0 8
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 26,841评论 0 195
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 35,648评论 2 274
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 35,550评论 2 270

推荐阅读更多精彩内容