第2章配对序列比对2bak

1、第四节 全局序列比对(Global Sequence Alignment)在介绍全局比对之前，我们得首先说明的是对任何比对方法（包括全局比对，局部比对）主要包括两部分：第一部分就是计算最大相似得分数，也就是这两个序列比对的最优化计算过程。第二部分是在根据这个最大相似得分数，通过“回溯”方式确定两个序列的比对结果即第一个序列中的各字符（在DNA中为对应的碱基，在蛋白质中为对应的氨基酸残基）与第二个序列中哪个字符相匹配。通过这个结果可以确定相应保守区，可变区，这样据此可推断这两个序列是否有同源关系以及其同源性程度。这部分工作在比对中一般称为“回溯”（traceback）。我们首先考虑两个序列 和

2、，现在要求解的问题是如何通过计算求出这两个序列的比对结果，以确定其相似程度。为使读者对比较容易理解全局比对的算法，我们这里通过一个具体例子的计算来说明。由于蛋白质有20个氨基酸残基，相应的得分表比较复杂，为此我们这里的例子采用一小段DNA序列来说明其具本计算过程。实例，两个序列S=ACCTGA 和T=CATGTAT，对应的得分矩阵假设为： 表2-3 四个碱基的得分矩阵ATCG*A 2-1-1-1-1T-1 2-1-1-1C-1-1 2-1-1G-1-1-1 2-1*-1-1-1-1 2 为此我们首先假设代表序列S中第a个字符即与序列T中第b个字符之间的相似分数，比如当a=3时，; b=4时，根

3、据得分矩阵，;同样的。这样我们可得到与序列S和T对应的表格即图2-7：i j01234T序列 5CATGT00-1-1-1-1-11A-1-12-1-1-12C-12-1-1-1-13G-1-1-1-12-14C-12-1-1-1-15T-1-1-12-126G-1-1-1-12-1S序列图2-20 序列S=ACCTGA 和T=CATGTAT之间的各碱基之间对应得分即得到对应的矩阵。然后我们根据比对的两大组成，首先就可在此基础上应用动态规划法计算各子序列之间最优的比对得分即相应的最佳相似性分数即比对的第一部分。然后应用回溯法给出比对结果。第一步，计算两个序列的最优（大）相似性分数。1，我们假设

4、是序列和 (0 i n, 0 j m)。m,和n分别为这两个序列的长度，在本例子中n=5;m=6。2，初始条件： 3，迭代，相应的迭代公式：对1 i n, 1 j m ，有这里的迭代公式对生物学工作者而言，可能比较复杂，为此，我们将它表示成相应的图：C BA图2-21 迭代过程示意图这样整个比对过程就可如下图所示： （0，0）（m,n）图2-22 整个动态规划计算示意图 4，寻找第m列（即最后一列）中所有元素的最大值: 5，寻找第n行（即最后一行）中所有元素的最大值：6，寻找最优比对值，即 现在我们以图2-20为例为说明上面公式的应用。根据初始条件，图2-20变为图2-23i j01234T序

5、列5CATGT00-1-2-3-4-51A-1-12-1-1-12C-22-1-1-1-13G-3-1-1-12-14C-42-1-1-1-15T-5-1-12-126G-6-1-1-12-1S序列图2-23 序列S=ACCTGA 和T=CATGTAT比对初始值设置图2-23中红色部分即是根据初条件计算得到（其它单元格中没有变）在迭代过程中，根据迭代公式计算V（1，1）（此时i=1,j=1）：据此可以得知的前一个单元是。 然后计算，此时i=1,j=2，即：据此可推断出的前一个单元为。此时图2-24变为（其它单元格中没有变）图2-24i j01234T序列5CATGT00-1-2-3-4-51A

6、-1-11-1-1-12C-22-1-1-1-13G-3-1-1-12-14C-42-1-1-1-15T-5-1-12-126G-6-1-1-12-1S序列图2-24 图2-23通过两步迭代计算后的结果不断应用上述迭代公式计算，便可得到全部的 即图2-26i j01234T序列5CATGT00-1-2-3-4-51A-1-110-122C-2100-1-23G-300-1214C-4-1-1-1115T-5-2-21036G-6-3-3032S序列图2-25 序列S=ACCTGA 和T=CATGTAT比对计算最终结果i j01234T序列5CATGT00-1-2-3-4-51A-1-110-1

7、22C-2100-1-23G-300-1214C-4-1-1-1115T-5-2-21036G-6-3-3032S序列图2-26 回溯图当两个序列每个碱基对的最优比对计算完成后，接下来的工作是在表2-7-3中找出最后一行和最后一列中最大的单元格，并比较这两个最大的单元格，以最终确定其最优比对结果。在这个例子中，最后一行的最优值（即最大值）是3，最后一列中最优值也为3。相应的单位格分别是（6，4）和（5，5）即图2-26中蓝色的两个方框。然后比较这两个方框中的值以确定最终比对结果（即S和T的Alignment）。在本例子中，由于它们相等，因此有两个比对结果，即有两个最优值。以上是全局比对算法的第

8、一步。第二步，回溯，给出两个序列的比对结果根据前面的讨论：当两个序列的最大相似得分数确定后，接下来的工作就是将第一个序列中的各字符与第二个序列中的各字符一一匹配起来。就是通常的回溯过程。由于我们在前面计算每对（第一个序列中的字符与第二个序列中的字符）字符时已记下它们对应的前一段字符，因此这个工作就相对比较容易。图2-26最终的回溯结果可以得到三个比对结果，见图2-27图2-27 三个比对结果得到上述三个比对结果，这里就涉及到生物大分子序列分析中的一个非常基本的概念，序列残基数相等（Sequence Identity）。它有很多方式的定义，但相比之下，我个人认为下面这个定义比较能反映比对的实际情

9、况：它的定义是通过比对，其残基数相同的配对数与全部配对数的百分比。比如在图2-8中，结果共有对残基相匹配即A-A，C-T，G-G，T-T。在这四对中有三对是相等的残基对，因此其序列残基数相等百分比为： ID%=3/4×100%=75% 这里应重点强调，序列残基数相等百分比的概念在生物信息学中尤其是生物大分子序列分析中非常有用，因此也非常重要。比如，现在在生物学界普遍认为当两个蛋白质通过比对后，如果序列残基数相等百分比大于20-30%，则就认为这两个蛋白质具有同源性，其三维构构象也比较相似。因而据此也可根据已知功能的蛋白质推断未知功能的新蛋白质。这正是生物信息学的主要宗旨之一：即通过计

10、算即生物信息学方法的处理，向实验工作者提供相应的有意义的结果，这里就是应用序列残基数相等百分比来推断新发现蛋白质的功能。显然有着非常现实的意义。 全局比对的计算时间是,(m,n分别为两个待比较序列的长度)级。显然如果两个较长的序列，其比对时间则会很快增加。近似算法FASTA和BLAST就是在尽可能不降低比对精度的情况下尽可能减少计算时间，以期增加比对速度。我们之所以用这么长的篇幅来介绍全局比对，主要原因有：1.全局比对是所有比对法的基础，如果全局比对方法搞清楚了，则其它比对方法就容易理解了。因为它们主要部分与全局比对方法是相同的。 2. 全局比对法涉及到的基本概念如序列残基数相等（Sequen

11、ce Identity），比对得分数等基本概念贯穿于整个生物大分子序列分析。对生物学工作理解应用某种软件分析生物大分子序列所得到的结果是必要的。如果这些基本概念没有搞清楚，则就谈不上应用生物信息学的基本方法来解决所面临的实际问题。第五节 局部序列比对（Local Alignment）顾名思义，局部比对就是要从两个序列中各自找出两个相似性最高的子序列片段。在许多情况下，如果从整体上看，两个序列并不怎么相似或者说它们通过全局比对得到的最大相似性得分数不大，但它们中的一些局部片段是很相似的，这就要求我们只要将这部分局分的片段找出来就达到目的，而不需要考查它们的整体相似性。这在现代生物学中是很有意义的

12、。比如，在蛋白质序列中，有一部分序列片段主要是与它的生物学功能有关，而另一部分片段则对该蛋白质的功能影响较小。因此，对两个类似这样的蛋白质，我们没有必要考虑它们整个蛋白质的比对而确定它们的相似性，我们只需要知道它们部分片段是否相似就可以了。这样就可以更加有针对性地比较两个蛋白质序列或DNA序列。典型的在DNA序列中，它有代码区（外显子 exons）和非代码区（内含子 introns）. 一般地，代码区能最终翻译成蛋白质，代表基因的功能，因此它不易突变，其主要目的是能保证该基因（或蛋白质）功能的稳定性，而非代码区即内含子则比较容易突变，因此如果应用局部比对技术比较两个不同种类的DNA，则找到的序

13、列最相似的那段片段就有可能是基因了。同样的，在蛋白质中，功能区与其它部分相对比较保守，因此应用局部比对技术比较两个不同来源的蛋白质，进而据此推断它们的功能。 对局部比对，我们在前面已有一个实例了。这里就不再重复了。局部比对的计算与全局比对相比，其主要不同点是：应用全局比对的计算方法进行计算，当计算到有分数值小于或等于零时，就将此分数值强制设定为零。而当回溯时碰到某个单元可为零时，就停止回溯，显然其比对结果中的子序列不一定要从原比对序列开始，也就是说可以是两个序列中的中间某个片段。这正体现了局部比对的精神。它的基本算法如下： 对两个待比对的序列S和T（与全局比对相同）1，初始条件：2，迭代，相应

14、的迭代公式：对1 i n, 1 j m ，有 3. 从V中寻找最大相似性得分即对应的对元格。从以止的计算公式中我们不难看出，局部比对计算与全局比对相比的不同点：1. 最初条件不一样，局部比对均为0，而局部比对则为字符与空位得分的累计。2. 在比对过程中，如果计算某个单元格小于零，就设置为零。3. 回溯过程中遇到零的单元格时就此停止。 我们仍以介绍全局比对法所用的例子来说明局部比对，主要是通过这个例子来说明局部比对与全局比对之间有什么不同。根据初始条件2-4和2-5，我们可得到如下的矩阵即图2-28.然后应用迭代公式2-6（A，B，C，D）计算相应的局部序列的得分总数，图2-28. 最后通过回溯

15、得到其比对结果即图2-29.由图2-29可知，所得到的比对结果中的两个序列均为原序列中的片段。这也就是前面多次提到的局部比对与全局比对的主要区别，也正是局部比对的意义之所在。 已有人证明：局部比对计算所需的时间与内存与全局比对基本相同。 局部比对计算因为对初始条件设置了0，加上在迭代过程中也设置了零，因此它所得到的比对结果可以在两个待比较的序列中任何地方开始，也可以在任何地方结束。这是局部比对区别于全局比对的最大特点。i j01234T序列 5CATGT00-1-1-1-1-11A-1-12-1-1-12C-12-1-1-1-13G-1-1-1-12-14C-12-1-1-1-15T-1-1-

16、12-126G-1-1-1-12-1S序列图2-28 局部序列比对中的值i j01234T序列 5CATGT00000001A0020002C0201003G0010304C0200025T0012006G000040S序列图2-29 局部序列比对的迭代结果然后根据回溯方法，我们得到图2-29的比对结果，将它写成相应的比对形式，得到如图2-30的局部比对结果（与全局比对结果的比较）图2-30 局部比对与全局比对的结果第六节 末端空位自由比对（End-Space Free Alignment）末端空位自由比对是比对的另一种变换形式。我们已知道，在全局比对中，其比对的开始和结束要在相同的地方，它允

17、许两个序列有所重叠。与前面介绍全局，局部比对一样，这里我们再用一个例子来说明：原始序列：S= C A C T G T A CT= G A C A C T T G与前面所讨论的一样，其四个碱基的得分矩阵与表2-1相同，这样我分得到其全局比对与末端空位自由比对结果如下： 图2-30 全局比对与末端空位比对结果比较从图2-30可以看出，末端空位自由比对两个序列的两端对应都有空位。那么末端空位自由比对的提出和建立到底有什么用呢。一个典型的例子就是在全基因组序列测定中“鸟枪序列聚集”。在实际生物学中，人们往往要从一个未知的DNA序列中获得许多序列碎片，这些碎片由于断列位置的差异和多样性，倒致它们之间的某

18、个片段相同，而两端不同，如图2-31：图2-31 鸟枪序列聚集示意图解决这类问题的就是如何将一个碎片序列与另一个碎片序列“粘合”，以达到重建原始序列的目的。由于这两个序列来源于一个原始序列的不同部分，而且，这两个重叠的子序列在原序列中不可能有相同的初始位置，也不可能有相同的末端位置。因此它们总的相似性是较底的。但在另一方面，如果我们应用“末端空位自由比对”法来比较这两个重叠的子序列，会得到很高的相似性分数。主要原因是它们重叠的部分是相同的。因此据此就可以将这些碎片（即子序列）“粘合”在一起，从而构建原始的DNA序列。这就是末端空位自由比对法的最大特点。接下来我们介绍这个方法的具体算法。 末端空

19、位自由比对法的算法是全局比对法算法的一个变形，而且与局部分比对在初始条件设置方面也相同。具体如下：1，初始条件：2，迭代，相应的迭代公式：对1 i n, 1 j m ，有3，寻找第m列（即最后一列）中所有元素的最大值： 4，寻找第n行（即最后一行）中所有元素的最大值：5，寻找最优比对值，即 以上就是末端空位自由比对算法。如果将上述计算应用下图中两个序列，则可得到如图2-11所示的结果：i j0A1C2T3G4T5T6A0 G00000001 T0-1212212 T0-1103433 A0-1002364 C02101255 T0-1432346 G0-1365437 A0-125876最佳值

20、 8图2-32 末端空位自由比对结果末端序列比对的时间计算量与空间复杂性同全局比对，局部比对相同。以上我们连续介绍了三种最常用的比对方法即全局比对，局部比对和末端空位自由比对。它们的基本算法几乎相同，但结果则不同，而且适用对象也不同：全局比对主要考虑两个序列总体的相似性；局部比对则是要找到两个序列中最相似的片段，比对蛋白质的功能区域；末端空位自由比对则是适用于两个序列有重叠部分，主要适用于得用原始序列的碎片重建原始序列。所以，尽管这些方法均建立在数学最优化法中的动态规划法基础上，而且算法基本相同，但在生物学中的适用对象是完全不同的。第七节 空段罚分比对（Gap Penalty）到目前为止，我们

21、在介绍全局比对，局部比对和末端自由空位比对时，分别来源于两个序列中的两个字符（如DNA序列中的碱基，蛋白质中的氨基酸残基）在比对时，可分为三类：匹配（match），不匹配（mismatch）和空位（Spaces）。在处理空位时，每个空位所给的得分数均相同。在这里我们再引进另一个概念：空段（gap）。它由连续几个空位组成，如下图：图2-33 空位与空段之间的关系空段与空位的不同在于空段将所含所有空位一起处理，而不是每个空位的得分相同。那么引进空段有什么生物学意义呢？在分子生物学中，DNA一次突变往往会引起多个碱基丢失或增加，而且每次的突变事件所丢失或增加的碱基数经常是不同的。因此一个空段是一个突

22、变事件的反映，因此将它作为一个整体来处理应该是比较合适的。引进“空段”概念在生物序列分析中的一个突出的作用就是DNA转录成RNA。我们知道，DNA序列主要包括外显子和内含子。当DNA通过转录（即A变为U，T变为A，C变为G，G变为C）变为转录RNA(即前信使RNA，前mRNA)时，DNA中包括内含子和外显子全部转录成mRNA，接下来便是“融合”过程。融合过程首先是确定外显子和内含子的边界，然后将RNA中内含子的部分除去，并将剩下的拼接连成mRNA，显然此时mRNA仅包含DNA中的外显子。生物体中的每个细胞一般都有含有该生物体的所有染色体，但对每个具体组织，器官的细胞如肝细胞，脑细胞，仅有其中的

23、一小部分基因被表达成最后的蛋白质。在分子生物学中，确定哪个蛋白质能在哪个器官，组织中被表达的一个标准方法是：首先找到与该蛋白质相对应的DNA即互补DNA（complementary DNA）,简称cDNA，cDNA所包含的碱基序列与原始的DNA中对应的能最终翻译成该蛋白质的外显子的序列是一致的，而且不包含其内含子，为了确定DNA中哪个基因能表达成该蛋白质，则需要比较cDNA与原始的DNA。 在这种情况下，如果采用一个很好的空段罚分模型，则可以阻止在比对时得到很低的得分，从而可以找到真正的基因。上述过程可以用图表示如下：图2-34 确定某个蛋白质的基因表达此外，减数分裂中染色体的不对称交换会产生

24、碱基的插入与删除；基因上某个单点突变；DNA在复制过程中融合时会产生一些复制品；逆转录病毒酶的插入；两个染色体之间的DNA转换。空段罚分方式主要有两种：一种是空段中每个空位的罚分均相同；另一种仿射空段罚分（Affine Gap Penalty），其数学表达式是上式中，为空段打开罚分，为空段延伸罚分，q为空段中的空位数。在介绍空段罚分模式的算法之间，我们先介绍如下符号： 类型1比对的最大得分值 类型2比对的最大得分值 类型3比对的最大得分值比对的最大得分值有了上述的符号，我们可得如下的罚分比对算法：初始条件：迭代过程：这里：其它的过程如回溯与前面所说的相同。这里不再重复。 到目前为止，我们比较详

25、细介绍了四种序列比对方法即全局比对方法、局部比对方法、末端空位自由比对方法和空段罚分比对方法。应该说，这些方法的思路基本相同，但在具体算法方面有或多或少的区别，现我们将它们总结如表2-4，读者可以从中比较。七总结 以上我们详细地介绍了生物大分子序列分析的几种精确比对方法。这几个方法非常重要，这也是我们详细介绍的原因。可以这样说，它们是整个现代生物信息学的基础。没有序列比对分析方法，生物信息学就只剩下蛋白质天然构象的从头预测算法，生物大分子模拟等内容。人类基因组计划所得到的基因组序列的分析也就无从下手。另外，我们从以上内容发现，生物信息学的基本概念都与生物学，分子生物学紧密相联，比如空位，空段等

26、均可在生物学中找到活生生的例子。此外，以上的几种比对方法均由其各自的应用对象，或者说，这些方法之所以建立，是因为生物学这门学科的实际需要。有关这方面我们在介绍各方法时均也有比较详细的介绍。表2-4 四种比对方法的比较初始条件迭代过程最大值的选取全局比对局部比对末端空位自由比对空段罚分比对可以是：也可以是：第八节 生物大分子序列分析的统计显著性从前面几种序列分析方法我们知道：当比对得分越高，则这两个序列越相似，表现在蛋白质上就是这两个蛋白质的同源性有可能越接近。但在实际情况中，由于不同的蛋白质序列长度是一不一样的。一般地，对序列较长的两对蛋白质或DNA，即使其实际相似性不是很高，其得分数也有可能

27、较高。因此仅凭得分数来推断两个待比较序列的相似性程度是不太符合实际情况的。比如对两段肽段S=AGPQSTM T=AGPQST，如果应用BLOSUM62得分矩阵及局部比对技术进行比对（其空段打开罚分为12，空位延伸为1），其最优比对的得分数为31，从这两个肽段来看，它们仅相差一个氨基酸残基，因此可以说它们实际上是很相似的，但我们用两个随机产生的序列，其长度分别为600个氨基酸残基，其最优比对的得分为58。如果我们仅用最优得分来推断其序列相似性，显然是两对随机序列的相似性要比上面两个肽段高。这显然与实际情况是不相符的。这就是我们通常所说的假阳性（false positive）。因此需要重新定义一种

28、新的得分标准来确定两个序列的相似程度。这里要介绍的序列分析统计显著性就是针对这种情况而提出的。但它也不是能够解决所有实际问题。当真正的比对的得分与随机序列比对的得分差异显著时，则采用统计显著性就能得到较好的结果；反之，当二者的分布有重叠时，则据此不易找到实际上真正相关的比对。图2-34与2-15就是两种典型的例子。所以，应用统计分析序列比对显著性也还需要看实际情况。图2-14说明有些序列统计显著性得到的结果一般为真阳性，而图2-15则说明应用统计显著性得到的结果可能会是假阳性。我们这里首先介绍一下序列比对统计显著性的基本研究情况。一开始，人们用一定数量的随机序列对进行比对，将相应的得分制成统计

29、分面图。对这个图，初看上去，象是正态分布，但仔细研究一下，则发现它符合Gumbel的极值分布。我们知道，在序列比对中有全局比对和局部比对，相比较而言，前者的统计显著性的建立要比较复杂，后者相对要比较简单。因为当两个序列应用全局比对方法进行比较时，既使是它们是随机序列或者这两个序列在实际情况中并不相关，但由于全局比对的计算思路是尽可能地多将两个序列中的字符（氨基酸残或碱基）排列起来，导致它们的得分仍会比较高。所以全局序列比对的统计显著性相对来说要显得比较复杂。一Z得分（Z-Score）最早引入的是Z-得分（Z-Score），以后Karlin, Altschul建立了局部比对的统计显著性。现在比较

30、通用的描述两个序列相似程度就是采用他们所建立的P值和E值。首先我们先介绍Z-得分。Z得分是比较老的一种，但在统计分析中非常有用，同样的，在比对分析中也非常有用。它的定义是离开随机序列比对得分的平均值的标准偏差数。如果某个比对得分为S，与这个比对相同序列长度的随机序列比对得平均值为u，对应的标准偏差为，则Z值的计算为： (2-19)在这个基础上，我们首先介绍一下全局比对显著性的基本性况。二、全局比对的统计显著性考虑到全局比对统计显著性的复杂性，许多人在这方面作了许多偿试性的工作。Abagyan建议在处理全局比对统计显著性时，采用值分面。但由于他在建立这个统计显著性时，采用的得分系统更象局部比对，

31、所以这个结论不被人们广泛采用。而Dayoff分别于1978年和1983年应用PAM250得分矩阵和Needleman-Wunsch方法，同时采用常数空段罚分方式，对一些不相关的实际蛋白质序列进行比对，发现其比对得分符合正态分布。在此基础上，他们建立了一个统计显著性方法，具体是：1. 对某对长度为m,n的两个序列A，B，先计算其相应的比对得分；2. 在保持与A，B的长度且氨基酸残基组成不变的情况下，产生100对或A不变，只变化B，或B不变只变化A的随机序列，对它们以同样的方式进行比对，得到相应的100个得分值，然后计算其平均值和标准偏差。然后根据上面Z值的公式计算相应的Z值。据此，Dayoff建

32、议有显著相关的序列的比对，其Z值一般为3-5，对应于相应的概率分别为和。但这个方有也有其潜在的问题：1.计算时间。因为这要做多个全局比对的计算，当在进行数据库搜索时，其计算时间可能要提高好几倍。但随着计算机性能提高，其这个问题不是很严重。2. 如果氨基酸残基的组成比较独特，或者其组成的复杂性较低，则这样的计算就显得过于简单。3. 当蛋白质序列相关性很高时，这样的序列模式不太适合于随机序列，而是适合于序列片段变化的序列。与此同时，Dayoff又设计了第二个统计显著性，具体做法是：先将一个序列所有可能的且长度一定的片段（segment）搜集起来，然后与另一个序列所有可能的相同的长度的片段一一比较，

33、然后计算这些得分的平均值和标准偏差，据此计算相应的Z值，根据Z值来判断序列A中的哪个片段与序列B中的哪个片段比较的得分是统计显著的。这个方法的优点是能找到两个序列中高度相似的片段，尤其是适用于保守区间小的活性位点的分析。三 随机DNA序列比对的统计显著性 首先将设计这样一个数学模型：在一个袋子内装满珠子，它们分别标记为A，T，C，G，换言之，在这个袋子内，共装有A，T，C，G的4种珠子。那么可以这样认为，如果取出10对相同类型的珠子（如A与A），则相当于在比对中的10个相等列。显然，取出一对是AA的珠子的概率,由于有4种可能性，所以取出相同类型的珠子的概率为。显然取出6个相等对的概率为。对应的

34、，取出一对不相匹配的概率为，连续取出6个这样的概率为。显然随机序列的比对则是这两种类型的混合体（匹配与不匹配）。对出现机率不相等的情况则稍有不同，但出入不会很大：取出相等对的概率为 上述模型严格来说是不太合适的。为此，人们引进了统计学中的“抛硬币试验”。我们知道，在序列比对中主要出现两种情况：一种是匹配，一种是不匹配，这类似于硬币的正面和反面。如果出现正面的概率为0.5，则它的对数的基为1/0.5=2。对抛100次硬币，连续为正面的最多的次数为。如果将序列比对的模式与抛硬币试验相结合的话，则有：图2-37 掷硬币试验与DNA序列比对这样，结合前面讨论，对长度分别为m,n的两个序列的比对，其连续

35、为匹配数的区间的长度为：比如，对两个长度为100，200的DNA序列比对，其R值为：为使结果更精确，人们建立了最长的匹配数数学期望和标准偏差：这里，.事实上，数学期望值也可简化为：这里有必要提一下序列校正的问题。比如对两个序列长度为100的比对，其为6.65,显然，当序列到93个时，其后面的字符不可能形成一个约有7个的连续的匹配对，所以其真实序列长度为100-7=93。这就是我们通常所说的序列校正应。值得一提的是，Altschul于1990年对值又作了改进，引入了一个常数，其值为,这样，值就变为根据前面的讨论，我们知道是两个相同字符匹配对的概率，在DNA中为0.25。可以根据得分矩阵中的平均概

36、率来计算。而相比之下，更具有普遍性。有关与的计算我们将在下面作介绍。二局部比对的统计显著性-极值分布（Extreme Value Distribution）对一个局部比对，其得分矩阵（或得分系统）为sij,背景概率为pi,则一个比对得分的希望值为 E=i,j=1m,npipjsij在给定一个比对得分，应用局部比对得到一个较高的得分，通常希望这个得分是一个偶然事件。人们通过进步研究发现，局部比对的得分符合统计学中的极值分布。即其概率密度函数为：（s）=EXP(-Z-EXP(-Z) (2-19)则对应的概率分布函数为： P（z>Z）=1-EXP(-EXP(-Z) (2-20)也可写成： P（

37、S<x）=EXP(-e-(x-) （2-21）相应的概率分布图见图2-16： 该图是=0, =1的标准极值分布图， 类似于标准准正态分布图。从图2-16可知，极值分布是一个不对称的分布，峰值的左边较右边陡。对一定的得分矩阵sij和组成已知的两个序列（各组成元素组成为p1,p2,pr, 对蛋白质序列r=20；对DNA序列r=4），我们希望其随机序列比对得分的数学期望即i,j=1rpipjsij为负值。其含义是当比对得分为正时，就说明这两个序列存在一定的相似（至少不是两个随机序列排列）。 在式2-21中，它有两个参数，一个是，它为特征化值，也是分布的中心，代表分布的衰减常数。图2-16就是=

38、0，=1的标准极值分布。其中是如下方程式的正解：i,j=1rpipjEXP(sijx)=1 (2-22)依赖于比对序列的长度，可由下式计算：=ln(Kmn)/ (2-23)这样我们将式（2-21）中的消除得：P（Sx）=1-EXP(-Kmne-x) (2-24)显然，对两个长度分别为m,n的两个序列的比对，如果知道了K和，就可以对任意比对得分为s的计算出相应的概率P（Ss）。问题是如何估算K和。已有许多方法对此作了深入的探讨，其中比较简单的方法是根据一定的背景概率pi应用计算机生成一定数量如10000对的随机序列，然后计算相应的比对得分，据此计算出相应的K与。这里比较详细介绍一下其计算方法。以

39、BLOSUM62作为得分矩阵，其罚分函数为（-12，-1），我们首先计算出10000对的长度均为600（即m=n=600）的得分值，据此可得到其概率密度分布,见图2-17根据这个概率密度分布图，我们可计算出其平均值u和标准偏差分别为36.449和4.427. 此外，我们可将Z写成Z=（S-u）/=S-ln(Kmn)=(S-ln(Kmn)/)/(1/)很明显， =1/ (2-25)ln(Kmn)/=u (2-26)根据2-25和2-26，我们可得到对应的=0.226 K=0.010为使读者能自己体验计算K和，我们这里给出各个氨基酸残基的背景概率：表2-2 20个氨基酸残基的频率表（% 总共氨基酸残基数为450431,