RepeatMasker基因组重复序列检测工具安装及使用

一. RepeatMasker简介:

基因组组装完成后,进行基因预测和注释。由于基因组中存在重复序列结构区,特别是高等真核生物,重复序列占了相当大的比例,会影响基因预测的质量,也会带来不必要的资源消耗,因此在基因预测前,首先要检测并屏蔽基因组中的重复序列。

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

那么,这时候就得使用检测重复序列的工具了,先来学习RepeatMasker。

RepeatMasker是重复序列检测的常用工具,通过与数据库相似性比对来准确识别或屏蔽基因组中的重复序列,属于同源预测注释的方式,依赖于RepBase数据库和第三方比对软件(RMBlast、Crossmatch、HMMER、ABBlast),是重复序列注释的重要软件。

二.安装:

首先是官网:

http://www.repeatmasker.org/RepeatMasker/

通过WSL在个人电脑上安装:

Conda install -c bioconda trf

Conda install -c bioconda rmblast

下载RepeatMasker-4.1.2-p1.tar.gz

tar -zxvf RepeatMasker-4.1.2-p1.tar.gz

./configure后,根据提示指定TRF和rmblast所在位置

在Libraries里会生成一个RepeatMaskerLib.h5文件,这就是Dfam+RepBase的整合。

三.运行:

首先确定数据库中是否收录了目标物种:

一些教程是利用./util/queryRepeatDatabase.pl -tree来查看,但我没有找到queryRepeatDatabase.pl文件;

通过“Libraries/taxonomy.dat”查看物种信息,所有已收录物种的名称都存储在该文件中,暂时找不到工具正确打开,按照文本文件去去打开后会有很多乱码,查询效果不佳,希望有了解的同学、老师看见后能够告诉我;

命令行(参考一篇Cell文章):

RepeatMasker -nolow -no_is -norna -parallel 4 -lib /Libraries/RepeatMasker.lib genomic.fna -dir /zhushi/repeatmasker > repeatmasker.log >&1 &

参数详解:

- nolow,不掩盖低复杂度的DNA或简单的重复序列;

-no_is,跳过细菌插入元素检验;

-norna,不掩盖小RNA(伪)基因;

-parallel 4,并行计算,多线程;

-lib [filename] Allows use of a custom library (e.g. from another species)

-dir, 输出结果目录;

参考的命令行并没有-species参数,且-species和-lib参数冲突,只能设置其中一个,理解查询了一下对这两个参数还是有些不明白;

四.结果解读:

大约32h后,程序结束,非常慢,下一次加线程数,结果生成如下四个文件:

.*cat

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

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

.*tbl文件, 重复序列统计文件:

包含了基因组长度、GC含量、重复区长度以及重复区各类别基本统计信息等。其中,“bases masked”就是重复序列的总长度和在基因组中的占比,视物种而定,一般都是比较可靠的。

*.out文件, 重复序列统计文件:

一共十五列,每一列含义如下图:

需要注意,部分注释会在ID列会有*标识,表明该区域有两条比对,*标记结果最佳,通常这部分大部分重叠或者是包含关系,可以选择其中一个;

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

生成gff3文件:

可以使用RepeatMasker提供的工具将结果转化为gff3, 但是需要注意,这里的gff3不是标准的gff3

RepeatMasker/perl rmOutToGFF3.pl *.fa.out > *.gff3

*.masked文件, 屏蔽后的序列:

相较于原始输入fasta文件中的序列,masked文件中将其中重复序列部分屏蔽为了N碱基,该文件默认生成。

注意区分,这里有的N碱基是屏蔽的重复序列,而有的N碱基则是原来这个基因组fasta文件中就有的(一般为gap)。

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

*.polyout文件:

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

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

当然,也有很多人不将“Simple_repeat”视作严格的重复序列类型,如果你也这么认为,直接在*.out中过滤掉所有注释为“Simple_repeat”的结果就可以了。

五.报错信息:

一开始我在服务器上安装,分别对以下几个软件、数据库进行配置:

1.Perl和Python3以及TRF非常好安装,基本上在官网上下载后解压就能用;

2.HMMER见我的公众号另一篇文章基因功能预测工具-HMMER的安装

3.Dfam和Repeat(https://www.girinst.org/repbase/update/)数据库:

Dfam数据库内容见我另一篇公众号文章,Dfam-转座元件TEs数据库, 其中RepeatMasker/Libraries目录下已经有了需要的Dfam.h5文件;

Repeat数据库进入DOWNLOAD页面下载,

tar -zxvf RepBaseRepeatMaskerEdition-20181026.tar.gz

将解压出来的:RMRBSeqs.embl和README.RMRBSeqs文件放到RepeatMasker/Libraries目录下;

4. h5py python 模块,pip3 不是pip,因为用pip3进行安装时报了错,就用conda install h5py安装在miniconda3/lib/python3.8/site-packages/目录下(后面证实安装错误)。

然后正式开始安装RepeatMasker:

tar -zxvf RepeatMasker-4.1.2-p1.tar.gz

./configure后,根据提示指定TRF和HMMER所在位置,然后开始报错了:

将miniconda3/lib/python3.8/site-packages/目录下的h5py文件夹复制到RepeatMasker文件夹下,再次perl ./configure,最后提示如下:

虽然有黄色的cannot提示,先接着往下运行,运行时报错:

运行命令:RepeatMasker -nolow -no_is -norna -parallel 4 -lib /home/ /RepeatMasker/Libraries/RepeatMasker.lib canu_removedup.fa

报错信息:

RepeatMasker version 4.1.2-p1

Search Engine: HMMER [ 3.3.2 (Nov 2020) ]

Taxonomy::new() needs a path for a famdb file!

at /public/home/lvqiang/software/RepeatMasker/RepeatMasker line 635.

后来我查了一下gitlib上issue栏目,作者这样回复:

还是没有正确安装h5py模块,这次想办法让软件去识别miniconda3/lib/python3.8/site-packages/目录下的h5py模块去,通过修改famdb.py文件,还是没有成功。

错误推测:因为我用的服务器系统安装的是python2, 我在我目录下安装了python3,很有可能是两个版本的python起了冲突,导致一些混乱,以后有网后,找个时间解决一下这个问题。

六.补充:

Repeat数据库下载链接:

https://www.girinst.org/server/RepBase/protected/repeatmaskerlibraries/RepBaseRepeatMaskerEdition-20181026.tar.gz;

参数:

参考:https://www.jianshu.com/p/ffdbedae80fa

七.遗留问题:

-species 和 -lib 参数设置;
Linux系统下h5py安装模块失败;
欢迎各位老师、同学针对我的遗留问题进行指导及讨论,在此谢过。

  • 2
    点赞
  • 8
    收藏
    觉得还不错? 一键收藏
  • 3
    评论
以下是一个基于Linux系统的脚本示例,用于使用RepeatMasker软件对鸭基因组进行重复序列检测并筛选出基序长度在2到6个碱基对之间的重复序列: ```bash #!/bin/bash # 设置RepeatMasker路径和其他参数 repeatmasker_path="/path/to/RepeatMasker" genome_file="/path/to/duck_genome.fasta" output_dir="/path/to/output_directory" # 创建输出目录 mkdir -p "$output_dir" # 运行RepeatMasker进行重复序列检测 $repeatmasker_path -species aves -dir "$output_dir" "$genome_file" # 过滤RepeatMasker结果,保留基序长度在2到6个碱基对之间的重复序列 awk -v min_length=2 -v max_length=6 '/^>/ {if (length(seq) >= min_length && length(seq) <= max_length) print header "\n" seq; seq=""; header=$0} !/^>/ {seq=seq $0} END {if (length(seq) >= min_length && length(seq) <= max_length) print header "\n" seq}' "$output_dir/$genome_file".out > "$output_dir/filtered_repeats.fa" # 删除生成的中间文件 rm "$output_dir/$genome_file".out ``` 请注意,上述脚本示例假设你已经正确安装RepeatMasker软件,并将其路径设置为`repeatmasker_path`变量。还需要将`genome_file`变量设置为鸭基因组文件的路径,将`output_dir`变量设置为你想要的输出目录路径。 脚本将运行RepeatMasker对鸭基因组进行重复序列检测,并将结果保存在指定的输出目录中。然后,它会使用awk命令过滤结果,只保留基序长度在2到6个碱基对之间的重复序列,并将结果保存在`filtered_repeats.fa`文件中。 请根据实际情况修改脚本中的路径和参数,并确保RepeatMasker软件正确安装并在系统路径中可用。

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 3
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值