卖萌哥

V1

2022/10/11阅读:92主题:橙心

用hifiasm组装基因组

用hifiasm组装基因组

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

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

可用资源

  • 文章地址:https://www.nature.com/articles/s41592-020-01056-5
  • GitHub地址:https://github.com/chhylp123/hifiasm
  • 官方教程:https://hifiasm.readthedocs.io/en/latest/index.html
  • 洲更学长的笔记:https://xuzhougeng.top/archives/assemble-hifi-pacbio-with-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)

分类:

工具介绍

标签:

开源软件

作者介绍

卖萌哥
V1