组装算法：为什么是k-mer?

Q:什么是组装算法？

A:无论是一代Sanger、二代Illumina还是三代Pacbio，其得到的测序数据（reads）相较于整个基因组而言都是极短的，基因组组装的任务就是将这些小片段连接起来，通过这些序列的关系构建Graph，然后根据算法从Graph中得到最优路径，从而得到最初的Contig序列。目前组装软件常用的两种算法：overlap-layout-consensus (OLC)和de-bruijn-graph(DBG)。

OLC算法

OLC是一种直观算法（intuitionistic assembly algorithm），主要用于长的低丰度序列的组装，特别是一代和三代测序数据，常见的软件有Arachne、Celera Assembler、CAP3、PCAP、TIGR、AMOS、Phrap、Phusion和Newbler。其方法主要分为三步，如下所示：

①Overlap：对所有reads进行两两比对，找到片段间的重叠overlap；

②Layout：根据得到的重叠信息将存在的重叠片段建立一种组合关系，形成重叠群，即Contig；

③根据构成每个Contig的片段的原始质量数据，在重叠群中寻找一条质量最佳的序列路径，并获得与路径对应的序列，即Consensus，如下图所示。

寻找路径的方法：将每条reads用一个节点代替，如果u的末端与w的首端存在overlap即创建一个有向连接directed-edge（u，w），这样一个重叠群的reads就形成一个网络，组装的过程就可以理解为在网络中寻找一条最短路径，使得可以经过所有的点并且每个点只经过一次，也即一个哈密顿路径（Hamiltonian path）问题。此算法的缺陷在于哈密顿路径本身所带来的NP难题（对于一个网络图中是否存在一条哈密顿路，没有可靠的充分必要条件），从而导致解决问题的时间过长。

DBG算法

DBG是一种非直观算法（anti-intuition algorithm），主要用于短的高丰度片段的组装，特别是二代测序数据，常见的软件有Velvet、ABySS、AllPath-LG、SOAPdenovo、EulerUSR、IDBA、SPAdes、Ray。与OLC算法不同，DBG算法将组装过程转换为一个在De Bruijn图中寻找欧拉路径（Eulerian path）的问题（从某点出发经过且只经过一次所有的边），而欧拉路径是P类问题，即有可靠的充要条件证明欧拉路径的存在，可以在确定的时间内解决。其方法如下所示：

①将reads分割为更短的长度统一的k-mers(长度小于k的reads将被舍弃)；

②寻找k-mer之间的重叠关系，建立De Bruijn图，即对于任意两个k-mer，如果u的后k-1个碱基序列与w的前k-1个碱基序列相同，则建立一条由u指向w的有向边；

③在De Bruijn图中寻找欧拉路径来获得结果序列Contigs。

由于二代测序得到的reads长度较短，包含的信息量较少，因此完成基因组拼接需要较高的覆盖度。OLC算法适用于读长较长的序列组装，通过构成的OLC图寻找Consensus sequence的过程，实际上是哈密顿通路寻找的问题，算法简单但是解决难度大。采用DBG算法，通过k-1的overlap关系，构建DBG图，通过寻找欧拉路径得到Contig序列，算法可靠性更高，两者的区别如下所示：

DBG算法将哈密顿路径问题转化为欧拉路径问题的关键在于De Bruijn图中每一条序列只比前一条序列多1个碱基，overlap为k-1个碱基，因此只需要知道起始节点和终止节点的序列以及中间所有边（overlap序列）便可得出contig。

组装流程

根据overlap构图并搜寻contigs仅仅是第一步，完整的基因组组装流程还包含Scaffold搭建、gap修补等过程。以经典的SOAPdenovo软件为例，其为一种新型的利用DBG算法的short read组装软件，设计服务于大型的植物和动物基因组，对细菌和真菌的基因组也有效。其组装分析过程如下所示：

⑴分割k-mers

将所有reads分割成固定长度为k的较短序列，形成等长的k-mers，长度短于k的reads将被舍弃。由于reads中仍有一些错误或者N的存在，造成一些错误k-mer或者低频k-mer，错误Kmer对后续组装会产生很大的困扰，因此在构建DBG图之前，需要先对数据进行纠错区分。SOAPdenovo采用绘制没有测序错误与错误k-mer的深度分布图的方法，如下所示：

Error free代表没有测序错误的Kmer频数分布，Error rate1%代表有1%错误率的Kmer频数分布。

⑵构建DBG图

根据k-mers之间的overlap关系划分Contigs，并对每一个Contig创建de bruijn graph。实际生成的DBG图非常复杂，例如由于测序错误、杂合或者高重复序列产生的tips翼尖（a）和bubbles气泡（c），由于，低覆盖率造成的链接low coverage links（b图）和（d图），微小重复造成的压缩(e图)，具体如下所示：

因此需要对DBG图进行简化，移除错误链接、删除低覆盖链接、解开短重复序列等，最终输出Contigs序列。

⑶构建Scaffold

双末端测序中，相同序列名字的read1和read2来自同一条序列（中间不一定有overlap），可以根据paired-end信息将不同的Contigs搭建成Scaffold，如下图所示：

⑷填补gap

得到的Scaffold中间会有较多的gap，为了使组装的序列更完整，需再次利用测序的双末端数据之间的配对关系连接Contigs，并利用测序数据与已经组装的Contig之间的覆盖关系对Contig之间空隙进行补洞，延长Contigs。

参考文献

[1] Z L, Y C, D M, et al. Comparison of the two major classes of assembly algorithms: overlap-layout-consensus and de-bruijn-graph[J]. Briefings in Functional Genomics, 2012, 11(1): 25.

[2] Luo R, Liu B, Xie Y, et al. SOAPdenovo2: an empirically improved memory-efficient short-read de novo assembler[J]. Gigascience, 2012, 1(1): 1-6.

END