1、第三章第三章 序列比对序列比对r寻找进化过程中的同源序列寻找进化过程中的同源序列;r基于同源物鉴定的功能预测基于同源物鉴定的功能预测;r基本假设:基本假设:序列的保守性序列的保守性 功能的保守性功能的保守性r注意:注意:1.蛋白质一般在三级结构的层面上执行功能;蛋白质一般在三级结构的层面上执行功能;2.蛋白质序列的保守性决定于其编码蛋白质序列的保守性决定于其编码DNA的保的保守性;守性;通常通常r第一节:数学基础:概率及概率模型第一节:数学基础:概率及概率模型r第二节:双序列比对算法的介绍第二节:双序列比对算法的介绍Dot matrix动态规划算法动态规划算法w(Needleman-Wunsc
2、h,Smith-Waterman算法算法)FASTA和和BLAST算法算法r第三节:打分矩阵及其含义第三节:打分矩阵及其含义r第四节:多序列比对第四节:多序列比对r从从N个物品中取出个物品中取出k个物品的个物品的排列排列数:数:r从从N个物品中取出个物品中取出k个物品的个物品的组合组合数:数:)!(!)1(.)1(kNNkNNNPkN)!(!kNkNkPkNCkNkNr概率模型概率模型:一个能够通过不同的概率产生不同结果一个能够通过不同的概率产生不同结果的模型。概率模型可以模拟或者仿真某一类型的的模型。概率模型可以模拟或者仿真某一类型的所有事件,并且对每个事件赋予一个概率。所有事件,并且对每个
3、事件赋予一个概率。161iipr色子模型:一个色子存在色子模型:一个色子存在6个概率值:个概率值:p1,p2,p6,其中,掷出其中,掷出i的概率为的概率为pi(i=1,2,6)。因此:。因此:pi0,且,且考虑三次连续的掷色子,结果为考虑三次连续的掷色子,结果为1,6,3,则总,则总概率为:概率为:p1p6p3r考虑连续变量考虑连续变量x,例如:物体的重量。则当,例如:物体的重量。则当重量确切为重量确切为1公斤时的概率,为公斤时的概率,为0。r变量的区间:变量的区间:P(x0 xx1)当区间无限小当区间无限小-0时,上式:时,上式:P(x-x/2 xx+x/2)=f(x)xf(x)称为概率密度
4、函数称为概率密度函数r因此:因此:10 xx10f(x)dx)XxP(X1)(dxxf且r1.事件只有两种可能出现的结果。例如掷硬事件只有两种可能出现的结果。例如掷硬币,正面记为币,正面记为“1”,反面记为,反面记为“0”。r2.则掷硬币则掷硬币N次,有次,有k次是次是1的概率为:的概率为:kNkppkNkP)1()(NpppkNkkkPNkkNkNk11)1()()1()1()()()(12122pNpppkNkkPkkNkNkNk期望值期望值 E(x)=方差方差 Var X=2r1.稀有事件发生的概率:在一个连续的时稀有事件发生的概率:在一个连续的时间或空间中,稀有离散变量出现的概率间或空
5、间中,稀有离散变量出现的概率r2.N-,E(x)=.2,1,0,!)()(xxexPxe=2.71828r对于大的对于大的N及小的及小的p值的二项分布,能够相值的二项分布,能够相当准确地用一个参数为当准确地用一个参数为=Np的泊松分布近的泊松分布近似。似。r当实验次数很多而概率很小时:当实验次数很多而概率很小时:二项分布二项分布泊松分布泊松分布r假设:需要测序的假设:需要测序的BAC长度长度200kbp;总共测序的序列数量:总共测序的序列数量:N;每次测序:每次测序:500bp;每次测序的覆盖率每次测序的覆盖率 p:500/200kbp=0.0025因此:每个点平均覆盖到的次数因此:每个点平均
6、覆盖到的次数:=N*p r k:测序能够覆盖到点测序能够覆盖到点X的次数。的次数。!)1()!(!)(keppkNkNkPkkNk点点X被覆盖被覆盖k次的概率:次的概率:(二项分布二项分布泊松分布泊松分布)当点当点X一次都不被覆盖时,一次都不被覆盖时,k=0;此时的概率为:此时的概率为:eP)0(rProf.Gene发现一种序列上的调控信号,在人的发现一种序列上的调控信号,在人的基因组上平均每基因组上平均每500kbp一个。那么,随机给一条一个。那么,随机给一条1mbp的序列,在上面发现的序列,在上面发现5个这样的信号,完全个这样的信号,完全是随机产生的概率是多少?是随机产生的概率是多少?r本
7、例中,本例中,E(x)=2(1mbp/500kbp)05.0036.0!52)5(52ePr统计显著性:统计显著性:p-value 0.05r与二项式分布的区别:不放回抽样。与二项式分布的区别:不放回抽样。r例:有例:有N个球,其中红球个球,其中红球M个,白球个,白球N-M个个,每次拿出一个球,每次拿出一个球再再放回,总共放回,总共n次,其中次,其中有有m个球是红球的概率为个球是红球的概率为(二项式分布二项式分布):mnmppmnmP)1()(p=M/Nr上例改为:有上例改为:有N个球,其中红球个球,其中红球M个,白球个,白球N-M个,每次拿出一个球个,每次拿出一个球不不放回,总共放回,总共n
8、次次,其中有,其中有m个球是红球的概率为:个球是红球的概率为:并且,并且,0mMNnNmnMNmMmP)(r上例再改为:有上例再改为:有N个球,其中红球个球,其中红球M个,白个,白球球N-M个,每次拿出一个球个,每次拿出一个球不不放回,总共放回,总共n次,其中次,其中至少有至少有m个球是红球的概率为:个球是红球的概率为:并且,并且,0mMNnmmnNmnMNmMmmPvaluep)(r上例再改为:有上例再改为:有N个球,其中红球个球,其中红球M个,白个,白球球N-M个,每次拿出一个球个,每次拿出一个球不不放回,总共放回,总共n次,其中次,其中最多有最多有m个球是红球的概率为:个球是红球的概率为
9、:并且,并且,0mMNmmnNmnMNmMmmPvaluep0)(r所有出现概率所有出现概率=观察表概率的概率之和观察表概率的概率之和)()(:)(mPiPiiPvaluepr超几何分布的精确概率计算。超几何分布的精确概率计算。r前提是固定边际分布,前提是固定边际分布,即即a+b、c+d、a+c与与b+d的值不变。的值不变。rRA Fisher,1935年文章示例:年文章示例:猜测先放的饮料牛奶 茶合计实际先放的饮料牛奶a=3b=1a+b=4茶c=1d=3c+d=4合计a+c=4b+d=4n=8r 计算公式:cancdcabaP=rP value:一种在原假设为真的前提下出现:一种在原假设为真
10、的前提下出现观察样本以及更极端情况的概率。观察样本以及更极端情况的概率。r显著性水平显著性水平A:认为预先设定的显著性水平:认为预先设定的显著性水平阈值,阈值,P A 为显著。为显著。r一般以一般以P 0.05 为显著,为显著,P O(log(n)O(n)O(n2)O(an)O(n!)r1.一般的,一般的,O(nk),当当k3 时,为多项式时时,为多项式时间,较为容易处理。间,较为容易处理。r2.当当O(an),则难以处理。,则难以处理。r3.NP完全问题(完全问题(NPC):无法找到能够在:无法找到能够在多项式时间复杂度内解决方法的问题;多项式时间复杂度内解决方法的问题;r4.近似算法近似算
11、法/优化算法,求近似解。优化算法,求近似解。r1900年,德国数学家年,德国数学家David Hilbert提出的提出的23个历史性数学难题。个历史性数学难题。r千禧年大奖难题美国克雷数学研究所(千禧年大奖难题美国克雷数学研究所(Clay Mathematics Institute,CMI)于)于2000年年5月月公布公布七个世界七个世界数学难题。数学难题。千禧年大奖难题P/NP问题:问题:P=NP?霍奇猜想庞加莱猜想(已证明)黎曼猜想杨-米尔斯存在性与质量间隙纳维-斯托克斯存在性与光滑性贝赫和斯维讷通-戴尔猜想rP问题:问题:Polynomial Problems可以在多项式可以在多项式(p
12、olynomial)时间内解决的问题时间内解决的问题;rNP:“Non-deterministic Polynomial”,并非并非“Non-Polynomial”可以在多项式的时间里验证一个解的问题可以在多项式的时间里验证一个解的问题;rNPC:NP-completer1.打分模型、替代矩阵以及空位罚分。打分模型、替代矩阵以及空位罚分。r2.比对算法:递归及动态规划算法;比对算法:递归及动态规划算法;r3.全局优化比对:全局优化比对:Needleman-Wunschr4.局部优化比对:局部优化比对:Smith-Watermanr5.工具资源:工具资源:http:/www.ludwig.edu
13、.au/course/lectures2005/Likic.pdfhttp:/ 2.修正的罚分:修正的罚分:d,第一次罚分的分数;第一次罚分的分数;g,空位数;,空位数;e,修正后的参数修正后的参数r 两条序列的比较,无空位:时间复杂度为两条序列的比较,无空位:时间复杂度为O(n2);r 两条序列比对,允许空位,时间复杂度为:两条序列比对,允许空位,时间复杂度为:r因此,有空位的双序列比对,时间复杂度为:因此,有空位的双序列比对,时间复杂度为:O(2O(22n2n),指数增加,指数增加,NPCNPC问题!问题!nnnnnn222)!()!2(2r 数学上保证提供最优解。数学上保证提供最优解。r
14、 动态规划算法:比较所有可能的字符对,动态规划算法:比较所有可能的字符对,考虑匹配、错配以及空位罚分,并且将比对考虑匹配、错配以及空位罚分,并且将比对次数控制在多项式时间内。次数控制在多项式时间内。r 替代矩阵:替代矩阵:BLOSUM62,空位罚分:空位罚分:11延伸的空位罚分:延伸的空位罚分:1(BLAST工具工具)r序列序列1:V D S C Yr序列序列2:V E S L C Yr替代矩阵中的分数:替代矩阵中的分数:4 2 4-11 9 7r 两序列比对的总分:两序列比对的总分:r Score=(AA pair scores)gap penalty=15GapVDSCYGap01gap
15、2gapV1gapE2gapSLCYgdgr)(本例:线性罚分GapVDSCYGap0-11-22-33-44-55V-11SijE-22S-33L-44C-55Y-66要求解要求解Sij的分数,我们必须先的分数,我们必须先知道知道Si-1,j-1,Si-1,j,以及以及Si,j-1的分的分数,这种方法叫做递归算法;数,这种方法叫做递归算法;采用这种方法,可以把大的问采用这种方法,可以把大的问题分割成小的问题逐一解决,题分割成小的问题逐一解决,即动态规划算法;需要存储如即动态规划算法;需要存储如何得到何得到Sij分数的过程。分数的过程。GapVDSCYGap0-11-22-33-44-55V-
16、11SijE-22S-33L-44C-55Y-66ijNeedleman-Wunsch算法;算法;Sij=max of Si-1,j-1+(xi,yj)Si-1,j-d(从左到右从左到右)Si,j-1-d(从上到下从上到下)GapVDSCYGap0-11-22-33-44-55V-11SijE-22S-33L-44C-55Y-66ij4-11-11Needleman-Wunsch算法;算法;Sij=max of Si-1,j-1+(xi,yj)Si-1,j-d(从左到右从左到右)Si,j-1-d(从上到下从上到下)GapVDSCYGap0-11-22-33-44-55V-114E-22S-33
17、L-44C-55Y-664-11-11GapVDSCYGap0-11-22-33-44-55V-114SijE-22S-33L-44C-55Y-66-3-11-11Needleman-Wunsch算法;算法;Sij=max of Si-1,j-1+(xi,yj)Si-1,j-d(从左到右从左到右)Si,j-1-d(从上到下从上到下)GapVDSCYGap0-11-22-33-44-55V-114-7E-22S-33L-44C-55Y-66-3-11-11GapVDSCYGap0-11-22-33-44-55V-114-7-18-29-40E-22-76-5-16-27S-33-18-510-1
18、-12L-44-29-16-19-3C-55-40-27-1287Y-66-51-38-23-31542GapVDSCYGap0-11-22-33-44-55V-114-7-18-29-40E-22-76-5-16-27S-33-18-510-1-12L-44-29-16-19-3C-55-40-27-1287Y-66-51-38-23-31542V D S C YV E S L C Yr下例:局部优化打分下例:局部优化打分r两条序列如下:两条序列如下:L D S C HG E S L C Kr 目标:使用局部优化算法寻找比对的结果目标:使用局部优化算法寻找比对的结果r1.Smith-Wate
19、rman算法;算法;r2.时间复杂度时间复杂度O(n2);r3.Sij=max of 0 Si-1,j-1+(xi,yj)Si-1,j-d(从左到右从左到右)Si,j-1-d(从上到下从上到下)r本例中:本例中:gap:12,线性罚分模型。,线性罚分模型。GapLDSCHGap000000G0SijE0S0L0C0K0Smith-Waterman算法;算法;Sij=max of Si-1,j-1+(xi,yj)Si-1,j-d(从左到右从左到右)Si,j-1-d(从上到下从上到下)0GapLDSCHGap000000G0SijE0S0L0C0K0-12-12-3GapLDSCHGap000000G000E0S0L0C0K0-12-12-4GapLDSCHGap000000G000000E002000S002600L040050C001092K000008-12-2GapLDSCHGap000000G000000E002000S002600L040050C001092K000008L D S C HG E S L C Kr1.局部优化比对:局部优化比对:9分分r2.全局优化比对:全局优化比对:2+4-12+9=3r3.为何不同?为何不同?L D S C HG E S L C Khttp:/zhanglab.ccmb.med.umich.edu/NW-align/