hifiasm(高质量组装PicBio HiFi数据,同样也是可
2023-07-29 本文已影响0人
GenomeStudy
hifiasm
如果你基因组是简单的二倍体,不复杂,不是高重复或者高杂合,测序数据还是PicBio HiFi数据,选择hifiasm进行组装是个不错的选择,对于处理PicBio HiFi数据,hifiasm有快、准、质量高的特点,还能分单倍型!
整体的流程
hifiasm.png1.从安装开始说起
# Install hifiasm (requiring g++ and zlib)
git clone https://github.com/chhylp123/hifiasm
cd hifiasm && make
记得添加环境变量就可以直接调用了。
2.使用(有多种模式可选择,提供三种常用模式)
1)在只有PicBio HiFi数据的情况下,也是最简单的
# Assemble inbred/homozygous genomes (-l0 disables duplication purging)
hifiasm -o CHM13.asm -t32 -l0 CHM13-HiFi.fa.gz 2> CHM13.asm.log
# Assemble heterozygous genomes with built-in duplication purging
hifiasm -o HG002.asm -t32 HG002-file1.fq.gz HG002-file2.fq.gz
-o 输出文件的前缀
-t 运行程序设置的线程
-l0 不进行purge
CHM13-HiFi.fa.gz hifi数据
模式不同,输出文件的后缀也会有所不同
在这个模式下输出的文件以前缀.bp.后缀
的形式输出。
#用awk转化gfa格式为fa格式,即得到组装的contig文件
#主文件
awk '/^S/{print ">"$2;print $3}' prefix.bp.p_ctg.gfa >prefix.bp.p_ctg.fa 2>2.log
#hap1文件,一型文件
awk '/^S/{print ">"$2;print $3}' prefix..bp.hap1.p_ctg.gfa > prefix..bp.hap1.p_ctg.fa 2>2.log
#hap2文件,二型文件
awk '/^S/{print ">"$2;print $3}' prefix..bp.hap2.p_ctg.gfa > prefix.bp.hap2.p_ctg.fa 2>2.log
2)有PicBio HiFi数据和HiC数据的情况
# Hi-C phasing with paired-end short reads in two FASTQ files
hifiasm -o HG002.asm --h1 read1.fq.gz --h2 read2.fq.gz HG002-HiFi.fq.gz
-h1 HiC一端的数据
-h2 HiC另一端的数据
在这个模式下输出的文件以前缀.hic.后缀
的形式输出
#同样的,用awk进行格式的转换
awk '/^S/{print ">"$2;print $3}' prefix.hic.p_ctg.gfa >prefix.hic.p_ctg.fa 2>2.log
3)有PicBio HiFi数据、HiC数据和ONT数据的情况,也是最容易达到T2T级别组装的模式
# Single-sample telomere-to-telomere assembly with HiFi, ultralong and Hi-C reads
hifiasm -o HG002.asm --h1 read1.fq.gz --h2 read2.fq.gz --ul ul.fq.gz HG002-HiFi.fq.gz
--ul ONT测序数据
3.还有其他复杂的参数,可以去学习,简单说几个
--hom-cov
--hom-cov INT homozygous read coverage [auto]
这个参数的使用,首先是你的基因组是一个杂合或者高重复的基因组,在进行上述的组装后,得到的组装结果比实际的偏大,或者分型结果,即hap1和hap2的结果相差较大。你可以看看日志文件,即*.log文件。在log文件中有一行[M::purge_dups] homozygous read coverage threshold: X.
的描述,这个X的值即可设置为--hom-cov X
-s
-s FLOAT similarity threshold for duplicate haplotigs in read-level [0.75 for -l1/-l2, 0.55 for -l3]
这个参数的使用,和--hom-cov
类似,在你得到的组装结果比实际的偏大的情况下,可以调整-s
,程序默认为0.5
,偏大的情况下可以往下调。
-n-hap
--n-hap INT number of haplotypes [2]
这个参数可以调整你需要分型的个数,即如果是四倍体材料,就有4个单倍型的基因组,可以试试--n-hap 4
,反正我组装三倍体的材料,设置了--n-hap 3
并没有成功。
--hg-size
--hg-size INT(k, m or g) estimated haploid genome size used for inferring read coverage [auto]
这个参数可以输入,预估的基因组大小,例如:--hg-size 500m
大多数的情况下,用默认的参数就是足够了的。
参考连接
https://hifiasm.readthedocs.io/en/latest/index.html
https://github.com/chhylp123/hifiasm