基因组结构预测|Braker3流程（基因组结构分析）-技术精选-玖叶教程网

BRAKER(https://github.com/Gaius-Augustus/BRAKER)是一个基因组注释流程，能够结合GeneMark，AUGUSTUS等软件BRAKER pipeline自动注释真核基因组，可以根据我们的数据：组装好的基因组、蛋白序列、转录组选择不同的流程。BRAKER 1 使用短read RNA-Seq在GeneMark-ET、AUGUSTUS进行训练和预测；BRAKER 2使用蛋白质数据在GeneMark-EP、AUGUSTUS中进行训练和预测；BRAKER 3则能够同时结合短读RNA-Seq与蛋白质数据库，在GeneMark-ETP和AUGUSTUS中训练并预测更加可靠的基因。了解一下它的使用吧！

一、软件安装

Braker3依赖的软件比较多，整个流程有十几个注释用的软件，手动安装会比较繁琐，可使用container来安装，以下是安装方法的介绍，仅供参考。

手动安装

https://github.com/Gaius-Augustus/BRAKER

1.1 新建conda环境

# 这是我使用的conda
wget -c https://repo.continuum.io/miniconda/Miniconda3-latest-Linux-x86_64.sh
# 新建braker3环境用于下面perl模块的安装
conda create -n braker3 
conda activate braker3

1.2 安装perl模块

conda install -c anaconda perl
conda install python=3.7
conda install -c anaconda biopython=1.76
conda install -c bioconda perl-app-cpanminus
conda install -c bioconda perl-file-spec
conda install -c bioconda perl-hash-merge
conda install -c bioconda perl-list-util
conda install -c bioconda perl-module-load-conditional
conda install -c bioconda perl-posix
conda install -c bioconda perl-file-homedir
conda install -c bioconda perl-parallel-forkmanager
conda install -c bioconda perl-scalar-util-numeric
conda install -c bioconda perl-yaml
conda install -c bioconda perl-class-data-inheritable
conda install -c bioconda perl-exception-class
conda install -c bioconda perl-test-pod
conda install -c bioconda perl-file-which # skip if you are not comparing to reference annotation
conda install -c bioconda perl-mce
conda install -c bioconda perl-threaded 
conda install -c bioconda perl-math-utils
conda install -c eumetsat perl-yaml-xs 
conda install -c bioconda perl-data-dumper
cpan Thread::Queue
perl -M模块名 -e1 检查是否安装好

1.3 安装Braker流程依赖的软件

Braker所需的samtools在usr/bin/中已装好
----------
diamond 0.9.24 (diamond和blast选择一个下载）
conda install diamond=0.9.24
----------
cdbtools 0.99
conda install -c bioconda cdbtools
----------
hisat2 2.2.1
conda install -c bioconda hisat2
----------
Bamtools 2.5.1 (Augustus依赖)
# 下载
wget https://github.com/pezmaster31/bamtools
# 解压
tar -zvxf bamtools-2.5.1.tar.gz 
# 创建安装文件夹 
mkdir build && cd build 
# 生成Makefile文件
cmake ../
# 编译和安装（为了使Augustus安装时找到bamtools路径，我把它装在bamtools-2.5.1/src/usr/下）
make && make DESTDIR=../src install
# 在~/.bashrc中加入环境变量
export BAMTOOLS_PATH=/home/yangkaiyan/biosoft/bamtools-2.5.1/src/usr/local/bin
----------
samtools (Augustus依赖)
# 下载
wget https://github.com/samtools/samtools
# 解压
tar -jxvf samtools-1.17.tar.bz2
cd samtools-1.17
./configure --prefix=/home/yangkaiyan/biosoft/samtools-1.17
# 编译和安装
make && make install 
# 在~/.bashrc中加入环境变量
export TOOLDIR=/home/yangkaiyan/biosoft/samtools-1.17
# 把samtools-1.17中htslib-1.17改成htslib
----------
MYSQL (Augustus依赖)
# 可以将common.mk中的ZIPINPUT=true注释掉或按照下面方安装载MYSQL
# 下载MYSQL
wget https://tangentsoft.com/mysqlpp/releases/mysql++-3.3.0.tar.gz
# 解压
tar -xvf mysql++-3.3.0.tar.gz
cd mysql++-3.3.0
pwd
# 复制当前路径/home/yangkaiyan/biosoft/mysql++-3.3.0
./configure --prefix=home/yangkaiyan/biosoft/mysql++-3.3.0/build_file
# 编译和安装
make && make install
cd build_file
pwd
# 复制当前路径/home/yangkaiyan/biosoft/mysql++-3.3.0/build_file（后面需要用到）
----------
Augustus 3.4.0 
开发者建议我们不要用conda直接安装Augustus，因为在这个conda环境中只能安装Augustus 3.3.3，
Braker不支持低版本的Augustus，看看如何手动安装吧！
# 所需库文件：gcc/boost/zlib/git(一般都有)
# 检查以下软件是否安装好
MYSQL(上一步已装好)
bamtools(上一步已装好)
htslib/bcftools/tabix/samtools(上一步已装好，都在samtools-1.17中)
# 下载
wget http://bioinf.uni-greifswald.de/augustus/binaries/augustus-3.4.0.tar.gz
# 解压
tar -zvxf augustus-3.4.0.tar.gz
cd Augustus-3.4.0
# 修改src中Makefile的路径（如果common.mk中的ZIPINPUT=true被注释掉，且运行make没有报MYSQL的相关错误，就不用做这一步）
vim ./src/Makefile
:set nu
找到第55行将-I/usr/include/mysql++改成
-I/home/yangkaiyan/biosoft/mysql++-3.3.0/build_file/include/mysql++
找到第36行LIBS += -l，在下一行空格到跟上一行LIBS对齐的位置，并新增
LIBS += -L/usr/include/lpsolve
LIBS += -L/home/yangkaiyan/biosoft/mysql++-3.3.0/build_file/lib
# 修改bam2hints中Makefile的路径
vim ./auxprogs/bam2hints/Makefile
INCLUDES = -I /home/yangkaiyan/biosoft/bamtools-2.5.1/src
LIBS = /home/yangkaiyan/biosoft/bamtools-2.5.1/src/usr/local/lib/libbamtools.a -lz
# 修改filterBam中src/Makefile的路径
vim ./auxprogs/filterBam/src/Makefile
BAMTOOLS = /home/yangkaiyan/biosoft/bamtools-2.5.1/src
LIBS = /home/yangkaiyan/biosoft/bamtools-2.5.1/src/usr/local/lib/libbamtools.a -lz
# 安装
make 
# 在bin目录下
augustus --help或./augustus --help
# 在~/.bashrc中加入环境变量
export $AUGUSTUS_CONFIG_PATH=/home/yangkaiyan/biosoft/augustus-3.4.0/config/
export $AUGUSTUS_BIN_PATH=/home/yangkaiyan/biosoft/augustus-3.4.0/bin/
export $AUGUSTUS_SCRIPTS_PATH=/home/yangkaiyan/biosoft/augustus-3.4.0/scripts/
总结安装方法就是缺少哪些文件就按照报错提示去下载安装，不是/usr/local/bin中的文件都要指定路径
----------
geneMark-ETP 1.0(for braker3)
# 安装perl模块
YAML::XS Data::Dumper Thread::Queue （上一步已装好）
# 下载：
wget https://github.com/gatech-genemark/GeneMark-ETP
# 查看INSTALL文档
# 解压
unzip *.gz
cd 解压后的文件
# 密钥拷贝为/home/yangkaiyan/下的隐藏文件
cp gm_key_64 ~/.gm_key_64
# 在~/.bashrc加入环境变量
export GENEMARK_PATH=/home/yangkaiyan/biosoft/GeneMark-ETP-main/bin
# 测试
perl gmetp.pl --cores 32 --cfg input.cfg
测试方法和流程：https://github.com/gatech-genemark/GeneMark-ETP-exp
基因组和蛋白质数据准备：https://github.com/gatech-genemark/GeneMark-ETP-exp/tree/main/Drosophila_melanogaster/data


如果使用的是braker2，就需要安装GeneMark-ES/ET/EP
GeneMark-ES/ET/EP(for braker2)
# 下载
wget http://exon.gatech.edu/GeneMark/license_download.cgi
# 解压
unzip *.gz
cd 解压后的文件
# 密钥拷贝为/home/yangkaiyan/下的隐藏文件
cp gm_key_64 ~/.gm_key_64
# 在~/.bashrc下加入环境变量
export PATH=$PATH:/home/yangkaiyan/biosoft/gmes_linux_64_4
# 测试
gmes_petap.pl --ES --sequence genome.fa
# 如果报错提示缺少一些perl模块，就按照提示一个一个安装
----------
ProtHint 2.6.0
# 下载
wget https://github.com/gatech-genemark/ProtHint
# 解压
tar -xzvf ProtHint-2.6.0.tar.gz
cd ProtHint-2.6.0
# 安装perl模块
cpan MCE::Mutex threads YAML Thread::Queue Math::Utils （上一步已装好）
# python版本查看 （Python 3.3 or higher is required）
Python --version
Python 3.3 or higher is required
# 下载蛋白质数据库
https://bioinf.uni-greifswald.de/bioinf/partitioned_odb11/
# 在~/.bashrc加入环境变量
export PROTHINT_PATH=/home/yangkaiyan/biosoft/ProtHint-2.6.0/bin
# 测试（一般装好都能正常运行）
bin/prothint.py genome.fasta proteins.fasta
-----------
Stringtie 2.2.1 （可以用conda安装）
https://ccb.jhu.edu/software/stringtie/#install
# 下载
wget https://github.com/gpertea/stringtie
cd stringtie
# 安装
make release
# 在~/.bashrc中加入环境变量
export PATH=$PATH:/home/yangkaiyan/biosoft/stringtie-2.2.1
----------
Bedtools 2.30.0（可以用conda安装）
# 下载
wget https://github.com/arq5x/bedtools2/releases/download/v2.30.0/bedtools-2.30.0.tar.gz
# 解压
tar -zxvf bedtools-2.30.0.tar.gz
cd bedtools2
# 安装
make
# 在~/.bashrc中加入环境变量
export PATH=$PATH:/home/yangkaiyan/biosoft/bedtools2/bin
----------
Gffread 0.12.7（可以用conda安装）
# 下载
wget https://github.com/gpertea/gffread/releases/download/v0.12.7/gffread-0.12.7.Linux_x86_64.tar.gz
# 解压
tar xzf gffread-0.12.7.Linux_x86_64.tar.gz
cd gffread-0.12.7.Linux_x86_64
# 安装
make 
# 在~/.bashrc中加入环境变量
export PATH=$PATH:/home/yangkaiyan/biosoft/gffread-0.12.7.Linux_x86_64
----------
TSEBRA 1.1.0 
TSEBRA是一个整合软件
# 环境要求
Python 3.5.2 or higher is required
# 下载
wget https://github.com/Gaius-Augustus/TSEBRA
# 解压
tar -zvxf TSEBRA-1.1.0.tar.gz
cd TSEBRA-1.1.0
# 测试
python ./bin/tsebra.py --help
# 在~/.bashrc中加入环境变量
export TSEBRA_PATH=/home/yangkaiyan/biosoft/TSEBRA-1.1.0/bin
----------
可选择安装的软件： 
GUSHR 1.0.0 （预测UTR时需要）
https://github.com/Gaius-Augustus/GUSHR
----------
twoBitInfo和faToTwoBit （预测UTR时需要）
http://hgdownload.soe.ucsc.edu/admin/exe
----------
SRA Toolkit https://github.com/ncbi/sra-tools/wiki（下载SRA转录组数据时需要）
----------
Spaln https://github.com/ogotoh/spaln

1.4 安装Braker3

终于到安装braker3了，截止到目前，Braker的最新版本是Braker3.0.3，该版本修复了Braker3.0.2中gtf和gff文件之间的转换问题

# 下载
wget https://github.com/Gaius-Augustus/BRAKER
# 解压
tar -xzvf BRAKER-3.0.3.tar.gz
cd /biosoft/BRAKER-3.0.3/scripts
# 修改权限
cd scripts
ls -l *.pl *.py
chmod a+x *.pl *.py
# 加入环境变量
export PATH=$PATH:/home/yangkaiyan/biosoft/BRAKER-3.0.3
# 测试
cd example
# 下载RNAseq.bam数据
wget http://topaz.gatech.edu/GeneMark/Braker/RNAseq.bam
# 运行test3.sh
perl /home/yangkaiyan/biosoft/BRAKER-3.0.3/scripts/braker.pl --genome genome.fa --prot_seq proteins.fa --bam RNAseq.bam --threads 8 --skipOptimize & 
# 对比自己的结果文件与官方results文件是否一致

2. 使用container安装

https://hub.docker.com/r/teambraker/braker3

# 创建singularity镜像
singularity build braker3.sif docker://teambraker/braker3:latest
singularity exec braker3.sif print_braker3_setup.py
# 把test脚本复制到当前目录
singularity exec -B $PWD:$PWD braker3.sif cp /opt/BRAKER/example/singularity-tests/test3.sh .
# 创建镜像文件的环境变量 
export BRAKER_SIF=/your/path/to/braker3.sif
# 测试
bash test3.sh

二、使用方法

为了能顺利运行Braker3，我们按照官网的建议做好准备工作

使用高质量基因组
染色体名称不能过长，最好是>contig01，不含%&!*(){})等特殊字符。
使用软屏蔽基因组（重复序列中大写的ATCG改为小写而不是N碱基）
检查物种是否具有进化分支特征
预测结果不要直接使用，先在UCSC Genome Browser中检查

重复序列屏蔽

除了一些转座酶基因外，多数重复序列是没有编码蛋白质的功能的，当我们在一个新的基因组中寻找蛋白质编码基因时，基因组应该做重复序列屏蔽，以此避免预测结果中偏高的假阳性比例。

软件准备：Repeatmasker Repeatmodeler
参考基因组文件：genome.fa
# 构建基因组索引
BuildDatabase -name riparia -engine ncbi genome.fa
# 从头注释 （这里我加上了-LTRStruct管道）
RepeatModeler -database riparia -engine ncbi -LTRStruct -pa 20 >& run.out & 
# 若要进一步注释unknown，可将得到的consensi.fa.classified用seqkit划分为known和unknown
# 去冗余后，用deepTE或PASTEClassifier注释unkonwn，能很大程度上降低unknown比例，这里不再赘述啦
# 构建重复序列库
cat knowns.fasta unknown.fasta repbase_Hexapoda.fa > custom_lib.fa
# Repeatmasker进行软屏蔽（使用-xsmall参数）
RepeatMasker -xsmall -nolow -norna -html -gff -dir ./output -lib custom_lib.fa genome.fa >& run_out &
得到文件：soft_genome.fa

接下来，我们进行转录组reads与基因组的比对

2. 转录组reads比对到基因组（若使用RNA-seq使用fastaq文件，可不用做这一步）

软件准备：Hisat2 
参考基因组文件：soft_genome.fa
# 构建索引 （最好进入到参考基因组的文件夹里去操作，得到8个索引文件）
hisat2-build -p 36 soft_genome.fa idx
# 执行比对 (在braker中涉及到有参转录本组装，所以一定要加--dta参数，还可以使用--max-intronlen --min-intronlen限制内含子大小）
hisat2 --dta -x idx -p 36 -1 ./T01_1.fq.gz -2 ./T01_2.fq.gz |samtools sort -@ 10 > T01.bam & 
最终得到每个样品对应的bam文件，Braker3的stringtie进行有参转录本组装时会用到

3. 蛋白数据库的构建

3.1 OrthoDB蛋白数据库下载

# 我下载的是后生动物的同源蛋白 （压缩文件有4.6G）
wget https://bioinf.uni-greifswald.de/bioinf/partitioned_odb11/Metazoa.fa.gz

同源蛋白的数据库比较大，运行时非常耗费时间，所以我选择了几个近缘物种构建蛋白库。

3.2 自定义蛋白数据库

# 下载近缘物种的蛋白质序列 （以黑腹果蝇为例）
wget https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/215/GCF_000001215.4_Release_6_plus_ISO1_MT/GCF_000001215.4_Release_6_plus_ISO1_MT_protein.faa.gz 
# 合并蛋白质序列
cat Dromel_GCF_000001215.4/Dromel_pr.fa \
Bicany_GCF_947172395.1/Bicany_pr.fa \
Schgre_GCF_023897955.1/Schgre_pr.fa \
Tricas_GCF_000002335.3/Tricas_pr.fa \
Acypis_GCF_005508785.2/Acypis_pr.fa > homo_pr.fa
# 将得到的homo_pr.fa去冗余

可以比较一下官网建议的OrthoDB蛋白数据库和自定义数据库哪一种注释效果更好。

4. Braker3进行基因组结构预测

终于来到了Braker3，我选择的是结合转录组数据以蛋白数据运行BRAKER，最后TSEBRA会整合两者的基因集，也就是下图的流程：

软件准备：本文第一部分软件安装的所有内容
参考基因组文件：soft_genome.fa
转录组数据：T01.bam, T02.bam, T03.bam, T04.bam, T05.bam, T06.bam
蛋白数据库：homo_pr.fa
perl ./scripts/braker.pl --genome /home/yangkaiyan/soft_genome.fa \
--prot_seq /home/yangkaiyan/homo_pr.fa \
--bam T01.bam,T02.bam,T03.bam,T04.bam,T05.bam,T06.bam \
--threads 8 --skipOptimize --gff3 &

4.1 Braker3的三种输入模式

可以直接输入SRA ID，SRA tookit可以直接从网站上把转录组数据下载下来得到fastq文件，再执行比对，也可以直接输入fastq文件或使用比对生成的bam文件。

Braker1、Braker2、Braker3都用的是同一个perl脚本braker.pl，如果输入文件中有蛋白和转录组数据，则调用GeneMark-ETP，如果输入文件中仅有转录组或蛋白数据，则调用GeneMark-ES/ET/EP

4.2 参数介绍

--genome 基因组文件位置（默认为软屏蔽）
--prot_seq 自定义蛋白库的文件位置
--bam 输入的bam文件位置
–-threads 线程数
--skipOptimize 省略AUGUSTUS中的optimize_augustus.pl优化，优化很耗费时间，
但是准确度提升不了多少，优化完成后调用etraining重新训练
--gff3 输出gff3文件
–-workingdir 工作目录位置
其他参数详见帮助文档

4.3 结果文件展示

成功运行后的结果文件如下图

braker.codingseq 基因集的fasta文件
braker.aa 基因集的蛋白序列的fasta文件
braker.gtf Braker预测的基因集最终文件
braker.gff3 Braker预测的基因集的gff3格式文件
Augustus文件夹 AUGUSTUS预测的基因集
GeneMark-ETP文件夹 GeneMark-ETP预测的基因集以及其他中间文件
hintsfile.gff 从RNA-Seq数据和蛋白库数据中提取的外部证据，有内含子等信息

这些结果中是不包含UTR的，需要调用GUSHR且在参数中添加--addUTR=on预测UTR信息

三、总结

Braker与以前相比，流程清晰了许多，根据文献中的实验结果，Braker3的特异性(specificity)在基因水平和外显子水平都有所提升，准确性也比其他流程更高，由此表明预测的基因中假阳性比例偏少。

以上就是Braker3使用方法的简单分享，可能会有疏漏之处，恳请小伙伴们批评指正！博观约取，厚积薄发，探索的路上，勇往直前吧！

参考文献：

https://www.researchgate.net/publication/367409816_The_BRAKER3_Genome_Annotation_Pipeline

玖叶教程网

前端编程开发入门

基因组结构预测|Braker3流程（基因组结构分析）