BUSCO - 组装质量评估

2018-08-28 本文已影响3人白菜代码小推车

busco.jpeg

BUSCO组装质量评估软件

BUSCO - Benchmarking Universal Single-Copy Orthologs

文献

《BUSCO: assessing genome assembly and annotation completeness with single-copy orthologs》

说明

BUSCO是一款使用python语言编写的对转录组和基因组组装质量进行评估的软件。在相近的物种之间总有一些保守的序列，而BUSCO就是使用这些保守序列与组装的结果进行比对，鉴定组装的结果是否包含这些序列，包含单条、多条还是部分或者没有等等。

BUSCO 软件根据OrthoDB 数据库，构建了几个大的进化分支的单拷贝基因集。将转录本拼接结果与该基因集进行比较，根据比对上的比例、完整性，来评价拼接结果的准确性和完整性。

它的流程是

genoem assemble : tBLASTn --> Augustus --> HMMER3
Transcriptome   :             Find ORF --> HMMER3
Gene set        :                          HMMER3

下载与安装【简单】

conda install busco

使用conda自动完整安装过程，但是有时候没有管理员权限或者其他一些原因无法这样安装，那么就可以根据下面的方法进行安装。

下载与安装【麻烦】

下载

下载软件包

#  ============ 下载BUSCO  ============ 
cd ~/Applications/download
wget -c https://gitlab.com/ezlab/busco/-/archive/master/busco-master.zip -O busco.zip

# ============ 下载依赖的工具 ============
# 下载Augustus
wget -c http://bioinf.uni-greifswald.de/augustus/binaries/augustus.current.tar.gz
# 下载HMMER
wget -c http://eddylab.org/software/hmmer/hmmer.tar.gz -O hmmer.tar.gz
# 下载Blast+
wget -c ftp://ftp.ncbi.nlm.nih.gov/blast/executables/LATEST/ncbi-blast-2.7.1+-x64-linux.tar.gz

下载数据库文件

# 浏览器打开下列网址
https://busco.ezlab.org/
# 查找下载需要评估的组装的物种对应的类群
- Bacteria
- Protists
- Metazoa
- Fungi
- Plant

# 例如下载植物的BUSCO的数据库
mkdir -p ~/database/BUSCO/
cd ~/database/BUSCO/
# 下载
wget -c https://busco.ezlab.org/datasets/embryophyta_odb9.tar.gz
# 解压文件
tar -xzvf embryophyta_odb9.tar.gz

安装

cd ~/stq/Applications

# 安装busco
unzip busco.zip
# 改名
mv busco-master busco 
mv busco ../
cd ../busco
python setup.py install

# 安装Augustus
tar -xzvf augustus.current.tar.gz
cd augustus-3.3.1
# 打开common.mk文件，将ZIPINPUT = true注释掉（即在最前面加一#号）
vim common.mk
make

# 安装HMMER
tar -xzvf hmmer.tar.gz
./configure
make

# 安装blast+
tar -xzvf ncbi-blast-2.7.1+-x64-linux.tar.gz
rm ncbi-blast-2.7.1+-x64-linux.tar.gz
mv ncbi-blast-2.7.1+ blast+-2.7.1-linux

配置

这里配置是必须的，在安装好软件之后 ~/Applications/busco/config/之中·并没有config.ini文件，只有一个config.ini_default文件，可以把里面的内容复制下来，然后把：

将out_path = /workdir前面加上 #
因为这个工具的输出路径有时候会出错，所以干脆将它注释掉，之后假如运行busco之后，输出的路径就是你之前cd到的路径
修改下面的工具的路径
[tblastn]
[makeblastdb]
[augustus]
[etraining]
[gff2gbSmallDNA.pl]
[new_species.pl]
[optimize_augustus.pl]
[hmmsearch]

也可以按照下面的说明进行

# 修改配置文件
vim ~/Applications/busco/config/config.ini

加入如下内容，里面的路径需要更改为你自己的工具的路径

# BUSCO specific configuration
# It overrides default values in code and dataset cfg, and is overridden by arguments in command line
# Uncomment lines when appropriate
[busco]
# Input file
;in = ./sample_data/target.fa
# Run name, used in output files and folder
;out = SAMPLE
# Where to store the output directory
# out_path = /workdir
# Path to the BUSCO dataset
;lineage_path = ./sample_data/example
# Which mode to run (genome / protein / transcriptome)
;mode = genome
# How many threads to use for multithreaded steps
;cpu = 1
# Domain for augustus retraining, eukaryota or prokaryota
# Do not change this unless you know exactly why !!!
;domain = eukaryota
# Force rewrite if files already exist (True/False)
;force = False
# Restart mode (True/False)
;restart = False
# Blast e-value
;evalue = 1e-3
# Species to use with augustus, for old datasets only
;species = fly
# Augustus extra parameters
# Use single quotes, like this: '--param1=1 --param2=2'
;augustus_parameters = ''
# Tmp folder
;tmp_path = ./tmp/
# How many candidate regions (contigs, scaffolds) to consider for each BUSCO
;limit = 3
# Augustus long mode for retraining (True/False)
;long = False
# Quiet mode (True/False)
;quiet = False
# Debug logs (True/False), it needs Quiet to be False
debug = True
# tar gzip output files (True/False)
;gzip = False
# Force single core for the tblastn step
;blast_single_core = True

[tblastn]
# path to tblastn
path = ~/Applications/blast+-2.7.1-linux/bin

[makeblastdb]
# path to makeblastdb
path = ~/Applications/blast+-2.7.1-linux/bin

[augustus]
# path to augustus
path = ~/Applications/augustus-3.3.1/bin

[etraining]
# path to augustus etraining
path = ~/Applications/augustus-3.3.1/bin

# path to augustus perl scripts, redeclare it for each new script        
[gff2gbSmallDNA.pl]                                                      
path = ~/Applications/augustus-3.3.1/scripts                         
[new_species.pl]                                                         
path = ~/Applications/augustus-3.3.1/scripts                         
[optimize_augustus.pl]                                                   
path = ~/Applications/augustus-3.3.1/scripts                         
                                                                         
[hmmsearch]                                                              
# path to HMMsearch executable                                           
path = ~/Applications/hmmer-3.2.1/src                                
                                                                         
[Rscript]                                                                
# path to Rscript, if you wish to use the plot tool                      
path = /usr/bin/

运行

导入环境变量

这里没有对软件进行完整的安装，所以需要导入临时环境变量
下面的路径都是刚才安装的程序的对应的路径

# augustus工具的执行文件所在文件夹
export PATH="~/Applications/augustus-3.3.1/bin:$PATH"
# augustus工具配置文件的所在位置 。AUGUSTUS_CONFIG_PATH 需要使用绝对路径
export AUGUSTUS_CONFIG_PATH="/share/home/ssd/Applications/augustus-3.3.1/config"
# hmmer工具的执行文件所在文件夹
export PATH="~/Applications/hmmer-3.2.1/src:$PATH"
# blast+工具的执行文件所在文件夹
export PATH="~/Applications/blast+-2.7.1-linux/bin:$PATH"

如果经常使用，建议加入永久的环境变量

# 打开.bashrc文件
vim ~/.bashrc

# 在末尾添加上面的导入环境变量的内容

开始评估

# 首先cd到对应的组装文件的文件夹
# 

# -i 输入文件
# -l BUSCO的数据库文件
# -o 输出的文件以及文件夹的后缀
# -m 分析类型（genome、transcriptome、proteins）
# --cpu 线程数
~/Applications/busco/scripts/run_BUSCO.py \
    -i contigs.fasta \
    -l ~/database/BUSCO/embryophyta_odb9 \
    -o suffix\
    -m genome \
    --cpu 8

结果解读

在文件夹下会有

run_suffix 文件夹：因为上面-o选项设置了suffix，所以文件夹名称加上了后缀。在这个文件夹里面，有一个文件最为重要。就是short_summary_suffix.txt,

下面是个范例

# Summarized benchmarking in BUSCO notation for file assembly/spades/contigs.fasta
# BUSCO was run in mode: genome
    C:98.6%[S:98.6%,D:0.0%],F:0.0%,M:1.4%,n:148

    146 Complete BUSCOs (C)
    146 Complete and single-copy BUSCOs (S)
    0   Complete and duplicated BUSCOs (D)
    0   Fragmented BUSCOs (F)
    2   Missing BUSCOs (M)
    148 Total BUSCO groups searched

C	S	D	F	M	Total
Complete	single-copy	duplicated	Fragment	Miss	Total

C : 多少个BUSCO测试基因被覆盖。
S : 多少个基因经过比对发现是单拷贝。
D : 多少个基因经过比对发现是多拷贝。
F : 多少个基因经过比对覆盖不完全，只是部分比对上。
M : 没有得到比对结果的基因数
Total : 总共测试的基因条目数

一般来看S + D的数值也就是C的·值越大越好，但是在文献中作者说如果D的数值太多的话可能意味着组装错误的可能性较大。因为一个基因被覆盖多次，那么可能就是说改基因所在的片段组装可能出现问题。
因为理论上

理论上：--------------------------------------------
实际上：--------   ----------    ----------    -----
错误的：--------                 -------------
           \\\\                    ///
          -----------     --------------

理论上组装之后各个片段之间应该前后有序，之间除了重复区域或者其他特殊片段之外不应该有可以重叠的地方。
例如

sequence1 : ..TAGTCGTGA                   GTGCATGCTGTAGC..
                       \                 /
                        AAAATTGGCGATGAAAA
                       /                 \
sequence2 : ..GGGTAGCGG                   TTGACTAGCTAGCT..

也就是说中间一段序列是两个序列的共同部分，除非这个序列存在多个拷贝，否则就很可能是拼接错误。通常一般这种拼接错误的序列的两端会出现重复序列。

那么也就是说一般来看，按照我个人的看法。S似乎越大越好，M越小越好，说明组装的越完整。但是D与F这两个数值越大不见得就是好的，因为组装饿错误可能会带来这两个值的增大。究竟应该如何评判，我觉得不能仅仅只是通过这一个软件来判定。比如还可以借助QUAST和常规指标N50、总的核酸量、对角线图等等多个评判标准来进行。

参考

我接触生信的时间不长，上面的部分内容可能有错误，望能见谅，还望能指出，多谢。