MECAT2: 三代高效组装工具
MECAT2是三代测序数据PacBio的高效组装工具，是之前MECAT的改进版, 修复了之前的很多bug, 使用基于string graph的fsa替换了之前的mecat2canu。
MECAT2由4个模块组成:
-
mecat2pw: SMART reads快速和准确地配对比对工具 -
mecat2ref: SMART reads的参考基因组比对工具 -
mecat2cns: 基于配对的重叠对存在噪音的reads进行纠错 -
fsa: 基于string graph的组装工具
目前MECAT2只支持PacBio数据，对于Nanopore数据，肖老师他们开发了NECAT
安装
MECAT2软件安装非常简单，不存在上一代MECAT的安装难问题了
mkdir -p ~/opt/biosoft
cd ~/opt/biosoft
git clone https://github.com/xiaochuanle/MECAT2.git
cd MECAT2
make -j 4
export PATH=~/opt/biosoft/MECAT2/Linux-amd64/bin:$PATH
唯一要注意的是尽量要用新版的perl，经测试，Perl 5.26 是可以正常运行脚本。
最后加入到环境变量PATH中
export PATH=~/opt/biosoft/MECAT2/Linux-amd64/bin:$PATH
软件使用
我们以 E. coli 数据集作为案例，介绍MECAT2应该如何使用。
第一步，我们要下载数据
mkdir -p ~/ecoli
cd ~/ecoli
wget http://gembox.cbcb.umd.edu/mhap/raw/ecoli_filtered.fastq.gz
gunzip ecoli_filtered.fastq.gz
注意: 目前MECAT2还不支持gz压缩文件，如果你直接用gz作为输入，会提示core dumped。
MECAT2不是根据文件名来确定输入的序列格式, 也就是如果你把 ecoli_filtered.fastq 命名为 ecoli_filtered.fasta, 组装也不会出错。
第二步, 准备模板的config文件
mecat.pl config ecoli_config_file.txt
使用vim修改config文件
PROJECT=ecoli
RAWREADS=ecoli_filtered.fastq
GENOME_SIZE=4800000
THREADS=4
MIN_READ_LENGTH=2000
CNS_OVLP_OPTIONS=""
CNS_OPTIONS="-r 0.6 -a 1000 -c 4 -l 2000"
TRIM_OVLP_OPTIONS="-B"
ASM_OVLP_OPTIONS="-n 100 -z 10 -b 2000 -e 0.5 -j 1 -u 0 -a 400"
FSA_OL_FILTER_OPTIONS="--max_overhang=-1 --min_identity=-1"
FSA_ASSEMBLE_OPTIONS=""
USE_GRID=false
CLEANUP=0
CNS_OUTPUT_COVERAGE=30
GRID_NODE=0
第三步: 原始数据纠错
mecat.pl correct ecoli_config_file.txt
第四步: 对纠错后的reads进行组装
mecat.pl assemble ecoli_config_file.txt
输出结果:
- 纠错后的reads：
ecoli/1-consensus/cns_reads.fasta. - 最长30X纠错后用于trimming的reads:
ecoli/1-consensus/cns_final.fasta. - trimmed reads:
ecoli/2-trim_bases/trimReads.fasta - 组装的contigs:
ecoli/4-fsa/contigs.fasta
最终的组装结果为
Count: 1
Tatal: 4630437
Max: 4630437
Min: 4630437
N25: 4630437
L25: 1
能够完美的就装出一条序列
配置文件介绍
MECAT2的组装过程的参数都在配置文件中, 也就是之前的ecoli_config_file.txt。我按照不同部分对这些参数进行介绍
基本参数
这些参数属于改起来不怎么纠结的参数
-
PROJECT=ecoli: 项目名，之后的所有输出文件都在该项目名的目录下 -
RAWREADS=: 原始数据的位置，可以是Fasta或者是Fastq， H5格式要先转成Fasta (参考 Input Format). -
GENOME_SIZE=: 基因组大小，单位是bp，如果是120M，那么就是120000000. -
THREADS=: 线程数 -
MIN_READ_LENGTH=: 用于纠错和trim的reads的最低长度. 数据质量好，就长一点，比如说1000到2000 -
USE_GRID=false: 是否有多个计算节点 -
CLEANUP=0: 运行结束后删除MECAT2的中间文件, 大基因组的临时文件很大，所以要设置为1. -
GRID_NODE=0, 当USE_GRID=1时，设置用到的计算节点数，如果是单节点服务器，不需要设置。
纠错和修剪阶段
纠错阶段(correct stage)和修剪阶段(trimming stage)，MECAT2调用的mecat2pw和mecat2cns， 与之相关的配置如下
-
CNS_OVLP_OPTIONS="", 在纠错阶段是检测候选重叠的参数, 会传给mecat2pw -
CNS_OPTIONS=""：原始reads纠错参数，会传递给mecat2cns, -
TRIM_OVLP_OPTIONS="": 在trim阶段，用于检测重叠的参数,会传给v2asmpm
但是我发现v2asmpm没有文档说明，和软件开发者讨论之后，给出的说明如下
V2asmpm -Pworkpath -Tt -Sx -Ey -B
-Pworkpath 工作目录是workpath
-Tt 用t个CPU线程
-Sx -Ey 这两个参数一起表示只计算第x到第y卷(包括第y卷), 通常工作目录下会有000001.fasta, 000002.fasta这样的分卷
-B 如果有这个选项, 就输出二进制格式的比对结果; 如果没有这个选项, 就输出文本格式的比对结果
因此用默认的-B就行了。
组装阶段
组装阶段相关的参数如下
-
CNS_OUTPUT_COVERAGE=30: 选择多少覆盖度的最长纠错后reads进行trim和组装。举个例子，4.8M的 E. coli 30X, 是144MB。 一般选择20~40之间就行， 会传递给mecat2elr -
ASM_OVLP_OPTIONS="": 在组装阶段，用于检测重叠的参数，传给v2asmpm.sh -
FSA_OL_FILTER_OPTIONS="", 过滤重叠的参数，传递给fsa_ol_filter -
FSA_ASSEMBLE_OPTIONS="", 组装trimm reads的参数, 传给v2asmpm
这些参数调整起来都很复杂，只能建议看函数帮助文档了。举个例子FSA_OL_FILTER_OPTIONS="--max_overhang=-1 --min_identity=-1"根据阅读fsa_ol_filter的帮助发现，-1表示让程序自己决定
而ASM_OVLP_OPTIONS="-n 100 -z 10 -b 2000 -e 0.5 -j 1 -u 0 -a 400"， 根据和软件开发者的讨论，是不需要额外调整，用默认的即可。
不客观的个人体验测评
MECAT2是基于比对纠错的组装工具，组装速度上比同类型的Canu和Falcon都快。仅看组装结果的N50参数，MECAT2结果和Canu, Falcon结果都是差不多的。
考虑到BUSCO完整度，Canu是三者最高。不过用三代数据和二代数据进行几轮Polish后会提高到Canu相同水平。三个软件上的BUSCO（embryophyta_odb10）值如下:
Canu: C:98.8%[S:96.3%,D:2.5%],F:0.5%,M:0.7%,n:1375
MECAT2: C:93.5%[S:91.6%,D:1.9%],F:4.6%,M:1.9%,n:1375
Falcon: C:96.8%[S:95.1%,D:1.7%],F:2.0%,M:1.2%,n:1375
因此，如果你原本你是用Falcon做组装，那么可以改成MECAT2了，效果差不多，速度还更快了。对于越大的物种，在有限的资源下，更推荐用MECAT2。
