动物随随便便几个G，植物动不动十多个G，而细菌基因组一般就几M，基因组结构简单，没有重复序列杂合度这些组装障碍。我现在很多文章细菌基因组开始用三代测序组装了，还有些文章是二代结合三代测序？？？纳税人的钱不是钱的感觉，这篇文章是测试原核生物的三代数据组装软件及评测。测试软件及数据类型持续更新~
三代reads一般长度是几K到几十K，组装植物动物基因组，三代数据的优势是很明显的。二代短测序平台存在较三代更明显的GC偏好性，以及重复序列的问题。
格式转化软件
#下载两个测试数据 hdf5 文件
wget -c -r -nd -np -k -L -p https://sra-download.ncbi.nlm.nih.gov/traces/era21/ERZ/002695/ERR2695069/ERR2695069_hdf5.tgz
wget -c -r -nd -np -k -L -p https://sra-download.ncbi.nlm.nih.gov/traces/era21/ERZ/002632/ERR2632363/ERR2632363_hdf5.tgz
# md5sum check
软件安装：
## 格式转化软件
pbh5tools
## 组装软件
# canu
数据预处理：
### hdf5文件解压
$ tar -zxf *tgz
#解压后会生成多个文件
m160607_042209_00127_c100967912550000001823217106101682_s1_p0.metadata.xml
m160607_042209_00127_c100967912550000001823217106101682_s1_p0.bas.h5
m160607_042209_00127_c100967912550000001823217106101682_s1_p0.1.bax.h5
m160607_042209_00127_c100967912550000001823217106101682_s1_p0.3.bax.h5
m160607_042209_00127_c100967912550000001823217106101682_s1_p0.2.bax.h5
### 转换格式
# ”--readType“ 可选ccs，subreads，unrolled。不了解可以看测序原理
$ python /XXX/pbh5tools/bin/bash5tools.py m160607_042209_00127_c100967912550000001823217106101682_s1_p0.bas.h5 --outFilePrefix ERR2632363.subreads --readType subreads --outType fasta
$ python /XXX/pbh5tools/bin/bash5tools.py m160706_123219_00127_c100991152550000001823221807191673_s1_p0.bas.h5 --outFilePrefix ERR2695069.subreads --readType subreads --outType fasta
hdf5解压文件：
m160607_042209_00127_c100967912550000001823217106101682_s1_p0.metadata.xml
m160607_042209_00127_c100967912550000001823217106101682_s1_p0.bas.h5
m160607_042209_00127_c100967912550000001823217106101682_s1_p0.1.bax.h5
m160607_042209_00127_c100967912550000001823217106101682_s1_p0.3.bax.h5
m160607_042209_00127_c100967912550000001823217106101682_s1_p0.2.bax.h5
待讲解。。。。。。。。
我简单统计了一下两个fasta文件, 基因组本身大小约1.8M，两个测试数据平均长度就三四K，深度分别是257X, 743X。。。。，既然数据这么多，我简单过滤一下,过滤掉小于2K的reads。
软件测评 ——基因组组装
1.canu
用过这款软件装过基因组的都知道。这款软件留给大家最深刻的印象就是：慢。由于数据量小，这一劣势在细菌基因组上没有任何表现。
组装分为三步correct，trim,assemble
###偷懒批处理
$ cat sample.list | while read a;do echo "canu -correct -p ${a} -d ${a}_cor genomeSize=1.8m -pacbio-raw ${a}.subreads.2k.fasta 1>canu.${a}.cor.log 2>canu.${a}.cor.err && \
canu -trim -p ${a} -d ${a}_trim genomeSize=1.8m -pacbio-corrected ${a}_cor/${a}.correctedReads.fasta.gz 1>canu.${a}.trim.log 2>canu.${a}.trim.err && \
canu -assemble -p ${a} -d ${a}_0.045 correctedErrorRate=0.045 genomeSize=1.8m -pacbio-corrected ${a}_trim/${a}.trimmedReads.fasta.gz 1>canu.${a}.assemble.log 2>canu.${a}.assemble.err";done >run.canu.sh
# 看一行运行脚本
$ cat run.canu.sh | head -n 1
canu -correct -p ERR2632363 -d ERR2632363_cor genomeSize=1.8m -pacbio-raw ERR2632363.subreads.2k.fasta 1>canu.ERR2632363.cor.log 2>canu.ERR2632363.cor.err && \
canu -trim -p ERR2632363 -d ERR2632363_trim genomeSize=1.8m -pacbio-corrected ERR2632363_cor/ERR2632363.correctedReads.fasta.gz 1>canu.ERR2632363.trim.log 2>canu.ERR2632363.trim.err && \
canu -assemble -p ERR2632363 -d ERR2632363_0.045 correctedErrorRate=0.045 genomeSize=1.8m -pacbio-corrected ERR2632363_trim/ERR2632363.trimmedReads.fasta.gz 1>canu.ERR2632363.assemble.log 2>canu.ERR2632363.assemble.err
### cor，trim大大减少了数据量
file format type num_seqs sum_len min_len avg_len max_len
ERR2632363_cor/ERR2632363.correctedReads.fasta.gz FASTA DNA 6,104 69,719,744 205 11,422 29,482
ERR2695069_cor/ERR2695069.correctedReads.fasta.gz FASTA DNA 41,000 146,577,915 3 3,575.1 27,783
ERR2632363_trim/ERR2632363.trimmedReads.fasta.gz FASTA DNA 6,085 63,644,442 1,005 10,459.2 27,739
ERR2695069_trim/ERR2695069.trimmedReads.fasta.gz FASTA DNA 39,487 142,538,939 1,000 3,609.8 27,740
组装结果：
file format type num_seqs sum_len min_len avg_len max_len
ERR2632363_0.045/ERR2632363.contigs.fasta FASTA DNA 1 1,760,152 1,760,152 1,760,152 1,760,152
ERR2695069_0.045/ERR2695069.contigs.fasta FASTA DNA 8 1,888,736 2,055 236,092 1,872,113
错误率：
The default is 0.045 for PacBio reads, and 0.144 for Nanopore reads.
For low coverage datasets (less than 30X), we recommend increasing [correctedErrorRate] slightly, by 1% or so.
For high-coverage datasets (more than 60X), we recommend decreasing [correctedErrorRate] slightly, by 1% or so.
ERR2695069组装结果有点不理想，在79X的情况下，尝试多个correctedErrorRate值，可得到更好组装结果
运行时间：
没有统计运行时间，在我小服务器满载的情况下，两个组装一起跑，三步大概分别耗时：约1h，几分钟，几分钟。。。
网友评论