介绍

上一节基因课的视频中,我们介绍了使用N50等指标评价转录本的拼接结果。

但是这些指标并不能直接反应拼接结果的质量,比如真实转录本中有多少比例被重构出来?是完整的重构出了整条还是只有部分序列?

要进行这方面的评价,可以使用 BUSCO软件。

BUSCO 的原理很简单,在近缘物种中,总有一些基因是保守的。那么,我们就来看这些保守的基因拼接的完整性和准确性如何。

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

BUSCO 还可以用于评价基因组的拼接、基因集。
基因集可以直接使用 HMMER3 与参考数据库比对;
转录本鉴定出最长读码框架之后,再使用HMMER3 比对;
而基因组需要先使用 Augustus 软件进行基因结构预测,在使用 HMMER3 比对。

安装

BUSCO 以来 NCBI-Blast+、HMMER 和 Augustus(评价基因组使用),手动安装比较麻烦。

本着能自动化安装一定不手动操作的原则,使用 Bioconda 来安装。Bioconda 的使用参考课程《生物信息软件安装》

安装 BUSCO

conda install busco

下载数据库

BUSCO 目前支持的数据库列表如下,应选择包含当前物种的数据库进行下载。选择不同的数据库,评价结果可能完全不同。

根据当前物种,选择下载合适的参考数据库

点击 “Download all datasets”
选择灵长类的数据库下载链接

mkdir -p ~/database/busco
cd ~/database/busco
wget http://busco.ezlab.org/v2/datasets/euarchontoglires_odb9.tar.gz
tar -zxvf eukaryota_odb9.tar.gz

使用

需要特别注意的是,输入文件是转录本序列是,每条基因下仅选择一条转录本,否则会被归位‘Duplicate’

mkdir BUSCO
cd BUSCO

#提取每个基因的最长转录本
perl $TRINITY_HOME/util/misc/get_longest_isoform_seq_per_trinity_gene.pl ../Assembly/trinity_out_dir/Trinity.fasta >longest_isoform.fasta

#运行 BUSCO
busco -i longest_isoform.fasta \
        -l /home/zxd/database/busco/euarchontoglires_odb9 \
        -o busco \
        -m tran \
        --cpu 2

参数说明:

-i|--in: 输入文件
-o|--out: 输出文件、文件夹前缀
-l|--lineage: 数据库的路径
-m|--mode: 运行模式,geno|tran|prot
-c|--cpu: 线程数,默认1.

结果解读

BUSCO 的结果一看就懂,主要看统计文件run_busco/short_summary_busco.txt 。我们的测试数据数据量太少,下面放一个真实项目的数据。

$ cat run_busco/short_summary_busco.txt 
# 省略若干行
# Summarized benchmarking in BUSCO notation for file longest_isoform.fasta
# BUSCO was run in mode: tran

    C:89.0%[S:49.9%,D:39.2%],F:4.4%,M:6.5%,n:429

    382    Complete BUSCOs (C)
    214    Complete and single-copy BUSCOs (S)
    168    Complete and duplicated BUSCOs (D)
    19    Fragmented BUSCOs (F)
    28    Missing BUSCOs (M)
    429    Total BUSCO groups searched

真实项目中,Complete BUSCOs (C) 的比例通畅都能达到 80% 以上。不过如果低于这个值,也未必有问题,还是要根据实际项目情况判断。除了跟经验值比较。更有意义的是使用不同软件、参数多组装几个版本,挑选最优版。

如果有多个版本,就可以画这样一张图,分析哪个版本更优。

mkdir my_summaries 
#拷贝多个版本的结果到同一个目录中,注意重命名 
cp run_SPEC1/short_summary_SPEC1.txt my_summaries/ 
cp run_SPEC2/short_summary_SPEC2.txt my_summaries/ 
cp run_SPEC3/short_summary_SPEC3.txt my_summaries/ 
cp run_SPEC4/short_summary_SPEC4.txt my_summaries/ 
cp run_SPEC5/short_summary_SPEC5.txt my_summaries/ 
python BUSCO_plot.py –wd my_summaries