基因组拼接原理

生信喵

一、拼接原理
      序列拼接也叫做基因组组装。高通量数据分析一共有两条路线，一条路线是短序列比对，而另一条路线则是序列拼接。
      基因组上完整的 DNA，经过提取建库，最终测序，得到散在的 reads。这个 reads 可以是等长的，例如 illumina 测序的双末端 300bp，也可以是三代测序 不等长的片段。将测序得到的reads 拼接起来，得到一个物种的全基因组序列了。这个过程就是基因组拼接。
      拼接的操作并不难，直接输入给软件，经过大量数据处理过程即可得到结果。但是序列拼接是一项非常复杂的工作。我们常说基因组是一部天书，比如人全基因组有 30 亿个碱基对，那么要是打印出来，真的是一部天书。我们准备这样十本书，也就是基因组的十倍，然后我们将这些书随机撕成 100 个字符长度的一张张纸条，然后我们需要再将这些字条拼接起来，组合成原来的天书。这将是一项非常复杂的工作。需要非常大的计算。


1.1 拼接算法
      目前基因组拼接的算法主要有两种，一种是适用于长读长测序的基于有 overlap 的共有（consensus）序列组装的方法。这种算法被称为 Overlap Layout Consensus，简写为 OLC算法，另一种是用于短读长测序的基于 De Bruijn 图的短序列拼接。
      无论是 OLC 还是 De Bruijn 图的算法，基因组拼接的基本原理都是序列之间存在 overlap 关系。

1.1.1 OLC 算法主要分为三步：
      （1）Overlap：，对所有 reads 进行两两比对，找到片段间的重叠信息；
      （2）Layout：根据得到的重叠信息将存在的重叠片段建立一种组合关系，形成重叠群，即Contig；
      （3）根据构成 Contig 的片段的原始质量数据，在重叠群中寻找一条质量最重的序列路径，并获得与路径对应的序列，即 Consensus。
      OLC 算法最初成功的用于 Sange 测序数据的组装，比如 Celera Assembler，Phrap，Newbler等均采用该算法进行拼接组装。不过nanopore的测序数据与传统sanger测序有明显的不同，测序质量明显低很多，这就需要的是在组装前进行测序 reads 的纠错。Nanopore 三代数据进行纠错后得到高准确性的数据。对 Clean Data 进行纠错后，基于纠错后的数据用 Canu、SMARTdenovo、WTDBG2、Miniasm、MECAT、SPAdes、Unicycler、MaSuRCA 等软件分别进行组装。


1.1.2 De Bruijn 图算法
      基于 De Bruijn 图算法主要分为四个步骤。首先是构图（Pregraph），将 reads 切割成 kmer。然后是利用 kmer 之间的 overlap 关系构建 contig，接下来是利用 reads 之间的 pairend 关系将 contig 连接成 scaffold，最后补洞这四个步骤。


1.2 名词解释
      1、reads：直接测序出来的片段就称为 reads，不同测序平台测序出来的我们都可以称为 reads。
      2、pairend reads：illumina 高通量测序的是双末端测序，reads 具有物理上的关系。小片段文库的这种物理关系叫做 pairend，大片段文库的这种关系称为 matepair 文库。
      3、insertsize：illumina 双末端测序中间物理距离叫做 insertsize，其实就是随机打断后用来测序片段的长度。在序列拼接中会用到 reads 之间的这种关系，这个 insertsize 的值也非常重要。
      4、overlap：reads 之间具有共有序列，这个共有序列关系成为 overlap；
      5、kmer：将 reads 切割成固定长度的小片段。这样原来的 reads 就变成了很多的更短的 kmer。
      6、contig：通过 kmer 之间的 overlap 连接成更长的片段，就是 contig。contig 是片段重叠群的意思。有时候我们也直接说 reads 通过 overlap 连接成 contig，都可以。准确来说是利用 kmer 连接成 contig；
      7、scaffold：congtig 之间我们可以通过 reads 之间的 pairend 关系再连接成更长的片段，这个就叫做 scaffold，
      8、gap：contig 连接 scaffold 过程中，中间区域使用 N 碱基来填充，称为 gap。


二、利用 kmer 估计基因组大小
      获取一个物种的基因组大小是一件非常重要的工作。在测序之前知道了基因组的大小，我们就可以根据这个值来决定测序的数据量。例如人的基因组是 3G，我们测 10 倍的数据量就刚好是 30G。而且也可以根据基因组大小，估计 denovo 基因组软件对机器内存的要求。如果基因组太大。例如 10G，序列拼接就会超出了机器内存的限制，就得想办法提高硬件。也需要在序列拼接之后将序列拼接值与真实值之间做比较，来评估序列拼接的效果。所以，获取基因组大小是非常重要的。
      那么如何来获取基因组的大小呢。一般物种的基因组大小可以从（http://www.genomesize.com/ ）这个数据库查到。如果没有搜录，需要考虑通过实验方法，例如利用流式细胞仪来估计基因组大小。也可以采用定量 pcr 估计基因组大小。
      不过利用实验的方法显然非常的复杂，需要很多的操作，这里面我们不采用实验的方法，而是基于现有测序数据，基于数据分析方法，利用 kmer 分析来估计基因组大小。也就是不通过序列拼接，就预测出基因组大小。


2.1 什么是 kmer？
      所谓 k-mer，即为一段短的 DNA 片段。K 为一个奇数，k 等于几，就称为几 mer。例如一段17bp 的序列，则可以称为 17mer。Kmer 这个概念在后面序列拼接中依然会用到。

      
      我们可以将一条 reads 切割成很多小的 kmer 片段，从第一个碱基开始，每隔固定距离的碱基开始提取碱基。例如一条 100bp 长的 reads，每隔一个位置取一个 17mer 的片段。也就是1-17 取一个 kmer 出来，2-18 取一个 kmer，3-19 取一个 kmer。以此类推，最终 84-100为 1 个 kmer。那么最终将会生成将会生成 100-17+1，也就是 84 个 kmer 片段。原来一条长片段，变成了很多短的片段，碱基的数目也增加了很多倍。而且，每次取 kmer 是同一条reads 正反取两次，也就是对这条 reads 的反向互补序列再取一次 kmer。


2.2 为什么二代测序必须要取kmer？
      
      像上图黄色的kmer的包含错误的，要是测序错误越多，黄色的kmer越多，接下来就是排除掉黄色的kmer，这样就没有错误的kmer了。
      
      理论支持kmer频数为1的就包含测序错误。