序列分析(一)一一序列比对生物信息学研究的三个层面初级层面:基于现有的生物信息数据库和资源,利用成熟的生物信息学工具(专业网站、软件)解决生物信息学问题——生物信息数据库(NCBI、EBI、DDBJ、UniProt等)——基因组序列分析、序列比对软件(BLAST、CLUSTAL等)——系统发育树构建软件的简单使用(PHYLIP、PALM等)——搜集、整理有特色的生物信息学数据库中级层面:利用数理统计方法和相关的工具,研究生物信息学问题——概率、数理统计基础——现有的数理统计和科学计算工具(EXCEL、SPSS等)高级层面:提出有重要意义的生物信息学问题;自主创新,发展新方法,开发新工具,引领生物信息学领域研究方向。——面向生物学领域,解决重要生物学问题——利用数学、物理、化学、计算科学等思想和方法——建立模型,发展算法——自行编程,开发软件序列分析内容与DNA和protein序列相关的研究都可称为序列分析。主要包括:1.序列比对2.基因组序列分析3.蛋白质序列分析4.综合序列分析为什么要进行序列比对???序列拼接数据库搜索方面进化方面功能方面与进化相关的几个概念同源性与相似性:极易被混淆的两个概念!!!同源性(homology):是指序列们是由共同祖先进化而来,讲两条序列的同源关系,只有两种情况:同源、不同源。相似性:指序列间的差别,是一个度量。同源与相似的关系:一般认为序列相似性达到一定程度,即可认为是同源,但不绝对。Ortholog(直系同源):两个基因通过物种形成的事件而产生,或源于不同物种的具有共同祖先的两个基因,或者两个物种中的同一基因,一般具有相同的功能。Paralog(旁系同源):指相同的基因组内因基因复制形成的多个具有不同功能的基因。Xenolog(异同源):由某一个水平基因转移事件而得到的同源序列。Convergence(趋同):序列的相似性是由随机因素产生。第一节、双序列比对点阵分析动态规划1.点阵分析用途:1.寻找两条序列间所有可能的比对;2.寻找蛋白质、DNA序列上正向或反向的重复序列;3.发现RNA上可能存在的互补区域。优点:1.可以找到两个序列间所有可能的残基匹配;2.简单、易懂3.直观、整体性强工具:例1:自身的比对AKGFKCADEA100000100K10010000G1000000F100000K10000C1000A100D10E1例2:重复序列AKGFDKGFEA100000000K10001000G1000100F100010D10000K11000G1100F110E1例3:反向重复/回文AUGCACGUCA100010000U10000010G1000100C101000A10000C11001G1100U110C1例4:RNAstem/loopAUGUAGCAUA100010010U10100001G1001000C000001U00001A0010C100A10U1例5:不同序列的比对PKDFCKALVP100000000K10001000F0100000T00000K11000A100I00V1PKDFCKALVPK-FTKAIV字符的关系:1.匹配2.删除或插入3.替换点阵法的序列比对Sequence1#1nSequence2#1m“-”Insertion“-”Insertion寻找两条序列的最佳比对,实际上就是寻找在矩阵标记图中找非重叠平行斜线最长的组合。→序列1→→序列2→实例对于较长的序列,有很多匹配的字符,点阵图变得非常复杂和模糊。点阵法的滑动窗口技术:使用滑动窗口代替一次一个位点的比较,是解决这个问题的有效方法。假设窗口大小为10,相似度阈值为8,则每次比较取10个连续的字符,如相同的字符超过8个,则标记,基于滑动窗口的点矩阵方法可以明显地降低点阵图的噪声,并且明确无误的指示出了两条序列间具有显著相似性的区域。(a)对人类(Homosapiens)与黑猩猩(Pongopygmaeus)的β球蛋白基因序列进行比较的完整点阵图。(b)利用滑动窗口对以上的两种球蛋白基因序列进行比较的点阵图,其中窗口大小为10个核苷酸,相似度阈值为8。(a)(b)滑动窗口和阈值的选择过于经验化,信噪比较低,不适合进行高通量的数据分析,对长序列,计算时间长。点阵序列比对的缺点作为双序列比对的第一步,点阵图提供了一个大尺度的轮廓信息序列相似性打分矩阵简单的得分函数:p(a,a)=1p(a,b)=0p(a,-)=p(-,b)=-1例5:不同序列的比对PKDFCKALVP100000000K10001000F0100000T00000K11000A100I00V1PKDFCKALVPK-FTKAIV字符的关系:1.匹配2.删除或插入3.替换PKDFCKALVPK–FTKAIV对于例五中的情况:11-1101101Score=1+1+(-1)+1+0+1+1+0+1=5由于序列长度不同,因此相对长度的得分更有意义:Sim(s,t)=2×Score/(m+n)=2×5/(9+8)=0.588注意:不同类型的字符替换,其代价是不同的。ATCGA5-4-4-4T-45-4-4C-4-45-4G-4-4-45BLAST矩阵目前最流行的序列比较程序BLAST使用的矩阵核酸矩阵一ATCGA1-5-5-1T-51-1-5C-5-11-5G-1-5-51转移矩阵基于颠换、转换原理。核酸矩阵二遗传密码矩阵:通过计算一个氨基酸变化为另一个氨基酸,所需密码子中碱基的变换数目而得到。例如:苯丙氨酸M=色氨酸W(ATG)(TGG)替换值为:2蛋白质矩阵二疏水矩阵:利用氨基酸的疏水性蛋白质矩阵一PAM矩阵:通过统计氨基酸的相互替换率得到的矩阵,最早是由Dayhoff等研究了71个相关蛋白家族的1572个突变。BLOSUM矩阵:也是通过统计氨基酸的相互替换率。蛋白质矩阵三二者区别:前者Dayhoff模型,假设蛋白质序列各部位进化的速率是均等的。但事实上并非如此,因为保守区的进化速率显然低于非保守区。后者由Henikoff算法得到,对不同家族蛋白质序列片段进行比对,不加入gaps,这些序列区间对应于高度保守的区域。氨基酸匹配率可通过各区间可能的匹配率得到。再将这些匹配率计入匹配率表。PAM矩阵(PointAcceptedMutation)序列相似度=14%-27%40%50%60%||||打分矩阵=PAM250PAM120PAM80PAM60BLOSUM:BLOSUM60:序列相似度为60%左右的序列使用;BLOSUM80:序列相似度为80%左右的序列使用。PAM1矩阵,除以10000ARA98672R19913ReplacementaminoacidOriginalaminoacid2.动态规划算法动态规划往往被用于一个复杂的空间中寻找一条最优路径。全局优化比对:Needleman-Wunsch局部优化比对:Smith-Waterman为什么要设计动态规划算法?直接的序列比对,要分别计算所有比对情况的分值,以求得最大分值(或最小分值)。但两序列比对数是序列长度的指数函数,计算量很大,因此必须设计高效的算法。BLOSUM62替代矩阵以两序列VDSCY和VESLCY为例空位罚分d=-11动态规划算法:全局比对(1)GapVDSCYGap01gap2gap…V1gapE2gapS…LCYgdgr)(本例:线性罚分全局比对(2)GapVDSCYGap0-11-22-33-44-55V-11SijE-22S-33L-44C-55Y-66要求解Sij的分数,我们必须先知道Si-1,j-1,Si-1,j,Si,j-1的分数,这种方法叫做递归算法;采用这种方法,可以把大的问题分割成小的问题逐一解决,即动态规划算法;需要存储如何得到Sij分数的过程。全局比对(3)ijGapVDSCYGap0-11-22-33-44-55V-11SijE-22S-33L-44C-55Y-66Needleman-Wunsch算法;Si-1,j-1+σ(xi,yj)Sij=maxofSi-1,j+d(从上到下)Si,j-1+d(从左到右)BLOSUM62替代矩阵全局比对(4)GapVDSCYGap0-11-22-33-44-55V-114E-22S-33L-44C-55Y-664-11-11Needleman-Wunsch算法;Si-1,j-1+σ(xi,yj)Sij=maxofSi-1,j+d(从上到下)Si,j-1+d(从左到右)BLOSUM62替代矩阵以两序列VDSCY和VESLCY为例空位罚分d=-11全局比对(5)GapVDSCYGap0-11-22-33-44-55V-114-7E-22S-33L-44C-55Y-66-3-11-11VD:-3全局比对(6)GapVDSCYGap0-11-22-33-44-55V-114-7-18-29-40E-22-76-5-16-27S-33-18-510-1-12L-44-29-16-19-2C-55-40-27-1287Y-66-51-38-23-315424-1-297比对结果:VDS–CYVESLCYGapVDSCYGap0-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-315基于蛋白质疏水矩阵的全局比对以两序列VDSCY和VESLCY为例空位线性罚分d=-20课后作业结论:1.比对结果与构建的打分矩阵和罚分函数有关2.构建合理的矩阵和罚分函数才能进行最佳比对3.目前较为公认的是BLOSUM和PAM矩阵局部优化比对下例:局部优化打分两条序列如下:LDSCHGESLCK目标:使用局部优化算法寻找最佳比对的结果对全局比对策略稍作修改可得到局部最优比对算法。比对的路径不需要到达搜索图的尽头,如果某种比对的分值不会因为增加比对的数量而增加时,这种比对就是最佳的。依赖于记分系统的性质:因为某种路径的记分会在不匹配的序列段减少,当分值降为零时,路径的延展将会终止,一个新的路径就会产生。局部优化比对(1)GapLDSCHGap000000G0SijE0S0L0C0K0Smith-Waterman算法;Si-1,j-1+σ(xi,yj)Sij=maxofSi-1,j+d(从上到下)Si,j-1+d(从左到右)0gap:-11,线性罚分模型。BLOSUM62替代矩阵局部优化比对(2)GapLDSCHGap000000G00E0S0L0C0K0-11-11-3局部优化比对(3)GapLDSCHGap000000G000E0S0L0C0K0-11-11-4比对结果:GapLDSCHGap000000G000000E002210S002610L040052C001092K000008LDS–CHGESLCK序列比对的分值1.Smith-waterman算法打分:9分2.直接打分:2+4-11+9=4LDS–CHGESLCK比对的统计检验1.用来判断两条序列的比对分数是否足够高,是否有统计意义;2.没有数学理论描述全序列比对的期望分布;3.通常做法:两序列分别打乱进行比对,获得正态分布进行检验;4.Karlin-Altschul公式5.在众多序列比对软件中,均能计算显著性。Karlin-Altschul公式在一定的序列长度m和n限定下,高比值片段对的统计值可由2个参数(k和λ)确定。最简单的形式,即不小于比较值为S的高比值片段对个数,可由下列公式算得其期望值:参数K和λ可分别被简单地视为搜索步长和计分系统的特征数第二节、多序列比对1.不同物种中,许多基因的功能保守,序列相似性较高,通过多条序列的比较,发现保守与变异的部分;2.构建进化的树的必须步骤;3.比较基因组学研