原核生物基因组三代数据(pacbio/nanopore)组装

动物随随便便几个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,几分钟,几分钟。。。

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 193,812评论 5 459
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 81,626评论 2 371
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 141,144评论 0 319
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 52,052评论 1 263
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 60,925评论 4 355
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 46,035评论 1 272
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 36,461评论 3 381
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 35,150评论 0 253
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 39,413评论 1 290
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 34,501评论 2 307
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 36,277评论 1 325
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 32,159评论 3 312
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 37,528评论 3 298
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 28,868评论 0 17
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 30,143评论 1 250
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 41,407评论 2 341
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 40,615评论 2 335

推荐阅读更多精彩内容