考虑到cnblog不适合基因组领域这种类型的文章,进过多番折腾,终于用jekyll+github搭了个独立博http://www.stbioinf.com/,如今博客已经搬迁至这里,很是欢迎你们订阅! html
考虑这样一个问题,“若是要保证基因组上95%的区域其覆盖深度在30x以上的话,那么最低的平均测序深度应该是多少?”。git
关于测序量的估计,对于作生物信息的人来说应算是屡见不鲜了,多数时候咱们都能直接根据以往项目的经验来得到,或是说的更具体些,在变异检测中通常要有25x以上的覆盖度才能获得一个比较靠谱的结果,因而以此为目的给出测序量的估计值;固然少数状况下也会有直接拍脑壳拍出一个值来的疯狂行为,不过嘛,虽然说是拍脑壳,但也不是随便拍的,拍脑壳的背景靠的是身后丰富的经验。相对更好一些的估计方式就是直接模拟数据,不过老是用模拟数据仍是让人以为麻烦,最好是能不用花多少时间,也不用作不少的计算就能脱口给出。我想在这里说一下这种状况下个人解法。固然了并不必定彻底准确,仅做交流,欢迎各位看官拍砖。github
闲话说完,回到上面的问题,在不经过数据模拟也不“拍脑壳”的状况下,要如何才能估算出一个合理的值呢?其实在做出任何推断以前咱们都应当要先有一个合理的前提假设,或者说是理论依据来做为后续分析的基础。咱们都知道短序列测序的一个特色是,在理论状况下位点被覆盖到的深度符合泊松分布(测序没什么问题的话,实际的情形也相差很少),但实际上在这种状况下用正态分布来考虑也是合理的,做为一个估计值,偏差也是可以接受的,这是咱们的基础。之因此想用正态分布来考虑,是由于正态分布有许多方便于计算的性质。其中一个颇有用的法则,就是“68-95-99”法则,意思就是距离均值一个标准差的区域围起来的面积大约是整体的68%,2个标准差的区域范围的面积是整体的95%,3个标准差区域范围占到了整体的99%,若是你本身想要验证这一法则也并不困难,只需作些积分就能算出来,但这里就不作计算了。以下图,均值用μ表示,标准差用σ表示。htm
如今事情就很简单了,从图中咱们能够看出,只要30x深度的位置在-2σ如下,那么就能达到理论的要求。要获得这一结果,问题就只剩下一个了,此时咱们只须要知道测序深度分布的标准差就能粗略估计出此时咱们所须要的最低平均测序深度。虽然这个标准差跟许多因素有关,这里以illumina公司的Hiseq系列测序仪为例子,依照以往基因组重测序的经验,σ约等于10x。那么,简单算一下,此刻,理论上咱们只须要测50x就可使得基因组上有97.7%的区域其覆盖深度在30x以上了,注意这里不是95%了,由于咱们的区域其实是[-2σ, +∞),而不是[-2σ,+2σ]! 再除掉一些边边角角的偏差,50x这个值在这里应当是合理的了。blog
以上计算都是以正态分布为基础而作出的估计。固然了,若是必定要用泊松分布去推算也能够,只是运算起来会麻烦不少。此外,若是是不一样系列或是不一样公司的测序仪,σ就不必定是10了。get