「三代组装」Pacbio组装后如何用自身数据进行polish（更

之前那我由于需要对PacBio的组装结果进行polish，于是写了「三代组装」Pacbio组装后如何用自身数据进行polish。最近发现自己又有了需求，于是重新回顾了我之前写的这篇文章，但是在实践的时候，却发现新版本pb-assembly(0.0.8)没有之前的pbalign和variantCaller这两个命令，也就是意味着我之前的文章随着版本更新而失效了。

此时我就面临着两个选择，一个用conda安装之前的版本，逃避可耻但有用。另一个就是探索新版本下的代码，勇敢面对惨烈的人生。由于时间不是特别的紧张，所以我选择了第二条路。

通过我一波检索，我发现第二条路是一条没有人走过的路。https://github.com/PacificBiosciences/pb-assembly只讲了falcon组装后，使用falcon-unzip进行polish，而没有介绍使用其他软件(如Canu, Mecat2) 应该怎么做。https://github.com/PacificBiosciences/GenomicConsensus提到的quiver命令并不存在于我们pb-assembly中。

我虽然把pb-assembly下的命令都看了个遍，但是也没有想到哪些命令可以组合在一起用来polish。但是我突然想到，既然falcon-unzip里面会有polish，那么的它运行日志里面一定会有它如何进行polish蛛丝马迹。经过一波的分析和追踪，终于被我找到了相关的代码。修改之后如下

set -e
REF=$1
BAM=$2
THREADS=80

# align
pbmm2 align --sort  $REF $BAM aln.bam
pbindex aln.bam
# arrow
gcpp --algorithm=arrow -x 5 -X 120 -q 0 -j $THREADS -r $REF  aln.bam -o cns.fa,cns.fq,cns.vcf

新的脚本不再用BLASR进行比对，而是拥抱了minimap2，提高了比对速度，并且避免了之前BLASR可能会出现的内存泄露问题。使用gcpp替代了之前variantCaller，但是参数几乎没有变动，你可以认为就是改了一个名字而已。

另外这是一个简化的脚本，如果是多个bam文件，则需要分别比对然后合并。原来的falcon流程会将参考基因组按照染色体拆分然后并行处理，但是我们的脚本设置了更高的线程数，最后结果速度也差不多。