写在前面

前面两篇文章De novo组装#01 | 测序数据质控（fasqc+fastp） | De novo组装#02 | 基因组survey (jellyfish+GenomeScope2)我们利用二代数据做了基因组的survey，大致了解了该物种的基因组的大小，然后根据预算和目的制定相应的的三代测序策略（测序平台和测序深度），拿到三代测序数据后，用相应的组装软件对三代测序reads经行拼接，我这边用的数据是PacBio Sequel II平台的CLR（Continuous Long Reads）数据，所以下面我分别用flye，wtdbg2，mecat2，canu这四款经典能够组装CLR数据的三代组装软件尝试对我的基因组进行初步组装，pacbio官方的falcon软件因为跟canu一样跑得比较慢，没那么多运算资源，加之测序公司给我的一版基因组就是falcon最新版的pb-assembly跑的我这边就不跑了，后面有空再试下吧

原始测序数据的转换 bam2fastx

测序公司返回的数据是.bam格式的，里面的信息较多，我们需要用bam2fastx将bam文件文件转换成适用于组装fasta或者fastq格式的reads文件，其他软件也行。
1#安装bam2fastx
官网：https://github.com/PacificBiosciences/bam2fastx

conda create -n bam2fastx
conda activate bam2fastx
conda install -c bioconda bam2fastx

2#bam转换成fasta

bam2fastx --fasta -o P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta  P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.bam

一、flye组装基因组

flye 软件支持CLR、HiFi 和nanopore 数据输入。
软件主页：https://github.com/fenderglass/Flye
1#安装flye

## 手动安装
git clone https://github.com/fenderglass/Flye
cd Flye
make

## 自动安装
conda install flye

2#用flye进行初步组装

/newlustre/home/jfgui/Wangtao/software/Flye/bin/flye  \
    --pacbio-raw ../P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta \ ##转换好的测序reads文件
    --out-dir ./  \ ## 结果输出路径
    --genome-size 771346780 \  ## 前面预估的基因组大小
    --threads 24   ## 线程数

3#查看组装结果

flye输出结果

用了24个cpu核心运行，大致跑了1700个小时cpu time实际跑了3天多，整个过程分为5步，具体的中间文件放在了分别放在了00-aasembly，10-consensus，20-repeat，30-cotigger，40-polishing文件夹里。

#安装assembly-stats查看primary assembly

## 手动安装,/foo/bar/为安装路径
git clone https://github.com/sanger-pathogens/assembly-stats.git
mkdir build
cd build
cmake -DINSTALL_DIR:PATH=/foo/bar/ ..
make
make test
make install

## 自动安装
conda install -c bioconda assembly-stats

最后组装出的primary assembly共680个共contig772M，最长的一条为25.18M，平均长度为1.13M，contigN50为9.21M。

$assembly-stats assembly.fasta
stats for assembly.fasta
sum = 771989320, n = 680, ave = 1135278.41, largest = 25186124
N50 = 9211426, n = 27
N60 = 7566198, n = 36
N70 = 6405554, n = 47
N80 = 5050486, n = 61
N90 = 2266136, n = 81
N100 = 493, n = 680
N_count = 0
Gaps = 0

其他相关好文推荐： 一文看懂三代组装软件——Flye - 简书 (jianshu.com)

二、wtdbg2组装基因组

wtdbg2 软件同样支持pacbio和nanopore的数据输入，特点就是快，占用运算资源少，但也组装得较碎。
软件主页：https://github.com/ruanjue/wtdbg2

1#安装wtdbg2

## 手动安装
git clone https://github.com/ruanjue/wtdbg2
cd wtdbg2 && make

## 自动安装
conda install -c bioconda wtdbg

2#用wtdbg2进行初步组装
分两步，第一步是组装，第二步是得到一致性序列。

## wtdbg2组装
/newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/wtdbg2   \
    -t 12  \   ## 线程数
    -x sq  \   ## 设置测序平台，rs/sq/ont
    -g 771M  \  ## 基因组大小
    -i ../P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta -o prefix   ## 输入测序数据

## wtpoa-cns得到一致性序列
/newlustre/home/jfgui/wy/anaconda3/envs/pengfang/bin/wtpoa-cns  \
    -t 12  \   ## 线程数
    -i prefix.ctg.lay.gz  \   ## 前一步输出的ctg.lay.gz文件
    -f  \   ## 覆盖之前结果
    -o wtdbg2_assembly.fa   ## 输出的.fa文件名

3#查看组装结果

wtdbg2输出结果

虽然给了最少的运算资源但还是最快跑完了，12个cpu核心，两步大致跑了750个小时cpu time实际跑了大概2.5天，lay文件跑了2天consensus跑了0.5天。确实是快！！
最后组装出的primary assembly共2697个共contig796M，最长的一条为24.25M，平均长度为0.29M，contigN50为8.4M。虽然快这结果确实组装得比较碎。

$assembly-stats wtdbg2_assembly.fa
stats for wtdbg2_assembly.fa
sum = 796796038, n = 2697, ave = 295437.91, largest = 24257430
N50 = 8405757, n = 33
N60 = 6518067, n = 44
N70 = 4827838, n = 58
N80 = 2584532, n = 80
N90 = 370995, n = 162
N100 = 2250, n = 2697
N_count = 0
Gaps = 0

其他相关好文推荐： wtdbg2：原理及用法 - 简书 (jianshu.com)

三、mecat2组装基因组

mecat2 只适用于pacbio 数据输入的软件，基于canu 改写，加速了其中的比对和纠错模块，其他是一样的。感觉也还是相对慢些
软件主页：https://github.com/xiaochuanle/MECAT2

1#安装mecat2

git clone https://github.com/xiaochuanle/MECAT2.git
cd MECAT2
make 

2#用mecat2进行初步组装
先生成一个配置文件config，根据自身的情况修改填充配置文件；随后分3步，纠错correct，修剪trim，组装assemble运行。

## 生成配置文件
/newlustre/home/jfgui/Wangtao/software/MECAT2/Linux-amd64/bin/mecat.pl config mecat2.cfg

## 修改填充.cfg文件
$vi mecat2.cfg
  1 PROJECT=sl_mecat2   # 工程名/输出目录
  2 RAWREADS=../P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta   # 输出测序数据
  3 GENOME_SIZE=771346780   # 基因组大小
  4 THREADS=16   # 线程数
  5 MIN_READ_LENGTH=2000
  6 CNS_OVLP_OPTIONS="-kmer_size 13"
  7 CNS_PCAN_OPTIONS="-p 100000 -k 100"
  8 CNS_OPTIONS=""
  9 CNS_OUTPUT_COVERAGE=30
 10 TRIM_OVLP_OPTIONS="-skip_overhang"
 11 TRIM_PM4_OPTIONS="-p 100000 -k 100"
 12 TRIM_LCR_OPTIONS=""
 13 TRIM_SR_OPTIONS=""
 14 ASM_OVLP_OPTIONS=""
 15 FSA_OL_FILTER_OPTIONS="--max_overhang=-1 --min_identity=-1"
 16 FSA_ASSEMBLE_OPTIONS=""
 17 CLEANUP=0


## reads纠错-correct
/newlustre/home/jfgui/Wangtao/software/MECAT2/Linux-amd64/bin/mecat.pl  correct  mecat2.cfg

## reads修剪-trim
/newlustre/home/jfgui/Wangtao/software/MECAT2/Linux-amd64/bin/mecat.pl  trim  mecat2.cfg

## 组装-assemble
/newlustre/home/jfgui/Wangtao/software/MECAT2/Linux-amd64/bin/mecat.pl  assemble  mecat2.cfg

3#查看组装结果

mecat2输出结果

最终结果4-fsa

16个cpu运行了大概2000个小时的cpu time实际大概运行了5天多，比flye要稍微慢点，还行。

最终输出5个目录，我们需要的primary assembly就在4-fsa里的contigs.fasta。用assembly-stats查看下，总共组装出1210个contig共782.4Mb，平均长度为0.65M，最大的一条contig为20.04M，contigN50为6.79M，单纯从contigN50这个参数上这个软件是不如前两款软件的。

$assembly-stats contigs.fasta
stats for contigs.fasta
sum = 782440466, n = 1210, ave = 646645.01, largest = 20039767
N50 = 6788542, n = 36
N60 = 5033148, n = 50
N70 = 3470086, n = 68
N80 = 1490280, n = 101
N90 = 380864, n = 200
N100 = 641, n = 1210
N_count = 0
Gaps = 0

其他相关好文推荐：
用mecat2组装基因组 - 简书 (jianshu.com)
MECAT2: 三代高效组装工具 - 知乎 (zhihu.com)

四、canu组装基因组

canu 同时适配CLR，HIFI和nanopore 数据输入，大家都说这软件组装的挺好的，但就是慢非常慢，非常消耗运算资源，所以也这边也尝试下，组装结果后续再更新把太久了😂
软件主页：https://github.com/marbl/canu

1#安装canu

## 手动安装（官方推荐二进制版）
wget https://github.com/marbl/canu/releases/download/v2.2/canu-2.2.Linux-amd64.tar.xz   # 获取二进制文件
tar -xJf canu-2.2.*.tar.xz  # 提取文件

## 自动安装
conda install -c bioconda -c defaults canu

2#用canu进行初步组装
和上面的mecat2软件一样，同样内置3个步骤correct校正 reads，然后trim修剪 reads，然后将 reads assemble组装成 contigs，但我们不要分步设置，软件默认依次执行这三个步骤。由于第一步校正所用的时间太久了，所以默认情况下canu会根据基因组大小选取40X的最长的reads进行校正然后继续进行后续的步骤。我们测序深度很高的话这软件会自动舍去一部分数据（短reads）,我们如果要用全部的数据去做correct则可以设置corOutCoverage=999，当然代价是运行时间。

## 直接一步运行
/newlustre/home/jfgui/Wangtao/software/canu-2.2/bin/canu  \
     -pacbio ../P01TYD20308306-1_r64030_20201110_065049_1_A02.subreads.fasta  \  # 输出测序数据
     -p slcanu  \ # 输出前缀
     -d ./sl_canu  \ # 输出目录
     genomeSize=771.3m  \ # 基因组大小
     useGrid=false   \ # 是否使用集群，默认是使用的，但我实际测试时不设置false的话会没法运行不知道为什么，所以后面是单节点上跑的？
     maxThreads=40  # 最大线程数

## 如果需要全部数据做correct需要加上：
     corOutCoverage=999

3#查看组装结果

canu输出结果

sl_canu目录里，我们用assembly-stats查看初始组装，总共组装出2535个contig共927.2Mb，平均长度为0.36M，最大的一条contig为18.2M左右，contigN50为5.36M，感觉这软件组装得最散，而且应该是重复序列比例太高导致组装结果偏大，后续再看看能否去冗余purge。

$assembly-stats slcanu.contigs.fasta
stats for slcanu.contigs.fasta
sum = 927204322, n = 2535, ave = 365761.07, largest = 18199312
N50 = 5368167, n = 49
N60 = 3542224, n = 70
N70 = 1513046, n = 109
N80 = 321576, n = 235
N90 = 83471, n = 949
N100 = 1478, n = 2535
N_count = 0
Gaps = 0

其他相关好文推荐：
序列组装软件—Canu - 知乎 (zhihu.com)
三代组装软件Canu简介 - 简书 (jianshu.com)