从今天开始我们正式move on to基因组拼接的内容了。
我也是在边看边总结
有说得不正确的地方，欢迎各位同仁指正。
首先需要提一下在基因组拼接中不得不说的一个词儿——k-mer
mer这个词根在生物中用的也很广泛，比如说dimer，trimer表示的是构成单元。在测序中k-mer指的是一条字符串中所有可能具有长度为k的子串。
举个例子就好理解了
如下图，对于序列为ACTCGATGCTCAATG的reads，长度L为15nt，若K取5，那么能够产生15-5+1的k-mer，也就是L-K+1条k-mer。实际上每次对read的正反链都会取k-mer，大大增加了计算机需要处理的碱基的数量。
1
那么这个东西有啥用途呢？
1）估计基因组的大小
如图所示，k-mer的总数为L-K+1。基因组估计的误差百分数为(L-N)/L。可以看到随着L增加，基因组估计的误差逐渐减小，所以对于更大的片段大小，k-mer 的总数提供了一个很好的近似于实际的基因组大小的估计。
2
在实际的测序中，我们知道的参数有测序所得到的所有子序列为K的总数N，K值以及reads的读长L，reads的数量C。
那么N=C*（L-K+1）
因为用每个read的k-mer数量n去估计基因组大小，而n=L-K+1，所以基因组大小G=N/C。
这个C也叫做测序中的覆盖范围。但是在实际生活中，PCR总有其偏好性， 虽然拷贝数为C，但是C大概率是有重叠的，不能体现真实的覆盖水平。因此在实际中，需要用到平均覆盖深度去代替这个理论的覆盖深度C值。那么平均覆盖深度又怎么计算呢，这又得用到一个词——k-mer频率。
k-mer频率指的是，整个测序过程中产生的长度为K的子序列总数N,统计每个长度为K的子序列在N中出现的次数，这个频数就构成了横坐标。再把这个频数所对应的子序列的数值作为纵坐标。
比如3个子序列出现了1次，5个子序列出现了2次，10个子序列出现了3次，6个子序列出现了4次，那么所画出来的k-mer分布图如下。也可以将频率（频数/N）
3
近似呈现泊松分布，以峰值所对应的横坐标作为平均覆盖度，代入到上面的公式中。
但是在实际的情况下，通常会出低频数的地方出现下降的情况。
4
这是因为测序的错误，很多的K-mer只有孤立的几条，但是测序的累计错误还是很高的，所以对应的频率（也就是纵坐标）也会很高。
2）估计基因的复杂度
这里说的复杂度主要是指的杂合性，倍性，重复序列等等。
理想情况下，k-mer分布图只有一个主峰，但是一些基因组复杂度较高，会出现多个主峰的情况。
首先，若是测序的深度是2（测序深度是指测序得到的总碱基数与待测基因组大小的比值，可以理解为基因组中每个碱基被测序到的平均次数。 测序深度= reads长度× 比对的reads数目/ 参考序列长度。 假设一个基因大小为2M，测序深度为10X，那么获得的总数据量为20M），没有杂合情况下，某一段read的17-mer深度为1，如果某个位点是杂合，那么包含这个位点的read就有两个17-mer。在主峰为C的1/2处，会出现一个小峰，并且杂合度越高，小峰的纵坐标越大。
如果基因组是多倍体，特别是同源多倍体或者是相近物种杂交的多倍体，许多长序列片段甚至会完全相同，因此在2C/3C处会出现一个小峰。
但是如果是基因组中重复序列比较高，会在高深度k-mer数量增加，出现明显的拖尾。
5
3）样品污染评估
物种样品的污染可以通过两种方式来评估：k-mer分布和GC含量图来分析。
首先来谈谈k-mer分布图，一般来说取不同的k-mer得到的图呈现的趋势是相同的，原先是单峰的还是单峰，双峰的还是会呈现双峰。如果次峰的位置不在1/2或者2的位置，那么可能存在物样品污染。
其实是GC含量分布图。每个物种的GC含量都是固定的，上面的峰图指的是GC的含量，对应下来大概是40%左右。右边的峰图是k-mer分布，可以看到有两个峰，分别在20和40的深度，表示存在部分的杂合，左下的图表示，两个深度的GC含量都在40%左右，表明是同一个物种的样本。
6
如果是下面这样分布的，表明存在不同的GC含量，那么大概率就是存在污染。
7
参考链接：
[算法学习 1] 基因组组装算法 De Bruijn Graph - 知乎 (zhihu.com)
基因组大小估计教程|计算生物学核心 (uconn.edu)
https://www.jianshu.com/p/56421071e8e8
网友评论