用hifiasm组装基因组

hifiasm大概是目前为止支持PacBio HiFi数据组装的所有软件中表现最优异的软件了。它不但能输出primary assembly result,还能分别组装出单倍体(haplotype)基因组的组装结果。

今天的分享就从hifiasm的组装开始,其他的组装软件之后也会介绍到。

可用资源

软件安装

hifiasm的安装也非常简单,conda一条命令就可以了

conda install hifiasm

运行hifiasm

纯PacBio HiFi组装

把命令写到一个名为run_hifiasm.sh的脚本里

vim run_hifiasm.sh

只要填入最基本的

hifiasm \
-o test.hifisam \
-t 24 \
/path/to/PacBio/test1.hifi_reads.fastq.gz \
/path/to/PacBio/test2.hifi_reads.fastq.gz \
/path/to/PacBio/test3.hifi_reads.fastq.gz

用nohup挂在后台运行就可以啦。

nohup bash run_hifiasm.sh &

输出结果如下:

test.hifisam.ovlp.source.bin
test.hifisam.ovlp.reverse.bin
test.hifisam.ec.bin
test.hifisam.bp.r_utg.noseq.gfa
test.hifisam.bp.r_utg.lowQ.bed
test.hifisam.bp.r_utg.gfa 
test.hifisam.bp.p_utg.noseq.gfa
test.hifisam.bp.p_utg.lowQ.bed
test.hifisam.bp.p_utg.gfa
test.hifisam.bp.p_ctg.noseq.gfa
test.hifisam.bp.p_ctg.lowQ.bed
★  test.hifisam.bp.p_ctg.gfa
test.hifisam.bp.hap2.p_ctg.noseq.gfa
test.hifisam.bp.hap2.p_ctg.lowQ.bed
★  test.hifisam.bp.hap2.p_ctg.gfa
test.hifisam.bp.hap1.p_ctg.noseq.gfa
test.hifisam.bp.hap1.p_ctg.lowQ.bed
★  test.hifisam.bp.hap1.p_ctg.gfa

上面打了五角星的三个文件是后续重点关注的文件。根据
https://hifiasm.readthedocs.io/en/latest/interpreting-output.html的描述:

Hifiasm generates the following assembly graphs only with HiFi reads in default:
prefix.bp.p_ctg.gfa: assembly graph of primary contigs.
prefix.bp.hap1.p_ctg.gfa: partially phased contig graph of haplotype1.
prefix.bp.hap2.p_ctg.gfa: partially phased contig graph of haplotype2.

HiC data integrated Integrated

如果你有HiC数据(测基因组的话大概率会有HiC数据的吧?)可以尝试基于HiC的证据的分型,primary的结果不会有什么变化,但是haplotype的结果会有提高。

hifiasm \
-o test.hifisam_hic \
-t 24 \
--h1 /path/to/HiC/1_R1.fastq.gz --h2 /path/to/HiC/1_R2.fastq.gz \
/path/to/PacBio/test1.hifi_reads.fastq.gz \
/path/to/PacBio/test2.hifi_reads.fastq.gz \
/path/to/PacBio/test3.hifi_reads.fastq.gz

格式转换

用gfatools或者自己写个脚本把gfa格式的文件转换成fasta就得到最终的组装结果啦。

如果写脚本转换

参考自:https://github.com/linsalrob/EdwardsLab/blob/master/bin/gfa2fasta.sh
(假设命名为gfa2fasta.sh

#!/usr/bin/env bash

if [ -z $1 ]; then 
    echo "Usage: $0 <gfa file>  >  <fasta file>"
    exit $E_BADARGS
fi

awk -v id=$1 '/^S/{print ">"id"_"$2"\n"$3}' $1

运行:

bash gfa2fasta.sh test.hifisam.bp.hap1.p_ctg.gfa > test.hifisam.bp.hap1.p_ctg.fasta

用gfatools转换

先用conda安装gfatools

conda install gfatools

接下来也写一个如上的脚本方便使用(假设命名为gfa2fa.sh

#!/usr/bin/env bash

if [ -z $1 ]; then 
    echo "$0 <gfa file>"
    echo "It will generate the fasta file like ${1%.*}.fasta"
    exit $E_BADARGS
fi
gfatools gfa2fa ${1} > ${1%.*}.fasta

使用起来就很简单啦:

bash gfa2fa.sh test.hifisam.bp.hap1.p_ctg.gfa

私货时间

  1. 我曾经尝试过把多个PacBio HiFi的gz压缩文件直接用cat给连接成一个文件然后去组装,但是组装结果和保持多个文件一起输入给hifiasm的结果是不一样的。cat到一起的结果的N50要略低于把三个文件分别输给hifiasm的版本。我和老板冥思苦想了好一会儿没想明白原因在哪儿,也就作罢了。所以建议如果有多个文件,就保持多个文件,不需要合并成一个。反正hifiasm是可以接受输入多个文件的。
  2. 建议有HiC的情况下尽量加上HiC数据,结果会更好一些。
  3. PacBio HiFi的组装结果不需要再经过其他的软件polish,不像传统的三代组装软件那样流程冗长而繁琐。

Hifiasm does not need polishing tools like pilon or racon, either.

引用自:(https://hifiasm.readthedocs.io/en/latest/index.html)

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

推荐阅读更多精彩内容