Canu: scalable and accurate long-read assembly via adaptive k-mer weighting and repeat separation

Abstract

对于 SMS read,准确高效的组装大的重复块和密切相关的单倍体仍然具有挑战性。CANU 包含新的 overlapping 和 assembly 策略,包含 tf-idf weighted MinHash 和 构建稀疏组装图可以避免折叠发散的重复区域和单倍体。

Introduction

有时基因的重复区域可能会比 overlap 区域的长度长,这会导致基因重建具有不确定性,并组装会碎片化。Parametric Complexity of Sequence Assembly: Theory and Applications to Next Generation Sequencing 有两个方法能克服这个基础的限制:增加有效的 read 长度或者基于拷贝的特异性变体分离不精确的重复。

本文的目标是重新构建整个基因组,着重于分层策略,它能够生成最具有连续性的从头组装。分层方法不需要其他的 read,而是在组装前使用多轮 read overlapping (alignment) 和校正来提高 SMS read 的质量

Results

Architecture

CANU 包含几个新特征:计算资源探测,自适应 k-mer 加权,自动估计错误率,稀疏图构建,图的片段组装输出

Canu 包含三个阶段:纠错(绿色),修剪(红色),组装(紫色)。Canu 各个阶段共享二进制磁盘存储和并行存储构造的接口。所有阶段的第一步都是,索引输入序列,生成 k-mer 直方图,构建 all-vs-all overlap 的索引,整理总结分析。
纠错步骤,首先选择最好的 overlap 用于纠错,估计纠错后 read 的长度,生成纠错后的 read
修剪步骤,确定输入不支持的区域,然后修剪或分割 read 到它们能支持的最长范围
组装步骤,确定序列的错误,构建最好的 overlap graph,输出 contigs,组装图和总结分析  

Adaptive MinHash k-mer weighting

处理 read 间的 overlap 一般有两个阶段,首先构建一系列具有相似性的 read pairs,然后对它们进行更直接的比较。在第一步通过确定公共的 k-mer 一般都能找到候选的 overlap,然而重复区域减少了 k-mer 的分布熵,并且一些 k-mer 的出现频率的大幅增加导致了候选 overlap 的增加,这些必须由更消耗资源的第二阶段进行处理。通用方法是掩码低复杂性序列或在索引的时候忽略高重复的 k-mer,但这可能会导致探测不到部分正确的 overlap。

Canu 使用一个弹性的方法去处理重复的 k-mer,它概率性的减少(但不去除)高重复性的 k-mer 被选择生成 overlap 的机会 。Canu 使用 MinHash Alignment Process (MHAP) 比较整个 read 的压缩草图,因为 MinHash 草图包含一个固定大小的k-mer( read 中选择的 ) 的子集,草图中包含特定 k-mer 的概率是可以调整的。例如,一个在基因中重复性很高的 k-mer 应该有更小的权重,因为它携带的关于原始 read 的相关信息较少。相反,在一条 read 中出现多次的 相对独特的 k-mer 应该含有更高的权重,因为他能代表更长的 read,这些 k-mer 的相对重要性在 NLP 中称为 tf-idf weight

在 read overlapping 上,权值就是,单个 read 中 k-mer 出现的次数和 k-mer 在所有 read 中的整体的稀有性。比如文档的相似性,在确定相似文档的时候,一个稀有的词语在文档中出现多次说明他是一个好的相似文档。对于overlap,重复的 k-mer 有小的权值,通过减少草图中重复的 k-mer 出现的次数,可以是 k-mer 的索引分布更加均匀,这减少了草图中的无意义信息。更重要的是,这是有概率实现的不是屏蔽重复,重复 read 之间真正的 overlap 仍然可以被探测到。

Best overlap graph

Canu 使用 “best overlap graph” (BOG) 算法 的变体,构建了一个稀疏的 overlap 图。The fragment assembly string graph 把整个 overlap graph 都加载到内存中了,非常消耗资源。相对来说,本文的贪婪算法只加载最好(长)的 overlap 到内存中去,当 read 足够长这种算法是最优的。然而,当重复区域的长度超过 overlap 的长度时,算法可能会被误导,导致出现错误的组装,Canu 的新算法 Bogart,通过统计性的过滤由于重复引起的重叠,并回顾性检查图形中的潜在错误来解决这个问题。

原始的 BOG 方法,在所有的 overlaps 中选择错误率低于用户指定值的 overlap,overlap 的错误率是编辑距离除以 overlap alignment 的长度。因此,阈值必须确定的足够低,避免引入错误的 overlap,因为序列的错误,为了探测正确的序列又必须设的足够高。在新的 Bogart 算法中,最优的错误率参数通过全局和局部数据自动估计,因为原始的 SMS read 错误率达 10%-20%,导致这具有挑战性,因此在构建图之前,Canu 进行多次纠错步骤,在这些纠错之后,剩余的读取错误是从最长的 overlap 中估计出来的。整个 overlap 集合然后过滤到只包含那些错误率小于全局错误率中位数的 overlap,然后会用这个子集重新计算最长的 overlap。低错误率的阈值有效的移除了大部分不对的 overlap,使贪婪算法能构建一个干净的最好的图。在图中,从最长的没有分支的路径中构建最初的 contigs。

overlap 错误率估计,重复区域确定和分割
A 所有最好的边的错误率
B 虚线表示全局错误率阈值(1.6%),此图代表这次组装中最大的 contig 计算出来的局部错误率
C 中间的双向黑线是 contig ,紫色代表的是重复区域,一条 read (蓝色)能横穿整个重复区域表明 contig 重建是正确的
D 这种情况下,没有单条 read 能横穿整个 重复区域,最初的 contig 是通过两条蓝色的 read 构建的,如果由两条蓝色 read 构建的 contig 不是明显的比任何一条蓝色到红色 read 构成的 contig 效果好,就需要把这个 contig 进行拆分

构建和评估后,Canu 以 Graphical Fragment Assembly (GFA) 格式提供一个组装图,这个图类似一个稀疏的 overlap 图,一些染色体组装成了单个的 contig,但这个图显示了那些组装中更复杂的,未解决的重复区域的结构,Canu 的输出定位了这些复杂的结构,但重复区域的长度阻碍了完全组装,将 Canu 组装图与补充的其他信息相结合,例如来自光学或染色质接触映射可以帮助识别正确的路径并解析此类结构。

Canu 的 GFA 输出,可以定位复杂的区域,并且能够提升 scaffold
A Bandage 绘制输出,Nodes 是 contigs,edges 表示 contigs 间没用到过的 overlap,最大的 contigs 用随机的颜色标注
B 唯一的 contig 用黑色阴影表示,红色阴影表示 重复区域的 contig

Methods

MinHash overlapping

  1. refer Assembling large genomes with single-molecule sequencing and locality-sensitive hashing
  2. 使用MHAP计算所有read之间的overlap,他有两个阶段,第一个阶段找到那些有可能有overlap的一对,
  3. TF-IDF算法 :TF(词频) x IDF(逆文档排序) ,简单快速容易理解,但重要性不够全面
  4. 第一阶段用了IF-IDF来加权k-mer,第二阶段用bottom sketch多重hash

Parallel overlap sort and index

  1. bucket sort

Read correction

Canu 使用 all-vs-all overlap 对单条read进行纠错,因为简单的使用所有 overlap 计算每个 read 的 overlap 会导致掩盖特定重复区域的拷贝变体。因此对纠错单条 read 使用的 overlap 进行过滤,有全局过滤和局部过滤,策略的目的是克服序列质量和重复区域所造成的偏差。比如,高于平均质量的read 更可能决定纠错结果,而不管它是否是处于正确的重复区域,为了避免这个,每个read只能被其他read进行纠错使用 C 次,C read 的期望深度。全局过滤器对每个 overlap 打分(length * identity),对每个 read 只留下最好的 C 个overlap,从而会将来自同意副本区域的相似的重复 read 聚集到一起。当错误是均匀分布的时候,因为来自同一重复区域的 read 的比来自不同重复拷贝你区域 read 的总体差异性更小,所以来自同一区域的更有可能会分到同一组。

在进行纠错计算之前,每对 overlap 都会用于计算每条read 纠错后的期望长度,根据这个估计,选择最长的且达到用户制定的覆盖深度的 read 来进行纠错,使用修改后的 falcon_sense 算法进行纠错,对于给定的进行纠错的 read,使用 Myers’ O(ND) 算法进行 align,然后根据 alignment 生成有向无环图,根据权重最高的路径生成纠错后的序列,当没有足够的证据进行纠错时,会分割原始 read

Overlap Based Trimming

对纠错后的 read 重新计算 overlap,移除那些不被其他read所支持的区域,

Overlap Error Adjustment

在修剪后和构建图之前,Canu再次计算overlap 对序列错误做最后一次检测,Holt et al. 算法的目的是提高真实有差别序列的分离性(不同的重复区域和单倍体)和随机测序c错误倒追的错误的差异。每条read 通过 overlapping alignment 的多数票进行纠错,只有当有其他充足的read支持这个不同的碱基时,才保存这个碱基。每个overlap的错误率是根据解决序列错误后生成的alignment进行调整的。算法需要使用两次 overlap,第一次用于检测可能的测序错误,第二次暂时应用这些修改,重新计算 alignment 更新计算的错误率。

Graph Construction

Bogart 使用基于的 The fragment assembly string graph 的 best overlap graph 策略,在 Bogart 中,最好的 overlap 是经过多重过滤后选择的,过滤包括不正常的高错误的 overlap,潜在的嵌合 reads 和 read 的overlap显示可能不是正常序列的,这会生成一个更精确和整洁的图。

在纠错、修剪、overlap修正后,所有的 overlap 被用来选择一个最好的边的集合,用这个集合计算 overlap 错误率的中值和中值绝对偏差,这个分布代表了之前纠错步骤的错误率剩余,低的平均 overlap 错误率,代表了好的序列数据和之前的纠错过程是成功的,然后自动计算 overlap 错误率的一个最大限制值,超过这个值的 overlap 不在这个图的构建中使用,这个是算法有能力区别紧密的相关read 和单倍型。

为了过滤不好的 overlap, Bogart 过滤那些没有经过合适修剪和纠错的可疑 read

剩下的 reads 和最好的 overlaps 集合构成了最好的 overlap graph,然后根据这个图构建出了最初的 contig

当由于单倍型的差异导致特定基因上某个特殊的区域不止一个重建时,就会出现组装气泡,小部分的不同,例如几十个碱基,在 overlap 中通常不会被探测到,因为相比于overlap的大小,这个不同就显得不是那么重要了,更大的差别,会导致形成两个几乎是冗余的 contigs 覆盖了相同的区域,具有更多 read 的单倍型通常会在跨越整个区域的的大 contig 中重建,而较少的仅作为变异区域。

尽管经过仔细的过滤,贪婪算法仍然易于出错,相较于全 string graph 所代表的,会缺少部分边,所以最后的步骤是增加缺少的边,并且打断那些组装好但不正确的 contigs。使用所有的成对的 overlap 信息,再次使用read放置机制和所有非气泡contig中的read,可以将每个组装好的 contig 使用相应的 read 位置注释。