11)利用Reads的序列。
12)将Reads进行比对,调整Contig分值。
由于Shotgun方法产生的数据大部分是小段的序列(500bp目标序列可能包含30 kbp,则有500-2000,每段的平均长度为m(n/m,2-W算法的每个比对的时间为O),则总的比对时间为O);而剩下的工作,大部分工作是在形成的Pairs(n2在Phrap-W算法,而一些LLR下面内容不包含在程序中:
多序列比对:
在一些情况下,是需要进行多序列比对的工作,比如要比较多个序列之间的相似性来构造系统发育树,也就是为了更好的理解蛋白质或DNA的方法来比对多序列,设有k,则使用类似于S算法的公式:F(i,j) =max{ F(i-1,j-1)+s(xi,yj), F(i-1, j-1)-d, F(i-1, j-1)–d,0}来计算,则每个元素的计算来自2k-1个元素,总共有nk个元素需要计算,则计算的复杂度为O(2k*nk),而占用的空间在没有优化的情况下是O(nk),设单位计算量为1ns,则6个长度为100的序列进行比对,需要的时间为26*1006/10-9=6.4*104(s)个序列,则时间为2.56*109(s)=81(yrs),实际上是没法用于实用中的。而使用的方法往往是一些带有启发信息的算法。
例如下图的比对显示了三个序列的最佳比对:
<![endif]--> 
而使用标准的动态规划的方法使用的每个节点的计算量为2k-1,参考下图:
<![endif]-->
1)利用Pairwise比对构造树来构造多序列比对:
设S={s1,s2,..,sm}为序列集,则进行Pairwise的序列比对,从S中删除得分值最高的两个序列,然后将它们的比对放入S中。则循环直到S只剩下一个序列为止。
则设为第k个循环,k=0..m-2,共有m-k个序列,则要计算的比对数为C(m, 2)(k=0)和m-k-1(只需计算新的序列和原有序列的比对分值),则总共需要计算的分值次数为:
C(m, 2)+(m-2)+(m-3)+…+1 = m*(m-1)/2 + (m-1)*(m-2)/2 = (m-1)2,而计算一个分值的时间是O(n2)。则总时间为O(n2*m2),则远小于O(nm)的普通方法。而且如果是做一个序列和数据库内序列进行比较,则数据库中可以保存已经存在的序列之间的比对值,则总的计算时间最坏情况是:每次产生了新的比对都不存在于数据库中,也就是每次计算中,新的序列都参与了新比对的生成,则总的计算次数为(m –1) + (m-2) + … + 1 = m*(m-1)/2次比对之间的比较。而最好的情况是新的比对直到最终时才参与新的比对的生成,则总计算次数为(n–1) + 1 + … + 1 = 2*n – 3次比对计算。
示意图为:
<![endif]--> 
则完成时的树形结构也可以看作是系统发育树,也就是从最初的分子(最顶层)是怎么逐步演化到现在的分子的(最低层)。
2)启发式算法:

假设我们得到了启发式的比对结果,然后又得知该结果和最优解十分接近,如上图,在启发式结果比对为中心的一个圆柱内,则搜索最优解的空间缩小了很多。而S算法中Banded之间的匹配来定义的,而搜索圆柱的宽度则是由Bandwidth序列比对的启发式算法:
:和下面的BLAST中的有些类似,也是先比较数据库中的String来寻找相似度的。设长度为n的两个序列s,则计算k-tuple长度的Word=1或2和t[j]相同,则两者之间的偏移为i-j到n-1,设置一个偏移量表,范围是-n+1->m-1,则检查s,并建立一张表存储着每个k-tuple中的k-tuple与s中相同的k-tuple之间的偏移量为多大,则在前面的偏移量表中相对应的值所在的单元增加一个,如图所示,设k=1中A的位置为2、7,则t和5所在位置的A和-3/1/2:
),然后以这个偏移量作为中心做比对,首先整合距离比较接近的k-tuples为一个区域(region)。然后将序列和整个数据库中所有序列的原始分进行比较,得到分布图和平均值,再依照原始分排序,计算排序靠前的数据库中序列和当前序列的比对值,得到优化分(optimal score)。算法中的k :
算法在比较相似性时先搜索当前序列和数据库中的序列之间的相似片断(HSPs然后再计算每个匹配的统计显著性,最后返回那些满足用户选择的显著性域值的匹配。计算HSP的Word的Word匹配上或者满足一定的正值的域值,然后将Word和FastA算法没有仔细看过,其中的内容来自一些国外的课程网页。
序列和结构
序列和结构这两大类不同性质的数据在数据量方面有天壤之别。对这一点必须有个明确的概念。蛋白质三维结构数据库显然难以与序列数据库的数据量相比,这是因为结构数据的采集、存储与处理远比序列数据复杂。从信息理论角度看,结构数据与序列数据之间数据量的巨大差异,反映了这两类既不相同、却又相关的数据之间信息量的差异。随着基因组计划的实施,序列数据大量积累,这种差距会越来越大。当然,结构数据也在快速增长。当然,这与序列数据每年翻番的增长速度相比,依然不可同日而语。
从目前的研究状况来看,下图分子生物信息数据库的构造过程中的上3/4部分产生的数据相对容易。但是作为生物大分子来说,只知道序列是不够的,蛋白质的功能往往是和分子的三维结构有关的,而且如果由于外界环境的改变导致蛋白质的三维结构展开,则蛋白质产生了变性,则再不能完成以前所承担的生物活性。蛋白质折叠问题所要解决的就是蛋白质一级结构中的氨基酸序列最终怎样折叠成三维空间结构。1973通过实验发现,变性的核糖核酸酶(ribonuclease)
一级结构: 二级结构: 螺旋,β转角,β超二级结构:由相邻的二级结构单元组合而成的松散结构体,可作为蛋白质三维结构的构件,如(βαβ折叠桶等);
三级结构:由二级结构和超二级结构缔结而成的完整蛋白质结构;
四级结构:一些独立的蛋白质经非共价键缔结而成的聚集体;
五级结构:由独立的生物大分子相互作用而成的聚集体,如蛋白质—核酸聚集体
在预测生物大分子的结构时,经典的而且也是数学家和物理学家所喜欢的约束只有一个:能量最小约束。蛋白质及其他大分子之所以呈现稳定的分子结构是因为具有该结构的大分子具有了最小的能量,分子本身不能释放能量而改变结构,除非吸收能量,因此从分子本身来说,这样的结构是稳定的。由于在三维空间中原子的自由度很高,原子之间有化合键及其他非化合键相关联,完全计算具有很高的时间和空间复杂度,而且计算时存在了很多的能量局部极小值,很容易陷入局部极小值中去。因此只凭着能量最小约束在目前是没有可能完成普通的蛋白质分子的结构预测,通常要增加一些约束或者一些启发信息,比如说是保守区域、二硫键的存在等。通常蛋白质二级结构预测的方法有三种。一是由已知结构统计各种氨基酸残基形成二级结构的构象趋势(其中最常用的是Chou 法);二是基于氨基酸的物理化学性质,包括堆积性(compactness)、电荷性、氢键形成能力等;三是通过序列比对,由已知三维结构的同源蛋白推断未知蛋白的二级结构。
虽然从每两年一届的蛋白质分子结构预测竞赛的结果来看,似乎预测的结果很鼓舞人心,但是很多的预测方法并不是普适的方法,都是另外加了约束或者和已知结构的蛋白质序列做比较后再预测,因此只能对某一类结构的蛋白质或者和已知蛋白质类似的蛋白质的预测结果较好。
根据与当前蛋白质结构库中的蛋白质的相似性或者同源性来进行结构的预测,要依赖与蛋白质结构库中的数据的多少。从代数上讲,就是蛋白质库中的数据是否包含了蛋白质结构空间的所有的基向量,而且这些基向量具体是哪些,一个新的蛋白质和这些基向量之间的关系系数是怎么确定的,以及是怎么由基向量来得到结构。而从概率上讲则是蛋白质库中的数据是否已经达到了足够大的样本量,能够对未知序列进行概率估计。如果确实是存在了足够的样本,则最后预测就变成了模式识别。Chou 法之所以现在准确性不高,也是因为用作统计二级结构构象趋势的蛋白质空间结构数据库中非同源蛋白的数量还不够多。而且由于结构测定速度的制约,这一数据库容量不足的问题将长期存在。
而蛋白质三级结构预测,特别是基于二级结构预测的三级结构预测,尽管已经由个别成功的例子,总的说来,还远远没有成熟。
模式识别和预测:
模式识别和预测是生物信息学中两种基本分析工具,虽然上文中有些说预测的内容实际上是模式识别,但是从它们所要解决的最基本的问题来看是存在着一定的区别。
模式识别的基本思想是利用存在于蛋白质序列或结构中的某些特征模式识别相关蛋白质的性质。如果某一蛋白质序列或结构中的一部分具有保守性,这种保守性或者与蛋白质的生物活性有关,或者与折叠方式有关;那么,这种特征模式可以用来识别该蛋白质家族中的新成员。也就是所如果将已知蛋白质的特征序列模式和特征结构模式搜集起来,构建成数据库,则可以用来确定新测定的蛋白质序列中是否具有某种特征模式,从而确定该未知蛋白属于哪个蛋白质家族。目前,利用序列模式和结构模板数据库查询确定蛋白质家族关系,从而推断该新序列的功能和结构,已经成了常用的方法。
显然,无论是序列模式识别,还是结构模式识别,都建立在已知序列和结构的基础上,这些已知序列和结构存放在一定的数据库中。应该说,序列模式识别比较容易,其结果也比较可靠。相比之下,结构识别亦即折叠模式识别要困难得多,往往需要有专门研究人员参与。即使如此,其准确性也只能达到40%相反,预测是生物信息学中的棘手问题,目前尚无行之有效的方法,预计在未来十年内也很难取得关键性突破。所谓预测,是指直接从氨基酸序列推断某一蛋白质的功能位点或预测其三维结构,它并不依赖于已知蛋白。因此,预测方法不需要建立序列模式或结构模式数据库,而需要研究开发解决蛋白质折叠问题的方法和软件。
序列分析的主要难点是如何把生物学信息和序列数据联系起来。如何从未知序列数据库搜索所得同源序列中提取功能信息,并非一件易事。由计算机程序自动生成的注释,其中不少是错误的;而这些错误结果已经整合到序列数据库中,并随之而扩散。目前,尚未找到对注释进行质量控制的有效方法。在信息爆炸的今天,计算机的使用无疑是十分必要的;然而,完全依赖计算机对数据进行自动处理,具有潜在的危险。
当序列相似性只局限于部分区域时,对分析结果的解释就更加复杂。这在蛋白质模块分析时尤其需要注意。所谓蛋白质模块,是指蛋白质序列中可形成独立折叠单元的连续片段。蛋白质模块可以看作是组成蛋白质结构域的基本单位,也是蛋白质的基本构造单位。从遗传学的观点看,这些模块在核酸序列中的分布,不是简单的基因复制和融合,而是由某种基因混排(gene shuffling)机制产生的。Jacob年认为不管具体过程如何,“自然”象一个巧裁缝,他把各式各样的补丁搜集起来,缝制成一件百衲衣。进化过程并不是完全从头开始,而是利用原有的材料,通过改造而产生新功能,或者是把几种不同的系统整合到一起,形成更好的新系统。
简单构件的再利用和整合,是更大、更复杂的系统中产生新的、未曾预料的新功能的关键。但是依照目前的水平,由已知的简单构件预测整体功能的正确率极为有限。根本原因在于一个系统的特征可以由组成该系统的部件的特征来解释,却不能由这些部件的特征进行推断,原因在于:(a) 已知的嵌合蛋白中的大部分单体都无法帮助我们预测缺失的模块(c一个经典Rube Goldberg联动机说明了预测复杂的生物系的难度。这种联动机的95%未知部件的功能仍然无法判断。生物系统就象这个联动机。Jacob 依照现在的分子生物学的发展看来,将逐渐从获取序列逐渐转向获取结构、获取功能、预测结构、预测功能、设计结构、设计功能方向发展,真正成为一门完全以计算为主的计算生物学学科。(汪冬)