hifiasm续：为什么我的组装结果不可重复？结果居然会受顺序影响？

前情提要

之前在介绍hifiasm的时候提到过：

我曾经尝试过把多个PacBio HiFi的gz压缩文件直接用cat给连接成一个文件然后去组装，但是组装结果和保持多个文件一起输入给hifiasm的结果是不一样的。cat到一起的结果的N50要略低于把三个文件分别输给hifiasm的版本。我和老板冥思苦想了好一会儿没想明白原因在哪儿，也就作罢了。所以建议如果有多个文件，就保持多个文件，不需要合并成一个。反正hifiasm是可以接受输入多个文件的。

推文发布之后我有一个师妹联系我说她也遇到了差不多的情况，多个cell测到的PacBio的HiFi reads在按照不同顺序输入给hifiasm之后得到的结果是不一样的。我表示非常震惊。原来不止我遇到了这样的问题呢？？

于是我就到hifiasm的GitHub上用order作为关键词看看有没有其他人也有类似的提问，果然有。

顺序会影响结果吗？会！

https://github.com/chhylp123/hifiasm/issues/291

在这个issue里使用者遇到了一样的问题，组装的结果不可重复，于是开了个issue。

下面李恒大神亲临现场：

This is expected.

李恒大大说这样的结果是预料之中的。即出现这样的差别并不奇怪。在后面hifiasm的仓库拥有者chhylp123也进行了回复：

The order matters. If the order is exactly the same, the assembly will be also exactly the same.

顺序是重要的，如果顺序不同，结果是会不一样的。

后面又进一步解释了一下：

There are millions of HiFi reads involved in genome assembly. In rare cases, selecting any of two or more reads for assembly is fine. In this case, hifiasm may select the reads with the smallest ID. I personally think the order does not matter too much. The difference will be very minor.

基因组组装涉及数百万个HiFi读数。在极少数情况下，选择两个或多个读数中的任何一个进行组装都可以。在这种情况下，hifiasm会选择ID最小的reads。我个人认为顺序没有太大关系。差别将非常小。

我自己做了个测试

于是乎我用我自己手头的HiFi数据在学校的HPC上做了个测试，保持其他的参数不变，仅改变HiFi数据输入的顺序：

hifiasm \
-o test1 \
-t 40 
--h1 hic_R1.fastq.gz \
--h2 hic_R2.fastq.gz \
hifi_reads1.fastq.gz \
hifi_reads2.fastq.gz \
hifi_reads3.fastq.gz

根据排列组合，把六个组合都运行了一遍。但是不知道啥原因，在得到初步结果的test1.hic.p_utg.gfa之后我的任务就被HPC杀掉了。但是好在每个测试组都生成了test1.hic.p_utg.gfa文件，把他们用gfatools转换成fasta数据之后统计：

cat *stat | grep "N50_length"
test1    N50_length    533,342    2514
test2    N50_length    538,212    2506
test3    N50_length    530,652    2541
test4    N50_length    530,928    2530
test5    N50_length    535,658    2523
test6    N50_length    538,754    2520

可以看到，确实每个版本的N50都是不一样的，总长度也不一样：

cat *stat | grep "Total_sequences"
test1    Total_sequences 4,772,977,416 27205
test2    Total_sequences 4,773,102,026 27204
test3    Total_sequences 4,775,409,060 27372
test4    Total_sequences 4,774,439,816 27281
test5    Total_sequences 4,772,800,836 27184
test6    Total_sequences 4,772,966,066 27202

看来如果当时在做初步组装的时候我换个输入顺序，可能得到的初步组装的结果又会不一样。啊，人生总是这样的无常。

萌哥碎碎念

1. 虽然作者说差别应该不会很大，但实际上我之前的结果，一个组装版本的N50是28 Mb，另一个是26 Mb，感觉也不算是很小的差别。
2. 根据作者的介绍，在选择任意一个reads都可以进行组装的情况下，hifiasm会选择ID号更小的那个，那是否意味着我们可以把所有的reads组成一个文件然后根据长度逆向排序后重新赋予ID呢？当然这只是我拍脑袋想到的，没有经过测试哦。
3. 一个瞎指挥想法：如果条件允许，有多个HiFi文件的情况下，可以跟我一样把各种排列组合都做一遍，看最终结果N50的高低来决定最终的组装版本可能也是种可行的方案。
4. 这个发现解答了我之前由于把多个文件结合成一个文件输入hifiasm组装得到较差结果的疑问，虽然不可能把现在的项目全部推倒重新再做一遍了，但是心里踏实多了。虽然或许只是增加了一个无用的小知识。