DNA特征提取之序列组装K-mer算法

基因就像一本天书，里面的字都是有A,T,C,G组成的。我们测序出来的大于几百bp的小的序列，我们叫做reads。我们测序的结果中有上百万条的reads，而这些reads的位置我们又不知道，我们只能根据他们重叠的部分来尽量还原他的原型。

k-mer是指将reads分成包含k个碱基的字符串，一般长短为m的reads可以分成m-k+1个k-mers。举个例子吧：为了简化，有这么个reads（当然实际比这个长）：AACTGACTGA.如果k-mer的k为3的话，我们可以将其切割为AAC  ACT  CTG  TGA GAC ACT CTG TGA。我们将这些k-mers放入计算机中拼接，假设第一个为TGA ,那么下一个应该为GA-,.……

基于这样的思路，我们很快就发现了问题，下一个点可能有很多的选择，或者没有选择

       我们需要找到Hamiltonian path，我们需要找到包含每个点的，但是只包含一次。下图像不像我们小学之前做过的游戏，遍历每个点，但是每个点只能经过一次。

这是俄国一个我一个伟大的科学家William Hamilton的一个发明。

但是又提出了一个假设，如果有两条或多条的Hamiltonian path呢？如何才能知道其中的一条是DNA的序列呢

我们上面提到的是3个碱基为一个node，现在我们就只要其中的2个来作图，然后将作出的图中相同的node合并，

这是两种算法思想，在Eulerian Path Problem中，visit every edge of the graph exactly once.
而在Hamiltonian Path Problem,visit every node exactly once.第一种算法更好实现，所以我们接下来讲关于Eulerian Path Problem的de Bruijin graphs。

在实际组装基因的时候，我们知道的是reads和k-mers，通过这个我们来基于Eulerian Path Problem来构建de Bruijin graphs，然后找Eulerian Path。可是会有很多的de Bruijin graphs，或者一个de Bruijin graphs有几个Eulerian Path。为了减少contigs,发明了read pair sequencing

把很多拷贝的相同基因，尺寸随意剪切为大的相同大小InsertLength片段。产生read-pairs:两个reads来自每个片段的末尾。A paired k-mer就是两个k-mer距离相隔d,

根据我的实际使用经验，如果你的read足够长，覆盖度足够高，kmer设的越高越好。
但是实际情况是，测序的覆盖度经常不够，或者用早期的GA平台测出来read长度只有35bp，或者为了节省成本，在mate-pair library(长片段insert的文库，一般>2kb)测序时双端只有70bp,甚至40bp之类的，情况比较复杂。

       一般来说，我尽量使用更高的kmer，如果我有100bp的pair-end,50bp的mate-pair,而且覆盖度挺高，我就用到kmer=45左右，如果mate-pair只有40bp，kmer=35左右。如果mate-pair更短，只有35bp，kmer值就再降一点。

转载：http://qinqianshan.com/k-mer-de-bruijin-graphs-based-on-algorithm-of-sequence-assembly/