1、第三章第三章 序列比对序列比对主讲教师：丁彦蕊主讲教师：丁彦蕊单位：信息工程学院单位：信息工程学院序列比对：序列比对：寻找两条或者两条以上的序寻找两条或者两条以上的序列中各个字符的一一对应关系。列中各个字符的一一对应关系。序列比对的根本任务：序列比对的根本任务：发现序列之间的相似性发现序列之间的相似性寻找共同区域寻找共同区域辨别序列之间的差异辨别序列之间的差异目的：目的：相似序列相似序列 相似的相似的结构，相似的功能结构，相似的功能 判别序列之间的同源性判别序列之间的同源性推测序列之间的进化关系推测序列之间的进化关系 第一节第一节 序列的相似性序列的相似性同源（同源（homology）- 具有
2、共同的祖先具有共同的祖先(P81)直向同源（直向同源（Orthologous ）:共同祖先的直接后代共同祖先的直接后代（没有发生基因复制事件）之间的同源基因称为直（没有发生基因复制事件）之间的同源基因称为直向同源。向同源。 共生同源（共生同源（paralogous ）:两个物种两个物种A和和B的同源基因，的同源基因，分别是共同祖先基因组中由复制事件而产生的不同分别是共同祖先基因组中由复制事件而产生的不同拷贝的后代，这被称为共生同源基因。拷贝的后代，这被称为共生同源基因。相似（相似（similarity）同源序列一般是相似的同源序列一般是相似的 相似序列不一定是同源的相似序列不一定是同源的 进化
3、趋同（同功能）进化趋同（同功能）序列比较的基本操作是比对（序列比较的基本操作是比对（Alignment） 两个序列的比对是指这两个序列中各个字符的一种两个序列的比对是指这两个序列中各个字符的一种一一对应关系，或字符的对比排列一一对应关系，或字符的对比排列 。设有两个序列：GACGGATTAG，GATCGGAATAGAlignment2: GA CGGATTAGGATCGGAATAGAlignment1:GACGGATTAG GATCGGAATAG1、字母表和序列、字母表和序列字母表字母表4字符DNA字母表：A, C, G, T扩展的遗传学字母表或IUPAC编码单字母氨基酸编码符 号含 义说 明
4、GGGuanine AAAdenine TTThymine CCCytosineRG or APurine YT or CPyrimidine MA or CAmino KG or TKeto SG or CStrong interaction (3 H bonds) WA or TWeak interaction (2 H bonds) HA or C or TNot-GBG or T or Cnot-AVG or C or Anot-T(not-U) DG or A or Tnot-C NG or A or T or CAny 扩展的遗传学字母表或IUPAC编码氨基酸名称氨基酸名称简简 写
5、写氨基酸名称氨基酸名称简简 写写甘氨酸G丝氨酸S丙氨酸A苏氨酸T缬氨酸V天冬酰胺N异亮氨酸I谷酰胺Q亮氨酸L酪氨酸Y苯丙氨酸F组氨酸H脯氨酸P天冬氨酸D甲硫氨酸M谷氨酸E色氨酸W赖氨酸K半胱氨酸C精氨酸R特定的符号 代表字母表 A* 代表由字母表A中字符所形成的一系列有限长度序列或字符串的集合 a、b、c代表单独的字符 s、t、u、v代表A*中的序列 |s|代表序列s的长度为了说明序列s的子序列和s中单个字符，在s中各字符之间用数字标明分割边界例如，设s=ACCACGTA，则s可表示为 0A1C2C3A4C5G6T7A8 i:s:j 指明第i位或第j位之间的子序列, 当然，0 i j |s|。
6、 子序列0:s: i 称为前缀，即prefix(s,i)子序列 i:s:|s|称为后缀，即suffix(s, |s|-i+1) i:s: i 为空序列j-1:s:j 表示s 中的第j 个字符，简记为sj子序列与子串（子序列与子串（p82）子序列：子序列：选取s中的某些字符（或删除s中的某些字符）而形成s的子序列例如： TTT 是 ATATAT的子序列。 s的子串的子串：是由s中相继的字符所组成。 例如：TAC是AGTACA的子串，但不是TTGAC的子串（是子序列）。 子串是子序列子串是子序列子序列不一定是子串子序列不一定是子串字符串操作字符串连接操作：两个序列s和t的连接： s + + t例如
7、：ACC+CTA = ACCCTA 字符串k操作 删除字符串两端的字符 其定义如下：prefix(s,l) = sk|s|-lsuffix(s,l) = k|s|-ls i:s:j = ki-1sk|s|-j 序列比较可以分为四种基本情况序列比较可以分为四种基本情况(P83)（1）两条长度相近的序列相似）两条长度相近的序列相似 找出序列的差别找出序列的差别（2）判断一条序列的前缀与另一条序列的后缀相似）判断一条序列的前缀与另一条序列的后缀相似（3）判断一条序列是否是另一条序列的子序列）判断一条序列是否是另一条序列的子序列（4）判断两条序列中是否有非常相似的子序列）判断两条序列中是否有非常相似的
8、子序列2、编辑距离（、编辑距离（Edit Distance)GCATGACGAATCAG TATGACAAACAGC GCATGACGAATCAG TATGAC-AAACAGC 说明两条序列的相似程度说明两条序列的相似程度 定量计算定量计算 两条序列的相似程度的定量计算相似度相似度，它是两个序列的函数，其值越大，表示两个序列越相似 两个序列之间的距离距离。距离越大，则两个序列的相似度就越小 字符编辑操作（Edit Operation）字符编辑操作可将一个序列转化为一个新序列 Match（a，a）Delete（a，-） Replace（a，b）Insert（-，b）直接距离计算的不足直接距离计算
9、的不足扩展的编辑操作ACCGACAATATGCATA ATAGGTATAACAGTCAACCGACAATATGCATA ACTGACAATATGGATA 第二条序列头尾颠倒CTAGTCGAGGCAATCTGAACAGCTTCGTTAGT ？反向互补序列反向互补序列RNA发夹式二级结构发夹式二级结构3、通过点矩阵进行序列比较、通过点矩阵进行序列比较“矩阵作图法矩阵作图法” 或或 “对角线作图对角线作图” 序列序列1 序列序列2 实实 例例 序列序列1 序列序列1 自我比较自我比较滑动窗口技术滑动窗口技术两条序列中有很多匹配的字符对，因而在点矩阵中两条序列中有很多匹配的字符对，因而在点矩阵中会形成
10、很多点标记。会形成很多点标记。滑动窗口技术滑动窗口技术 使用滑动窗口代替一次一个位点的比较是解决这使用滑动窗口代替一次一个位点的比较是解决这个问题的有效方法。个问题的有效方法。 假设窗口大小为假设窗口大小为10，相似度阈值为，相似度阈值为8，则每次比较，则每次比较取取10个连续的字符，如相同的字符超过个连续的字符，如相同的字符超过8个，则标个，则标记记 基于滑动窗口的点矩阵方法可以明显地降低点阵基于滑动窗口的点矩阵方法可以明显地降低点阵图的噪声，并且明确无误的指示出了两条序列间具图的噪声，并且明确无误的指示出了两条序列间具有显著相似性的区域。有显著相似性的区域。 （a）对人类（）对人类（Hom
11、o sapiens）与黑猩猩（）与黑猩猩（Pongo pygmaeus）的）的球蛋白基因球蛋白基因序列进行比较的完整点阵图。（序列进行比较的完整点阵图。（b）利用滑动窗口对以上的两种球蛋白基因序列进）利用滑动窗口对以上的两种球蛋白基因序列进行比较的点阵图，其中窗口大小为行比较的点阵图，其中窗口大小为10个核苷酸，相似度阈值为个核苷酸，相似度阈值为8。 (a) (b) 具有连续相似区域的两条具有连续相似区域的两条DNA序列的简单点阵图序列的简单点阵图4、 序列的两两比对序列的两两比对序列的两两比对序列的两两比对(Pairwise Sequence Alignment）:通过字符匹配和替换，或者通
12、过字符匹配和替换，或者插入和删除字符，使两条序列达到一样插入和删除字符，使两条序列达到一样的长度，并使两条序列中相同的字符尽的长度，并使两条序列中相同的字符尽可能一一对应。可能一一对应。 s:AGCACACAAGCACACA t:ACACACTAACACACTA Match(A, A)Match(A, A)Delete(G, - )Replace(G, C)Match(C, C)Insert( -, A)Match(A, A)Match(C, C)Match(C, C)Match(A, A)Match(A, A)Match(C, C)Match(C, C)Replace(A, T)Insert
13、( -, T)Delete(C, -)Match(A, A)Match(A, A)图3.6 序列AGCACACA和ACACACTA的两种比对结果Alignment -1 Alignment -2不同编辑操作的代价不同不同编辑操作的代价不同编辑操作定义函数w，它表示“代价（cost）”或“权重（weight）”。对字母表中的任意字符a、b，定义 w (a, a) = 0 w (a, b) = 1 a b w (a, -) = w ( -, b) = 1 也可以使用得分（score）函数来评价编辑操作 p (a, a) = 1 p (a, b) = 0 a b p (a, -) = w ( -,
14、b) = -1 概念概念两条序列两条序列s 和和 t 的比对的得分（或代价）等于将的比对的得分（或代价）等于将s 转化为转化为t 所用的所有编辑操作的得分（或代价）总所用的所有编辑操作的得分（或代价）总和；和；s 和和t 的最优比对是所有可能的比对中得分最高的最优比对是所有可能的比对中得分最高（或代价最小）的一个比对；（或代价最小）的一个比对；s 和和t 的真实距离应该是在得分函数的真实距离应该是在得分函数p值（或代价值（或代价函数函数w值）最优时的距离。值）最优时的距离。 例如：例如：s:AGCACAC At:A CACACTA cost=2 s:AGCACAC A t:A CACACTA
15、score (s，t）= 5序列比对的目的是寻找一个得分最大（或代价序列比对的目的是寻找一个得分最大（或代价最小）的比对。最小）的比对。5、打分矩阵（、打分矩阵（Weight Matrices）(P87) （1）核酸打分矩阵设DNA序列所用的字母表为 = A，C，G，T a. 等价矩阵（相同核苷酸得分为1，不同核苷酸替换得分为0） b. BLAST矩阵（相同核苷酸得分为+5，不同核苷酸得分为-4）c. 转移矩阵（transition，transversion） （嘌呤：腺嘌呤A，鸟嘌呤G；嘧啶：胞嘧啶C，胸腺嘧啶T） ATCGA1000T0100C0010G0001ATCGA5-4-4-4T-
16、45-4-4C-4-45-4G-4-4-45ATCGA1-5-5-1T-51-1-5C-5-11-5G-1-5-51表3.1 等价矩阵表表3.3 转移矩阵表3.2 BLAST矩阵（2）蛋白质打分矩阵（i）等价矩阵等价矩阵jijiRij01其中Rij代表打分矩阵元素i、j分别代表字母表第i和第j个字符。（ii） 氨基酸突变代价矩阵氨基酸突变代价矩阵GCM GCM矩阵通过计算一个氨基酸残基转变到另外一个氨基酸残基所需矩阵通过计算一个氨基酸残基转变到另外一个氨基酸残基所需的密码子变化数目而得到，矩阵元素的值对应于代价。如果变化一的密码子变化数目而得到，矩阵元素的值对应于代价。如果变化一个碱基就可以使
17、一个氨基酸的密码子改变为另一个氨基酸的密码子，个碱基就可以使一个氨基酸的密码子改变为另一个氨基酸的密码子，则这两个氨基酸的替换代价为则这两个氨基酸的替换代价为1，如果需要两个碱基的改变，则替换，如果需要两个碱基的改变，则替换代价为代价为2，以此类推。，以此类推。 GCM矩阵常用于进化距离的计算，其优点是计算结果可以直接用于矩阵常用于进化距离的计算，其优点是计算结果可以直接用于绘制进化树，但在蛋白质序列比对尤其是相似程度较低的序列比对绘制进化树，但在蛋白质序列比对尤其是相似程度较低的序列比对中很少使用。中很少使用。（iii）疏水矩阵）疏水矩阵根据氨基酸残基替换前后疏水性的变化得到的矩阵。如果氨根
18、据氨基酸残基替换前后疏水性的变化得到的矩阵。如果氨基酸基酸A被氨基酸被氨基酸B替换后，疏水性变化不大则替换得分高，替换后，疏水性变化不大则替换得分高，否则替换得分低。否则替换得分低。（iv）PAM矩阵（矩阵（Point Accepted Mutation）统计自然界中各种氨基酸残基的相互替换率。如果两种特定统计自然界中各种氨基酸残基的相互替换率。如果两种特定的氨基酸之间替换发生得比较频繁，则这一对氨基酸在得分的氨基酸之间替换发生得比较频繁，则这一对氨基酸在得分矩阵中的互换得分就高。矩阵中的互换得分就高。PAM矩阵基于进化原理，建立在进化的点接受突变模型基矩阵基于进化原理，建立在进化的点接受突变
19、模型基础上，通过统计相似序列中的各种氨基酸替换发生率而得到础上，通过统计相似序列中的各种氨基酸替换发生率而得到的矩阵。的矩阵。PAM矩阵（矩阵（Point Accepted Mutation） 基于进化的点突变模型基于进化的点突变模型 一个一个PAM就是一个进化的变异单位就是一个进化的变异单位, 即即1%的氨基酸改变的氨基酸改变 这类矩阵里列出同源蛋白质在进化过程中氨基酸变化的可能性。这类矩阵里列出同源蛋白质在进化过程中氨基酸变化的可能性。 这类矩阵式基于进化原理的这类矩阵式基于进化原理的 证据：证据： 编码相同蛋白质的基因随着进化发生分歧，相似度降低。编码相同蛋白质的基因随着进化发生分歧，相
20、似度降低。 科学科学 用得多用得多 矩阵集合矩阵集合- PAM-N如，如，PAM120矩阵用于比较相距矩阵用于比较相距120个个PAM单位的序列。单位的序列。一个一个PAM-N矩阵矩阵元素（元素（i，j）的值：的值： 反应两个相距反应两个相距N个个PAM单位的序列中第单位的序列中第i种氨基酸种氨基酸替替换换第第j种氨基酸的频率。种氨基酸的频率。针对不同的进化距离采用针对不同的进化距离采用PAM 矩阵矩阵序列相似度序列相似度 = 40% 50% 60% | | |打分矩阵打分矩阵 = PAM120 PAM80 PAM 60PAM250 14% - 27% （v） BLOSUM矩阵矩阵 （Bloc
21、ks Amino Acid Substitution Matrices）通过统计相似蛋白质序列的替换率得到的。PAM矩阵是从蛋白质序列的全局比对结果推导出来的，而BLOSUM矩阵是从蛋白质序列块比对而推导出来的。BLOSUM 62第二节第二节 两两比对算法两两比对算法1、序列两两比对基本算法、序列两两比对基本算法直接方法直接方法 生成两个序列所有可能的比对，分别计算代生成两个序列所有可能的比对，分别计算代价函数，然后挑选一个代价最小的比对作为最终结果。价函数，然后挑选一个代价最小的比对作为最终结果。本质问题：优化本质问题：优化动态规划寻优策略动态规划寻优策略动态规划算法（动态规划算法（Dyna
22、mic Programming）（）（P93）最短路经问题最短路经问题起点起点终点终点C1 C2 W1 W2路径1：C1 + w1 ？路径2：C2 + w2 ？ 取最小值！取最小值！算法求解算法求解: 从起点到终点逐层计算从起点到终点逐层计算利用利用动态规划动态规划方法求解序列的两两比对方法求解序列的两两比对起点起点终点终点ATTCCGAAGA AGTCGAAGGTATTCCGAAG AGTCGAAGGAT+（1）ATTCCGAAGA AGTCGAAGG-T+（2）ATTCCGAAG AGTCGAAGGTA-+（3）求解过程求解过程起点起点终点终点ATTCCGAAGA AGTCGAAGGT 从
23、两个序列前端开始从两个序列前端开始 逐步推进逐步推进 直到两个序列的末端。直到两个序列的末端。序列S: i-1 i i+1序列t: j-1 j j+1序列S: i-1 i i+1序列t: j-1 j j+1Case1:匹配（si，tj )中间过程：比对中间过程：比对0:S:i 与与 0:T:j 序列S: i-1 i i+1序列t: j-1 j j+1序列S: i-1 i i+1序列t: j-1 j j+1Case2:删除（si， -） 序列S: i-1 i i+1序列t: j-1 j j+1序列S: i-1 i i+1序列t: j-1 j j+1Case3:插入（ -，tj ） 设序列设序列s
24、、t的长度分别为的长度分别为m和和n。考虑两个前缀考虑两个前缀 0:s:i 0:t:j 假如已知序列假如已知序列0:s:i 和和0:t:j 所有较短子列的最优比对，即已知：所有较短子列的最优比对，即已知：（1）0:s:(i-1) 和和 0:t:(j-1) 的最优比对的最优比对（2） 0:s:(i-1) 和和 0:t:j 的最优比对的最优比对（3） 0:s:i 和和 0:t:(j-1) 的最优比对的最优比对则则0:s:i和和 0:t:j 的最优比对一定是上述三种情况之一的扩展的最优比对一定是上述三种情况之一的扩展（1）替换（）替换（si，tj）或匹配（）或匹配（si，tj ) ，这取决于，这取决
25、于si 是否等于是否等于tj ；（2）删除（）删除（si， -）；）；（3）插入（）插入（ -，tj ）。）。 ):, :(00jitsS令：令：为序列为序列0:s:i和与序列和与序列 0:t:j 比对的得分比对的得分按下述方法求解按下述方法求解 其初值为：其初值为：for i=1 , 2 ,., mfor j=1 , 2 ,., n),():, :(),():,:(),():,:(max):, :()1(000)1(0)1(0)1(000jjiijijijijitptsSsptsStsptsStsS),():,:():,:(),():,:():, :(0):,:() 1(00000000)
26、1(00000000jjjiiitptsStsSsptsStsStsS距离矩阵距离矩阵按照上述方法，对于给定的得分函数按照上述方法，对于给定的得分函数p，两，两个序列所有前缀的得分定义了一个个序列所有前缀的得分定义了一个(m+1) (n+1)的距离矩阵的距离矩阵D = ( d i , j )其中其中d i , j = S (0:s:i , 0:t:j ) d i , j的计算公式如下：的计算公式如下：),(),(),(max1, 11, 1,jjiijijijijitpdspdtspddd i , j 最小值的三种选择决定了各矩阵元素之间的关系，最小值的三种选择决定了各矩阵元素之间的关系，用下
27、图表示：用下图表示：di,jdi,j-1di-1,jdi-1,j-1 距离矩阵元素d i , j 的计算S (0:s:i , 0:t:j )S (0:s:i-1 , 0:t:j )S (0:s:i-1 , 0:t:j-1 )S (0:s:i , 0:t:j-1 )动态规划算法计算过程动态规划算法计算过程： 计算过程从计算过程从d 0 , 0开始开始 可以是按行计算，每行从左到右，也可以是按列计算，每可以是按行计算，每行从左到右，也可以是按列计算，每列从上到下。列从上到下。 当然，任何计算过程，只要满足在计算当然，任何计算过程，只要满足在计算d i , j时时 d i-1 , j、d i-1 ,
28、 j-1、和、和d i, j-1都已经被计算这个条件即可。都已经被计算这个条件即可。 在计算在计算d i , j后，需要保存后，需要保存d i , j是从是从d i-1 , j、d i-1 , j-1、或、或d i, j-1中的哪一个推进的，或保存计算的路径，以便于后续处理。中的哪一个推进的，或保存计算的路径，以便于后续处理。上述计算过程到上述计算过程到d m , n结束。结束。 最优路径求解：最优路径求解：与计算过程相反与计算过程相反 从从d m , n开始，反向前推。开始，反向前推。 假设在反推时到达假设在反推时到达d i ,j，根据保存的计算路径判断，根据保存的计算路径判断d i , j
29、究竟是根据究竟是根据d i-1 , j、d i-1 , j-1、和、和d i, j-1中的那一个中的那一个计算而得到的。找到这个点以后，再从此点出发，计算而得到的。找到这个点以后，再从此点出发，一直到一直到d 0 , 0为止。为止。 走过的这条路径就是最优路径（即代价最小路走过的这条路径就是最优路径（即代价最小路径），其对应于两个序列的最优比对。径），其对应于两个序列的最优比对。 计算过程：计算过程：（1）初始化）初始化计算过程：计算过程：（2）反复计算）反复计算按列计算按列计算计算过程：计算过程：（2）反复计算）反复计算按行计算按行计算其他方式其他方式计算过程：计算过程：（3）求最佳路径）求
30、最佳路径t ts s ACACACTAAGCACACA例：例：s = AGCACACAt = ACACACTA 得分矩阵得分矩阵D （99 9）t ts s ACACACTA0-1-2-3-4-5-6-7-8A-1G-2C-3A-4C-5A-6C-7A-8初始化初始化计算计算d（2，2）t ts s ACACACTA0-1-2-3-4-5-6-7-8A-110-1-2-3-4-5-6G-201C-3A-4C-5A-6C-7A-8最终的得分矩阵最终的得分矩阵及序列比对及序列比对t ts s ACACACTA0-1-2-3-4-5-6-7-8A-110-1-2-3-4-5-6G-2010-1-2-
31、3-4-5C-3-11110-1-2-3A-4-2021210-1C-5-3-1132321A-6-4-2024333C-7-5-3-113543A-8-6-4-202455AGCACAC A| | |A CACACTA举例例 1: Si 和 Tj对齐 S: C A T T C A C T: C - T T C A Gi - 1ijj -1 S: C A T T C A - C T: C - T T C A G - 例 2: Si 中加入空位i - 1ij S: C A T T C A C - T: C - T T C A - G例 3: Tj 是入空位i jj -1计算过程C(n,m)C(0
32、,0)C(i,j) 1j, i (C,) j, 1i (C),T,S(w) 1j, 1i (Cmax) j, i (CjiC(i-1,j)C(i-1,j-1)C(i,j-1) C T C G C A G CACTTCAC+10 匹配, -2 不匹配, -5 空位 0 -5 -10 -15 -20 -25 -30 -35 -40-5-10-15-20-25-30-3510 50-5-10-15-20-25-30-35-40-51050-5-10-15-20-25-10583-2-70-5-10-150151050-5-2-7-20-5101383-2-7-4-25-1052015181383-3
33、0-150151813282318-35-20-5101328232633 C T C G C A G CACTTCAC回溯得到最佳的比对*C - T C G C A G CC A T - T C A - CC - T C G C A G CC A T T - C A - C第一种比对方式第二种比对方式2 2、子序列与完整序列的比对、子序列与完整序列的比对-AGCT-ATGCAGCTGCTT目标：目标：使使S(s, i:t:j ) 最大最大序列S:序列t: i j不计前缀0:t:i 的得分, 也不计删除后缀的j+1:t:|t|得分局部序列比对给定两条序列0:s:m和0:t:n，从t中寻找一个子
34、序列i:t:j使得S(s,i:t:j)最大.不计前缀不计前缀0:t:i 的得分的得分处理第一行处理第一行t ts s ACACACTA000000000AGCACACA0):,:(000itsS不计删除后缀的不计删除后缀的j+1:t:|t|得分得分 处理最后一行处理最后一行):,:(),():,:(),():,:(max):,:() 1(000) 1(0) 1(0) 1(000jmmjmjmjmjmtsSsptsStsptsStsSdm,jdm,j-1dm-1,jdm-1,j-1S (0:s:i , 0:t:j )S (0:s:i-1 , 0:t:j )S (0:s:i-1 , 0:t:j-1
35、 )S (0:s:i , 0:t:j-1 )不计代价不计代价距离矩阵初始化时，对第一行进行如下处理：距离矩阵初始化时，对第一行进行如下处理：d0,j = 0 for 0 j n 最后一行的计算应该是：最后一行的计算应该是：同样，同样，d m, n依然是最优局部比对的得分，而匹配的子列依然是最优局部比对的得分，而匹配的子列i:t:j 按如下方式寻找：按如下方式寻找： （1） j = min k d m ,k = d m ,n （2）由（）由（m,j）反推比对路径，最终通过斜线（非空位）到达（）反推比对路径，最终通过斜线（非空位）到达（0，i）。）。 （3-10）（3-11） 1, 11, 1,)
36、,(),(maxjmmjmjmjmjmdspdtspdd3、寻找最大的相似子序列、寻找最大的相似子序列目标：目标：使使dw (i:s :j, i:t:j ) 最大最大序列S:序列t: i ji j数据结构：数据结构：（m+1） （n+1）的矩阵）的矩阵 D但是，对数组元素含义解释与基本算法有所不同但是，对数组元素含义解释与基本算法有所不同每个元素的值代表序列每个元素的值代表序列0:s:i 某个某个后缀后缀和序列和序列0:t:j 某个某个后缀后缀的最佳比对。的最佳比对。 这种局部比对不计前缀的得分，所以新的边界条件是：d0,j = 0 for 0 j n （3-12）di,0 = 0for 1
37、i m另外，由于0:s:i 和0:t:j 总有一个得分为“0”的空后缀比对，因此矩阵D中的所有元素大于或等于“0”，于是，新的递归计算公式为：(3-13)0),(),(),(max1, 11, 1,jjiijijijijitpdspdtspdd寻找最佳比对的子序列寻找最佳比对的子序列在矩阵中找最大值在矩阵中找最大值该值就是最优的局部比对得分该值就是最优的局部比对得分该值对应的点为序列局部比对的末点该值对应的点为序列局部比对的末点然后反向推演前面的最优路径，直到局部比对的起点。然后反向推演前面的最优路径，直到局部比对的起点。 TATA|TATA4、准全局比较、准全局比较所谓准全局比较就是在评价序
38、列比对时不计终端“空缺”（end space，或空位）的得分或代价 序列序列1 长度为长度为8序列序列2 长度为长度为18（a）6个匹配，个匹配，1个失配，个失配，1个空位个空位 （b）8个匹配个匹配情况1 ：不记s后面的空位与 t 后缀比对的得分在矩阵di,j中取最后一行的最大值.序列S:序列t: i ji j空位后缀情况2 ：不记s前面的空位与 t 前缀比对的得分将矩阵di,j中的第一行各元素值置为“0”序列S:序列t: i ji j空位前缀情况3：情况4：半全局比较算法与基本算法在计算半全局比较算法与基本算法在计算di,j时的区别归纳为时的区别归纳为下列四个方面：下列四个方面：（1）第一
39、行初始值为“0”，表示不计第一个序列的前端空位；（2）寻找最后一行的最大值，表示不计第一个序列的末端空位；（3）第一列初始值为“0”，表示不计第二个序列的前端空位；（4）寻找最后一列的最大值，表示不计第二个序列的末端空位。对于最后一行和最后一列的另一种处理办法是：对于最后一行和最后一列的另一种处理办法是：最后一行的横向移动不被空位罚分最后一行的横向移动不被空位罚分最后一列的纵向移动也不被罚分最后一列的纵向移动也不被罚分这样，就可以允许在两条序列终端自由存在空位。这样，就可以允许在两条序列终端自由存在空位。 当矩阵当矩阵D所有元素计算完以后，其右下角得值即为两条序列最所有元素计算完以后，其右下角
40、得值即为两条序列最终准全局比对的得分。终准全局比对的得分。 ACACTGATCG|ACACTG5、关于连续空位的问题、关于连续空位的问题K 阶空位 K个连续的空位字符 “-”ATG-A-T-C-A-GATG-ATCAGATGCAGTGCAATGATGTTTTTATCAG生物学意义 “插入” 或“删除” 突变突变次数连续空位可能对应于一次突变非连续空位对应于 多次突变对于连续空位的代价是一个线性的函数。设p（k）代表空位得分函数，其中k是连续空位的个数，则： p（k）= -bk 这里b（0）是单个“空位”得分的绝对值。处理方法：任何一个比对可以被唯一地分为若干个相继的块。有三类块：（1）两个字符
41、的比对（2）与序列s空位进行比对的t的最大连续字符序列（3） 与序列t空位进行比对的s的最大连续字符序列上述算法的时间复杂度为上述算法的时间复杂度为O(n3)。比起标准算法，其多花的时间主要用于处理连续的空位。那么，是比起标准算法，其多花的时间主要用于处理连续的空位。那么，是否可以改进连续空位的得分函数，而使得算法的时间复杂度降低为否可以改进连续空位的得分函数，而使得算法的时间复杂度降低为O(n2)呢？呢？如果认为k个连续空位比k个孤立空位出现的可能性更大，则p(k) kp(1) （3-22）或更一般地，p(k1 + k2 + + kn ) p(k1) + p(k2) + +p(kn) (3-
42、23)可以用下式重新计算连续“空位”的得分：p(0)=0（3-24） p(k) = h g(k-1), k1 （3-25） h0，g0，hg。 6、比较相似序列、比较相似序列相似序列快速比较算法 例如，有两个序列： s=GCGCATGGATTGAGCGA t=TGCGCCATGGATGAGCA 最优比对所对应的路径偏离主对角线，经过一段以后重新返回主对角线。经验法则（针对蛋白质序列）：经验法则（针对蛋白质序列）： 如果两个序列的长度都大于如果两个序列的长度都大于100，在适当地加入空，在适当地加入空位之后，它们配对的相同率达到位之后，它们配对的相同率达到25%以上，则两个以上，则两个序列相关；
43、序列相关； 如果配对的相同率小于如果配对的相同率小于15%，则不管两个序列的，则不管两个序列的长度如何，它们都不可能相关；长度如何，它们都不可能相关； 如果两个序列的相同率在如果两个序列的相同率在15% 25%之间，它们可之间，它们可能是相关的。能是相关的。第三节第三节 序列多重比对序列多重比对多序列比对的定义设有k个序列s1, s2, . ,sk，每个序列由同一个字母表中的字符组成，k大于2。通过插入操作，使得各序列达到一样的长度。多序列比对的目的发现多个序列的共性发现与结构和功能相关的保守序列片段多序列比对的应用研究蛋白质之间的关系研究一个家族中的相关蛋白质研究相关蛋白质的保守区域研究蛋白
44、质的结构和功能进行蛋白质结构预测推测各个序列的进化历史1、SP（Sum-of-Pairs）模型）模型评价多重序列比对的结果评价多重序列比对的结果按照每个对比的列进行打分，然后加和按照每个对比的列进行打分，然后加和(P106)处理每一列：处理每一列： k个变量的打分函数个变量的打分函数 用一个用一个k维数组来表示该显式函数（类似于打分矩阵）维数组来表示该显式函数（类似于打分矩阵）期望：期望：函数在形式上应该简单函数在形式上应该简单具有统一的形式具有统一的形式不随序列的个数而发生形式变化不随序列的个数而发生形式变化 11121),(),.,(kikijjikccpcccscoreSP其中，c1,c
45、2,ck是一列中的k个字符，p是关于一对字符相似性的打分函数。逐对加和逐对加和SP（sum-of-pairs）函数）函数(P106) 26GSGPALLscoreSP逐对计算逐对计算p(1,2)，p (1,3)，.，p(1,8)，p (2,3)，p(2,4)，.， p (2,8)，.，p (7,8) 的所有得分的所有得分（-7-6-5-4-3-2-1）+2 = -26 另一种计算方式：先处理每一个序列对另一种计算方式：先处理每一个序列对在处理序列对时，逐个计算字符对，最后加和在处理序列对时，逐个计算字符对，最后加和则则SP得分模型的计算公式如下：得分模型的计算公式如下：jiijscoreSP)
46、( 是一个多重比对是一个多重比对 ij是由是由 推演出来的序列推演出来的序列s i 和和s j的两两比对的两两比对 2、多重比对的动态规划算法、多重比对的动态规划算法多重序列比对的最终目标是得到一个得分最高多重序列比对的最终目标是得到一个得分最高（或代价最小）的序列对比排列，从而分析各（或代价最小）的序列对比排列，从而分析各序列之间的相似性和差异序列之间的相似性和差异。前趋节点的个数等于前趋节点的个数等于2k - 1 )(),(jjjkjjisccbisColumn假设以k维数组A存放超晶格，则计算过程如下： a 0, 0, ,0 = 0 a i = max a i - b + SP-scor
47、e(Column(s, i, b) (3-37)(3-38) if bj = 1if bj = 0多序列比对的过程实际是一个递推过程，在计算每个晶格节点得分的时候，将其各前趋节点的值分别加上从前趋节点到当前点的SP得分，然后取最大值作为当前节点的值。计算量问题对于k条序列的比对，动态规划算法需要处理k维空间里的每个节点，计算量与晶格中的节点数成正比，而节点数等于各序列长度的乘积，况且计算的每个节点依赖于前趋节点的个数，因此用动态规划进行多重序列比对的时间复杂度为O(2k i=1,.,k si )即即O(2kNk) 利用SP模型寻找最优多重序列比对是一个NP-完全类问题完全类问题。NP-完全问题
48、通常被认为是一些人们难以在有限的时间、空间内对问题求出最佳解得问题，几乎所有专家都认为不可能在多项式时间内准确求解的问题。解决办法1. 只求解规模比较小的问题2. 利用动态规划、分支约束等技术减小搜索空间，提高求解问题的效率。3. 针对具体问题的特点，根据实际情况，设计实用求解算法。4. 采用近似算法或者启发式方法，如局部搜索、模拟退火、遗传算法等。3、 优化计算方法（优化计算方法（p110)标准动态规划算法存在的问题：标准动态规划算法存在的问题： 搜索搜索空间大空间大对于两两序列比对，最优的路径常常处对于两两序列比对，最优的路径常常处于对角线附近，而对于多序列比对，最于对角线附近，而对于多序
49、列比对，最优的路径不可能在某个平面上，而是在优的路径不可能在某个平面上，而是在某一个区域范围内。某一个区域范围内。利用人工智能空间搜索策略的剪枝技术，根据问题本身的特殊性将搜索空间限定在一个较小的区域范围内。若问题是搜索一条得分最高（或代价最小）的路径，则在搜索时如果当前路径的得分低于某个下限（或累积代价已经超过某个上限），则对当前路径进行剪枝，即不再搜索当前路径的后续空间。4、星形比对、星形比对(P112)星形比对的基本思想是：在给定的若干序列中，选择星形比对的基本思想是：在给定的若干序列中，选择一个核心序列，通过该序列与其它序列的两两比对形一个核心序列，通过该序列与其它序列的两两比对形成所
50、有序列的多重比对成所有序列的多重比对 ，从而使得，从而使得 在核心序列和任在核心序列和任何一个其它序列方向的投影是最优的两两比对。何一个其它序列方向的投影是最优的两两比对。利用标准的动态规划方法求出所有利用标准的动态规划方法求出所有si和和sc的最优两两比的最优两两比对对时间为时间为O（2knk）将这些两两比对聚集起来将这些两两比对聚集起来并采用并采用“只要是空位只要是空位 ， 则永远是空位则永远是空位”的原则。的原则。scs1s2sk(sc, s1) (sc, s2) (sc, sk)两两比对两两比对 多重比对多重比对如何选择核心序列？尝试将每一个序列分别作为核心序列，进行星形多重序列比对，
