第三章序列比较3.3序列多重比对与序列两两比对不一样，序列多重比对（MultipleAlignment）的目标是发现多条序列的共性。
如果说序列两两比对主要用于建立两条序列的同源关系和推测它们的结构、功能，那么，同时比对一组序列对于研究分子结构、功能及进化关系更为有用。
例如，某些在生物学上有重要意义的相似性只能通过将多个序列对比排列起来才能识别。
对于一系列同源蛋白质，人们希望研究隐含在蛋白质序列中的系统发育的关系，以便更好地理解这些蛋白质的进化。
序列两两比对往往不能满足这样的需要，难以发现多个序列的共性，必须同时比对多条同源序列。
图3.14是从多条免疫球蛋白序列中提取的8个片段的多重比对。
这8个片段的多重比对揭示了保守的残基（一个是来自于二硫桥的半胱氨酸，另一个是色氨酸）、保守区域（特别是前4个片段末端的Q-PG）和其他更复杂的模式，如1位和3位的疏水残基。
实际上，多重序列比对在蛋白质结构的预测中非常有用。
多重比对也能用来推测各个序列的进化历史。
从图3.14可以看出，前4条序列与后4条序列可能是从两个不同祖先演化而来，而这两个祖先又是由一个最原始的祖先演化得到。
实际上，其中的4个片段是从免疫球蛋白的可变区域取出的，而另4个片段则从免疫球蛋白的恒定区域取出。
当然，如果要详细研究进化关系，还必须取更长的序列进行比对分析。
对于多重序列比对的定义，实际上是两个序列的推广。
设有k个序列s1,s2,...,sk，每个序列由同一个字母表中的字符组成，k大于2；通过插入操作，使得各序列s1,s2,...,sk的长度一样，从而形成这些序列的多重比对。
如果将各序列在垂直方向排列起来，则可以根据每一列观察各序列中字符的对应关系，如图3.14。
通过序列的多重比对，可以得到一个序列家族的序列特征。
当给定一个新序列时，根据序列特征，可以判断这个序列是否属于该家族。
对于多序列比对，现有的大多数算法都基于渐进比对的思想，在序列两两比对的基础上逐步优化多序列比对的结果。
进行多序列比对后，可以对比对结果进行进一步处理，例如构建序列的特征模式，将序列聚类，构建分子进化树等。
3.3.1SP模型SP模型（Sum-of-Pairs，逐对加和）是一种多重序列比对的评价模型。
在多重比对中，首先要对所得到的比对进行评价，以确定其优劣。
例如，对图3.14中的8条序列进行比对，可以得到另外两种结果，如图3.15所示。
那么，这样的三个多重比对，哪一个更好呢？这就需要有一种方法来评价一个多重比对。
评价一个多重序列比对比评价序列两两比对结果更复杂。
这里，我们假设得分（代价）函数具有加和性，即多重比对的得分是各列得分总和。
那么，我们首先考虑如何给比对的每一列打分，然后将各列的和加起来，成为一个总得分。
在处理每一列时，自然的处理方式是寻找一个具有k个变量的打分函数（k是参与多重比对的序列的个数），而每一个变量或者是一个来自特定字母表中的字符，或者是一个空位。
我们很难得到这样一种具有k个变量的表达式函数。
另一方面，这种隐式函数不具有统一的形式，随着k的变化，函数的表现形式也发生变化，不利于计算机处理。
可以考虑使用显式函数，在具体实现显式函数时，用一个k维数组来表示该显式函数（类似于打分矩阵），指定对应于k个变量各种组合的函数值。
这带来一个新问题，即所需的数组空间很大，而且随着k的变化，数据结构也要随之动态变化。
我们所期望的函数在形式上应该简单，具有统一的形式，不随序列的个数而发生形式变化。
根据得分函数的意义，函数值应独立于各参数的顺序，即与待比较的序列先后次序无关。
满足上述条件的一个函数就是常用的逐对加和SP函数。
SP函数定义为一列中所有字符对得分之和：(3-34)其中，c1,c2,…,ck是一列中的k个字符，p是关于一对字符相似性的打分函数。
对于p可采用不同的定义。
下面一个是SP函数计算例子。
总得分根据字符两两比较得分计算演化而来，即逐对计算p(1,2)，p(1,3)，...，p(1,8)，p(2,3)，p(2,4)，...，p(2,8)，...，p(7,8)的所有得分，再加和得到结果。
在上述计算中，我们假定：如果两个对比的字符相同，则得分为0，否则，得分为-1。
所以，上述一列的得分为：（-7-6-5-4-3-2-1）+2=-26。
将一个多重比对所有列的得分全部加起来，其和即为该多重序列比对的得分。
注意：在进行多重序列比对时，可能会出现两个空位字符的比对，因此需要扩充函数p的定义域，即增加p(-，-)的定义。
通常的定义是：p(-,-)=（3-35）这似乎有点意外，因为一般只要是遇到空位字符，其得分就应该是负的，所以两个空位字符的比对应得到更多的负分。
但是，这样处理是有道理的。
在序列多重比对中，我们往往在得到整体比对的基础上进一步分析两条序列的对应关系。
例如，根据图3.14的比对结果，取出最后两条比对的序列，见图3.16（a）。
这里存在空位字符对比的列，相当于这两条序列都进行了插入操作。
但是由于插入位置相同[如图3.16（a）箭头所指位置]，这两条序列本身在此位置上是完全相同的，所以，此位置上的编辑代价为零，或者得分为0。
因而在分析这两条序列时，可以去掉这些空位字符对比的列，得到由多重序列比对结果推演的序列两两比对，如图3.16（b）所示。
其结果又称为多重序列比对在两条特定序列上的投影（projection）。
若先处理每一个序列对，而在处理序列对时，逐个计算字符对，最后加和，则SP得分模型的计算公式如下：(3-36)其中，α是一个多重比对，αi,j是由α推演出来的序列si和sj的两两比对。
具体计算SP时，我们可以先对多重序列比对a的每一列进行计算，然后将每一列的得分值相加；也可以先计算每一对推演出来的两两序列比对的得分，然后再加和。
但是注意，上述两种计算过程仅在p(-,-)=0这一条件成立下才等价。
3.3.2多重比对的动态规划算法多重序列比对的最终目标是通过处理得到一个得分最高（或代价最小）的序列对比排列，从而分析各序列之间的相似性和差异。
如同处理序列两两比对一样，依然可以用动态规划算法。
回想前一节中穿越得分矩阵的路径（见图3.9）。
得分矩阵相当于二维平面，而对于三条序列，每一种可能的比对可类似地用三维晶格中的一条路径表示，而每一维对应于一条序列，如图3.17所示。
路径的起点为晶格的左下角，路径的终点为晶格的右上后角。
图3.17中的路径记为(0,0,0)，(1,0,0)，(2,1,0)，(3,2,0)，(3,3,1)，(4,3,2)。
如果存在多条序列，则形成的空间是超晶格（Hyperlattice）。
假设在三维晶格的顶端、前面、后面各有一平行光源，这些光源将代表多重序列比对的路径投影到相对的侧面，即投影到一对序列所在的平面（图3.17仅画出了右边的光源）。
将多重比对投影到两个序列所在的平面在加速优化计算中将起着重要的作用。
一个多重序列比对投影的结果可能比原来要短，例如，下面的多重序列比对G---SNSGN----SGNAVSNS如果在前两条序列方向进行投影，则投影结果为：G-SNSGN--S如果路径的某一段沿投影方向进行，那么该段路径就不可能产生投影。
图3.17中的第一段沿着第一条序列向前进，推进方向垂直于其他两条序列，则平行光源对这一段不产生直线投影，而仅仅产生一个投影点。
图3.17三条序列所对应的三维晶格序列两两比对的动态规划算法经改进后可直接用于序列的多重比对。
就二重序列而言，将动态规划算法的计算过程看成是在二维平面上按一定顺序访问每一个节点，访问节点的先后顺序取决于节点之间的关系，如图3.7所示。
二维平面（或二维晶格）与前面所介绍的超晶格相似，上一节得分矩阵中的位置与每一个晶格节点相对应。
在计算过程中，对每个节点逐一计算其得分（或代价）。
每个节点的得分代表两个序列前缀的最优比对的得分，而晶格最后一个节点（右下角节点）的得分即为两条完整序列的比对得分。
在超晶格中，序列的放置与上一节有所不同，计算是从左下角进行的，而不是得分矩阵中的左上角。
如图3.18所示，计算过程从超晶格的坐标点（0,0,…,0）开始，按节点之间的依赖关系向右上后方推进，直到计算完最后一个节点。
实际计算时，可以采用类似二维的计算方式，即先低维后高维（对应于先行后列）。
在图3.18中，当前点的得分计算取决于与它相邻的7条边，分别对应于匹配、替换或引入空位等三种编辑操作。
计算各操作的得分（包括前趋节点的得分），选择一个得分最大的操作，并将得分和存放于该节点。
在三维或超晶格中，计算公式与上一节相似，但计算一个节点依赖于更多的前趋节点。
在三维情况下要考虑7个前趋节点，在k维情况下要考虑2k-1个前趋节点。
假设以k维数组A存放超晶格，则计算过程如下：①a[0,0,…,0]=0②a[i]=max{a[i-b]+SP-score[Column(s,i,b)]}这里i是一个向量，代表当前点，b是具有k个元素的非零二进向量，代表i与前一个点的相对位置差（例如，在二维的情况下，b=（1，1）、（1，0）或（0，1）），s代表待比对的序列集合，而其中，sj[ij]表示第j条序列在第ij位的字符，SP-score(Column(s,i,b))代表SP模型的得分值。
计算过程是一个递推的过程，在计算每个晶格节点得分的时候，将其各前趋节点的值分别加上从前趋节点到当前点的SP得分，然后，取最大值作为当前节点的值。
随着待比对的序列数目增加，计算量和所要求的计算空间猛增。
对于k条序列的比对，动态规划算法需要处理k维空间里的每一个节点，计算量自然与晶格中的节点数成正比，而节点数等于各序列长度的乘积。
另外，计算每个节点依赖于其前趋节点的个数。
在二维情况下，前趋节点个数等于3，三维等于7，四维等于15。
对于k维超晶格，前趋节点的个数等于2k-1。
这个计算量是巨大的，而动态规划算法对计算空间的要求也是很大的。
可以证明，利用SP模型寻找最优多重序列比对是一个NP-完全问题。
1971年，Cook给出NP-完全问题的定义。
P类问题为多项式界的问题；而NP-完全问题是这样一类问题，如果其中的某个问题存在多项式界的算法，则NP类中的每一个问题都存在一个多项式界算法。
目前，人们对NP-完全问题通常采用以下二种方法进行求解：舍去寻找最优解的要求，寻找对一般问题比较接近最优解的近似解；利用非常规的求解技术求解，例如，利用神经网络、遗传算法等方法进行问题求解。
对于生物信息学的问题，有许多是NP-完全问题。
在实际应用中，有一些变通的方法可以用来求解NP-完全问题。
下面分别介绍几种实用的多重序列比对算法。
3.3.3优化计算方法如果待比对的序列很多，多重序列比对所形成的超晶格空间将会非常大，算法不可能访问所有的节点，而且也没有必要这样做。
正如序列两两比对一样，可以想象代表最优的路径处于主对角线附近，至少不会在超晶格中的某个平面上。
假设已知一个试探性的路径与最优的路径靠近，并且两条路径最多相距r个单位，则应当将动态规划算法所要计算的节点限制在以试探性路径为中心、半径为r的区域内，这样做无疑可以提高算法的效率。
我们不能盲目地搜索整个空间，应该利用人工智能空间搜索策略的剪枝技术，根据问题本身的特殊性将搜索空间限定在一个较小的区域范围内。
若问题是搜索一条得分最高（或代价最小）的路径，那么，在搜索时，如果当前路径的得分低于某个下限（或累积代价已经超过某个上限），则对当前路径进行剪枝，即不再搜索当前路径的后续空间。
在序列两两比对中,Fickett和Ukkonen设计了一种称为定界约束的方法来缩小搜索空间、减少计算量，其中得分矩阵的上界和下界可以预先确定或动态变化。
为了在多维空间上使用动态规划算法，Carrillo和Lipman将这种思想引入到多重序列比对，即先进行初步的序列双重比对，以便限制进一步作多重序列全面比对所需要的多维空间，从而减少计算量。
这样可以克服多序列的维数、空间以及运算量之间的矛盾。
这里介绍Carrillo-Lipman的优化计算方法，该算法是基于下述关系的，即多重序列比对与向两个序列投影之间的关系，也就是通过公式（3-36）将SP得分与序列两两比对的得分联系起来。
设k条序列s1,s2,...,sk的长度分别为n1,n2,...,nk，按照SP得分模型计算这些序列的最优比对。
假设α是关于k条序列s1,s2,...,sk的最优多重比对。
注意：这并不意味着α向某两条序列的投影是这两条序列的最优比对。
在说明具体的优化计算方法之前，先介绍一种计算两条序列经过特定断点的最优比对算法。
设有两条序列s、t，已知它们的两个断点分别是i、j（即s和t分别在i和j处一分为二），则s、t对于经过特定断点（i、j）的最优比对可分为两个部分，一是对应于两条序列前缀0:s:i与0:t:j的最优比对，另一个是对应于两条序列后缀i:s:m与j:t:n的最优比对,m、n分别为两条序列的长度。
实际上，经过特定断点的最优比对是两个比对。
为了得到特定断点的最优比对，用两个矩阵A和B，其值为：ai,j=sim(0:s:i,0:t:j)bi,j=sim(i:s:m,j:t:n)对于矩阵A的计算与标准算法一样，而对于矩阵B的计算则是反方向的，即先对B的最后一行和最后一列进行初始化，然后反向推进到（0，0）。
这样，矩阵A与B的和（A+B=C）包含了经过任意特定断点（i、j）的最优比对得分。
我们称C矩阵为总得分矩阵，而A、B分别是前缀和后缀的得分矩阵。
根据矩阵C，我们可以迅速求出两个序列的最优比对，而不需要像上一节那样用反向跟踪的方法求出最优比对所对应的路径。
图3.19列出了根据某个打分模型计算得到的矩阵A和C。
可以看出，根据C的最大值，可非常容易地找出最优比对所对应的路径。
最优路径通过一系列断点，经过这一系列断点的所有最优比对得分值相同，实际上这个得分值就是两条序列比对的得分。
换句话说，这些断点都处于最优路径上。
虽然最优多重比对的投影不一定是两两最优比对，但是我们可以为投影的得分设立一个下限。
设α是关于s1,s2,...,sk的最优比对，如果SP-score(α)≥L，则可以证明(3-39)其中score(αi,j)是α在si,和sj所在平面投影的得分，(3-40)这里，L实际上是最优多重比对的一个下限，而Li,j是序列si和序列sj比对得分的一个下限。
当然，上述条件只是一个必要条件，但不是充分条件。
满足条件只是说明i可能处于最优路径，但不一定处于最优路径。
条件的作用是限制搜索空间，提高算法的实施效率。
这里，Cx，y是序列sx和sy的总得分矩阵，Cx，y[ix,iy]表示在点[ix,iy]处的值。
这个具体的条件是根据前面的叙述推演出来的。
假设最优比对α所对应的路径通过节点（i1,i2,…,ix,..,iy,…,ik），则Cx，y[ix,iy]≥score(αi,j)≥Lx，y。
为了得到一个合理的下限L，我们可以任选一个包含所有序列的多重比对，计算该多重比对的得分，以此作为L。
若选取的L接近于最优值，算法速度将大大提高。
注意：L的值不能超过最优值，否则算法出错。
在实现上述优化计算方法时，必须非常仔细。
逐步扩展节点，直到终点(n1,…,nk)。
以数组a[]保存各节点的计算结果。
称一个节点i影响另一个节点j，如果在计算a[j]时用到i。
这又称为j依赖于i，每一个节点依赖于其它2k-1个节点。
首先将0放入其中。
当一个节点i进入缓冲区时，其对应的值a[i]被初始化，然后a[i]的值在随后的阶段中被更新。
当该节点离开缓冲区时，其值即为该点真正的值，并用于其他节点（依赖于此节点）的计算。
幸运的是，在实际工作中，许多生物分子序列的比较满足这两个条件。
所以，必须考虑其他的方法，首选的就是所谓启发式方法。
目前所用的算法大部分将序列多重比对转化为序列两两比对，逐渐将两两比对组合起来，最终形成完整的多序列比对。
这种方式又称为渐进法。
可以按照星形结构或者树形结构组合两两序列比对。
星形比对是一种启发式方法，首先由Gusfield提出。
星形比对的基本思想是：在给定的若干条序列中，选择一个核心序列，通过该序列与其他序列的两两比对，形成所有序列的多重比对α，从而使得α在核心序列和任何一个其他序列方向的投影是最优的两两比对。
设s1,s2,...,sk是k条待比对的序列（暂时不考虑如何选择核心序列）。
假设已知核心序列是sc，c介于1到k之间，则可以利用标准的动态规划算法求出所有si和sc的最优两两比对。
将这些序列的两两比对聚集起来，并采用“只要是空位，则永远是空位”的原则。
聚集过程从某一个两两比对开始，如s1和sc，然后逐步加上其他的两两比对。
在这个过程中，逐步增加sc中的空位字符，以适应其他的比对，但决不删除sc中已存在的空位字符。
假设在上述过程中的某一时刻，有一个由sc指导的、已经建立好的部分序列的多重比对，接下来就是加入一个新的、与sc两两比对的序列，如果需要，则插入新的空位字符。
这个过程一直进行到所有的两两比对都加入以后结束，所需的计算量为O(k2n)。
那么，如何选择核心序列呢？一种方法就是尝试将每一个序列分别作为核心序列，按上述过程进行，取结果最好的一个。
另一种方法是计算所有的两两比对，取下式值最大的一个：(3-42)例如，有5个序列，s1=ATTGCCATTs2=ATGGCCATTs3=ATCCAATTTTs4=ATCTTCTTs5=ACTGACC根据它们之间的相似性度量建立如图3.20所示的表格，得分函数见公式（3-2）。
从表格中可以看出，选sc=s1时，式（3-42）取最大值。
接下来，计算s1与其他序列的最优两两比对。
假设得到下述各两两比对：ATTGCCATTATTGCCATT--ATTGCCATTATTGCCATTATGGCCATTATC-CAATTTTATCTTC-TTACTGACC--s1s2s3s4s5s17-20-3s27-20-4s3-2-20-7s4000-3s5-3-4-7-3图3.20两两比对的得分在这个例子中，核心序列s1中的空位字符是由于s3序列在末端引入的两个空位字符。
最后多重序列比对的结果如下：ATTGCCATT--ATGGCCATT--ATC-CAATTTTATCTTC-TT--ACTGACC----星形比对是一种多重序列比对近似的方法，然而是一种比较好的近似方法。
如果用代价来评判来多重序列的比对结果，则可以证明，用该方法所得到的多重序列比对的代价不会大于最优多重序列比对代价的两倍。
3.3.5树形比对假设有k条待比对的序列，同时有一棵具有k个叶节点的树，其中每个叶节点对应一条序列。
如果将序列也赋予树的内部节点，则可以计算树中每个分支的权值。
权值代表对应分支连接的两条序列之间的相似度。
所有权值的和就是这棵树（对应于将序列赋予内部节点的特定方式）的得分。
树形比对的问题就是寻找一种树的内部节点序列赋予方式，使得树的得分最大。
星形比对可以认为是树形比对的一种特例。
例如，图3.21中有4条简单序列，将CT、CG、CT分别赋予节点x、y、z，则树的得分为8。
这里，假设a=b时，p(a,b)=1，否则p(a,b)=0，p(a,-)=p(-，b)=-1。