基因组重复序列检测:RepeatMasker的安装及使用

RepeatMasker是重复序列检测的常用工具,通过与参考数据库的相似性比对来准确识别或屏蔽基因组中的重复序列,属于同源预测注释的方式。下文除了RepeatMasker的使用,也提到了一些在安装和使用过程中可能会碰到的问题。

基因组组装完成后,进行基因预测和注释。由于基因组中存在重复序列结构区,特别是高等真核生物,重复序列占了相当大的比例,会影响基因预测的质量,也会带来不必要的资源消耗。因此在基因预测前,首先要检测并屏蔽基因组中的重复序列(真核生物很必要,原核生物可忽略)。对于这种情况,只是为了屏蔽重复序列降低后续基因预测的压力,重复序列的检测可以不必太纠结。

但有时候我们注释基因组的重复序列结构,也可能是专注于某些特定研究,例如,某些重复元件可能参与了重要功能,我们期望定位它们的位置。这种情况下,需要识别精准。

RepeatMasker提供了在线服务器,将基因组序列上传后即可在线运行。
RepeatMasker在线链接:http://www.repeatmasker.org/cgi-bin/WEBRepeatMasker/

RepeatMasker本地安装及使用:在正式安装RepeatMasker主程序之前,需要提前配置好其依赖的工具。
这些依赖的工具可以手动安装,也可以用conda 自动安装。

RepeatMaskr本地配置:

1、TRF
TRF(Tandem Repeats Finder)在RepeatMasker安装前必需要提前配置好,这个可以使用conda直接安装。

conda install TRF

2、序列搜索引擎(4选1即可)
RepeatMasker需要序列搜索引擎,以实现输入基因组序列和参考库中序列的比对,在RepeatMasker安装前必需要提前配置好。RepeatMasker主要通过以下4种工具实现序列比对功能,CrossMatch、RMBlast、ABBlast、HMMER,因此我们需要从中至少选1个。
选择了RMBlast,这个也可以使用conda直接安装。

conda install RMBlast

3、RepeatMasker
前两个工具配置好后,正式安装RepeatMasker,conda搞不定了,需要源码编译。

wget http://www.repeatmasker.org/RepeatMasker-open-4-0-6.tar.gz
tar xzvf RepeatMasker-open-4-0-6.tar.gz
cd RepeatMasker
chmod -R 755 *
./configure   #执行后,根据提示信息一步步来。

首先是perl环境,默认自动检测,或者手动更改perl主程序路径后回车继续。

1.png

然后是TRF,默认自动检测,或者手动更改TRF主程序所在路径后(可通过which trf获取,比方说我conda安装的trf程序所在“/home/my/software/Miniconda3/bin”),回车继续。

2.png

再然后是序列搜索引擎,因为刚才安装的是RMBlast,所以我们这里选择2。然后在接下来的新界面中,将RMBlast主程序所在路径输入后(可通过which rpsblast获取,比方说我conda安装的rpsblast程序所在“/home/my/software/Miniconda3/bin”),回车返回主界面后,再选择5,就完成了。

你也可以指定多种序列搜索引擎后,再选择5,不过实际运行时,一次只能选择一种序列比对方式。

3.png
4.png

RepeatMasker提示安装完成后,最后配置环境变量。

#例如,我的RepeatMasker安装路径是在“/home/my/software/RepeatMasker”
export PATH=/home/my/software/RepeatMasker:$PATH
 
#这时候没啥问题的话应该可以看到帮助界面了
RepeatMasker -h
5.png

如果提示有perl模块未配置好,通过cpan命令安装相应的perl模块即可。

如果反复提示以下关于“Text::Soundex module”模块的错误(即便你这个perl模块确实安装好了),建议将perl环境更改为/usr/bin下的perl(对,不建议使用conda中的perl),并sudo cpan重新安装该模块后,再重新安装RepeatMasker(安装时指定系统/usr/bin下的perl),就解决了。

4、数据库文件(Repbase)

Repbase,是遗传信息研究所发布的重复DNA数据库,收录了非常多的物种基因组的重复序列信息。在这里,我们将通过它对自己的物种基因组的重复序列进行检测。

登录REPBASE数据库(https://www.girinst.org/server/RepBase/index.php),需注册才能下载,非营利性组织可以免费使用,人工审批,需要等待1-2天时间。在官网上直接下载的肯定是最新版的。

6.png

当然如果你懒得注册等待,随便在网上一搜,现成的Repbase资源也有很多可直接下载,基本上都是个人分享的。缺点就是这些资源的时间可能比较久远,版本较老数据库不如现在的版本全是肯定的。

Repbase:百度网盘分享
链接:https://pan.baidu.com/s/1gNKFvcNLgGkafX1xv_PrMw
提取码:9v62

在RepeatMasker安装完毕后,将Repbase库解压放置到RepeatMasker安装路径下的“Libraries”中就可以了。

RepeatMasker使用:

接下来,我们使用某真菌基因组(真核生物,且基因组小,方便演示),测试RepeatMasker使用。
基因组获取自NCBI,Assembly:GCA_900382705.2,物种Fusarium tricinctum。

首先确定数据库中是否收录了目标物种
(一)Repbase配置好后,直接通过安装目录中的“RepeatMasker/util/queryRepeatDatabase.pl”,即可查看当前库中已存在的物种。

#查看物种列表
cd RepeatMasker
./util/queryRepeatDatabase.pl -tree > species.txt
less -S species.txt
 
#目标物种,这里为 Fusarium tricinctum,我们可直接查看 Fusarium 属所有的
grep 'Fusarium' species.txt

(二)或者参看“RepeatMasker/Libraries/taxonomy.dat”中的物种信息,所有已收录物种的名称都存储在该文件中。

#查看物种列表
cd RepeatMasker
./util/queryRepeatDatabase.pl -tree > species.txt
less -S species.txt
 
#目标物种,这里为 Fusarium tricinctum,我们可直接查看 Fusarium 属所有的
grep 'Fusarium' species.txt

(三)还有一种方式,我们大可以直接在运行的时候,让软件自己判断了。
也就是在正式使用RepeatMasker时(见下文使用示例),通过“-species”参数指定物种。对于常规的模式物种或类群,可以直接输入通用的名称;更保险的方式是在“-species”参数后接物种的拉丁名称,由于属名和种名之间有空格,所以需要用引号括起来。如果RepeatMasker没有提示异常,即库中存在目标物种,则将会正常运行,等待结果就可以了。

7.png

若无法在现有的Repbase库中找到一个合适的参考集,除了使用近缘物种代替外,还可能需要考虑基于当前的全基因组序列,执行重复序列的denovo预测,如RepeatModeler等工具,训练重复序列集构建本地的repeat library,运行RepeatMasker时通过“-lib”参数指定本地的repeat library。

此外,也并不是说找到参考物种了,就是完全可行了。参考库中的信息可能并不全,没有将目标物种基因组覆盖完全,由此会导致基因组重复序列识别的较少。可见下文的示例输出结果部分,以一个示例提到了这个问题。

RepeatMasker使用示例

我们的物种数据库中存在,就可以直接拿它的重复序列作参考了。

#查看帮助,各参数详情不再细说了,大家自行了解
RepeatMasker -h
 
#一个简单的示例
#基因组获取自:https://www.ncbi.nlm.nih.gov/genome/68823?genome_assembly_id=380814
 
RepeatMasker -pa 4 -species "Fusarium tricinctum" -poly -html -gff -dir repeat1 GCA_900382705.2_FTRI.INRA104.GCA2018.2_genomic.fna 1>log.o.txt 2>log.e.txt
 
#需要注意的地方
#-dir 指定的输出结果路径,必须提前建立好,否则无结果
#一定要通过 -species 指定物种,否则默认比对的是人类重复序列数据库
#如果使用本地的参考库,通过 -lib 指定,替代 -species
#-s、-q、-qq 等参数可控制序列比对的灵敏度,如果你的目标物种和参考物种不是很近,可能需要提升灵敏度

这个真菌基因组不大,所以很快就完事了。
其它常见问题,可参阅官方文档:http://www.repeatmasker.org/webrepeatmaskerhelp.html

主要结果

输出结果细节,可参阅官方文档:http://www.repeatmasker.org/webrepeatmaskerhelp.html
以下简要展示上述的示例命令所得的7个结果文件。

8.png

*.tbl
本次RepeatMasker运行的结果报告,该文件默认生成,包含了基因组长度、GC含量、重复区长度以及重复区各类别基本统计信息等。其中,“bases masked”就是重复序列的总长度和在基因组中的占比,视物种而定,一般都是比较可靠的。

9.png

然而这儿就碰到问题了。Fusarium tricinctum是一种霉菌,重复序列的比例不可能这么少(仅0.87%)。主要的原因可能是对该物种的研究相对较少,或者使用的数据库版本较老,数据库中该物种的重复序列信息不全(实际上,我使用了一个非常旧的Repbase库)。可以尝试使用最新版本的Repbase库,同时多找几个近缘物种的参考库作比对,更改RepeatMasker参数把灵敏度调高(或者将RMBlast替换为更灵敏的工具)等,尽可能找到更多的重复序列。如果最大原因仍是数据库本身不全所致,其实更改参数之类的提升效果也甚微了。此时可能需要考虑执行重复序列的denovo预测,如RepeatScout、RepeatModeler等工具,训练重复序列集构建repeat library,再通过RepeatMasker注释重复序列。

*.cat

记录了输入的基因组序列和数据库中参考重复序列的比对详情,该文件默认生成。

会存在个别碱基的差异,其中“i”和“v”分别代表了碱基转换(transitions)和颠换(transversions),“-”表示该位点存在碱基插入/删除。

10.png

*.out****和*****.out.html

这两个文件中的信息是一致的,将基因组中预测得到的重复序列和参考序列相比的碱基替换频率、插入/删除率,以及重复序列的位置、结构、类型等信息展示出。其中,.out(默认生成)是纯文本样式,.out.html(-html参数生成)以网页列表展示。

11.png

通过这个*.out文件(或下文的gff文件),就可以去定位你期望关注的特殊类型的重复序列元件在基因组中的位置了,后续可再自写脚本根据位置信息将这段序列提取出来,或者更进一步研究它们的功能等。

*.out.gff

将*.out中的内容,整理为标准gff文件的结构类型(-gff参数生成)。主要包含重复序列的位置、结构等信息。

12.png

*.polyout

命令行中通过-poly参数,可额外将预测结果.out中的微卫星注释识别出来,单独整理为一张表,文件结构同.out。

13.png

如果你不想将微卫星视作严格的重复序列类型,可通过.polyout文件中的注释位置,将.out中的微卫星去除。那么,为什么不直接在.out中根据注释作筛选呢?因为微卫星属于“Simple_repeat”的一种,但.out中并非所有“Simple_repeat”都是微卫星,所以直接去筛选很难操作。当然,也有很多人不将“Simple_repeat”视作严格的重复序列类型,如果你也这么认为,直接在*.out中过滤掉所有注释为“Simple_repeat”的结果就可以了。

.masked
相较于原始输入fasta文件中的序列,
..masked中将其中重复序列部分屏蔽为了N碱基,该文件默认生成。
注意区分,这里有的N碱基是屏蔽的重复序列,而有的N碱基则是原来这个基因组fasta文件中就有的(一般为gap)。

14.png

那么这个文件有什么作用呢?可以将该文件作为后续基因序列预测的输入文件。这样,基因预测时就不会再考虑这些重复序列区域(因为这些区域的碱基已经屏蔽为了N碱基,不会被识别),缩小了基因组范围,大大减少了资源消耗,提升准确度。

本文转载自刘尧科学网博客。
链接地址:http://blog.sciencenet.cn/blog-3406804-1202555.html

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

推荐阅读更多精彩内容