序列分析(一)一一序列比对

1、生物信息学研究的三个层面生物信息学研究的三个层面初级层面：初级层面： 基于现有的生物信息数据库和资源，利用成熟基于现有的生物信息数据库和资源，利用成熟的生物信息学工具（专业网站、软件）解决生物信的生物信息学工具（专业网站、软件）解决生物信息学问题息学问题 生物信息数据库（生物信息数据库（NCBINCBI、EBIEBI、DDBJDDBJ、UniProtUniProt等）等） 基因组序列分析、序列比对软件（基因组序列分析、序列比对软件（BLASTBLAST、CLUSTALCLUSTAL等）等） 系统发育树构建软件的简单使用（系统发育树构建软件的简单使用（PHYLIPPHYLIP、PALMPALM等

2、）等） 搜集、整理有特色的生物信息学数据库搜集、整理有特色的生物信息学数据库中级层面：中级层面： 利用数理统计方法和相关的工具，研究生物信息利用数理统计方法和相关的工具，研究生物信息学问题学问题 概率、数理统计基础概率、数理统计基础 现有的数理统计和科学计算工具（现有的数理统计和科学计算工具（EXCELEXCEL、SPSSSPSS等）等）高级层面：高级层面： 提出有重要意义的生物信息学问题；自主创新，提出有重要意义的生物信息学问题；自主创新，发展新方法，开发新工具，引领生物信息学领域研发展新方法，开发新工具，引领生物信息学领域研究方向。究方向。 面向生物学领域，解决重要生物学问题面向生物学领域

3、，解决重要生物学问题 利用数学、物理、化学、计算科学等思想和方法利用数学、物理、化学、计算科学等思想和方法 建立模型，发展算法建立模型，发展算法 自行编程，开发软件自行编程，开发软件序列分析内容序列分析内容 与与DNADNA和和proteinprotein序列相关的研究都可称为序列相关的研究都可称为序列分析。序列分析。 主要包括：主要包括： 1.1.序列比对序列比对 2.2.基因组序列分析基因组序列分析 3.3.蛋白质序列分析蛋白质序列分析 4.4.综合序列分析综合序列分析为什么要进行序列比对？为什么要进行序列比对？序列拼接序列拼接数据库搜索方面数据库搜索方面进化方面进化方面功能方面功能方面与

4、进化相关的几个概念与进化相关的几个概念同源性同源性(homology)： 是指序列们是由共同祖先进化而来，讲两条序列是指序列们是由共同祖先进化而来，讲两条序列的同源关系，只有两种情况：的同源关系，只有两种情况：同源、不同源同源、不同源。相似性：相似性：指序列间的差别，是一个度量。指序列间的差别，是一个度量。 同源与相似的关系：一般认为序列相似性达到一同源与相似的关系：一般认为序列相似性达到一定程度，即可认为是同源，但不绝对。定程度，即可认为是同源，但不绝对。Ortholog (直系同源直系同源)： 两个基因通过两个基因通过 物种形成物种形成 的事件而产生，或源于不同物的事件而产生，或源于不同物

5、种的具有共同祖先的两个基因，或者两个物种中的同一基因，种的具有共同祖先的两个基因，或者两个物种中的同一基因，一般具有相同的功能。一般具有相同的功能。Paralog (旁系同源旁系同源)： 指相同的基因组内因指相同的基因组内因 基因复制基因复制 形成的多个具有不同功形成的多个具有不同功能的基因。能的基因。Xenolog (异同源异同源)： 由某一个由某一个 水平基因转移水平基因转移 事件而得到的同源序列。事件而得到的同源序列。Convergence(趋同趋同)： 序列的相似性是由序列的相似性是由 随机因素随机因素 产生。产生。第一节、双序列比对第一节、双序列比对点阵分析点阵分析动态规划动态规划1

6、. 点阵分析点阵分析用途：用途： 1. 寻找两条序列间所有可能的比对；寻找两条序列间所有可能的比对； 2. 寻找蛋白质、寻找蛋白质、DNA序列上正向或反向的重复序列；序列上正向或反向的重复序列； 3. 发现发现RNA上可能存在的互补区域。上可能存在的互补区域。优点：优点： 1. 可以找到两个序列间所有可能的残基匹配；可以找到两个序列间所有可能的残基匹配； 2. 简单、易懂简单、易懂 3. 直观、整体性强直观、整体性强工具：工具： http:/www.ebi.ac.uk/Tools/emboss/index.html例例1：自身的比对：自身的比对A K G F K C A D EA 100000

7、100K10010000G1000000F100000K10000C1000A100D10E1例例2：重复序列：重复序列 A K G F D K G F EA 100000000K10001000G1000100F100010D10000K11000G1100F110E1例例3 3：反向重复：反向重复/ /回文回文A U G C A C G U CA100010000U10000010G1000100C101000A10000C11001G1100U110C1例例4：RNA stem/loopA U G U A G C A UA100010010U10100001G1001000C000001

8、U00001A0010C100A10U1例例5：不同序列的比对：不同序列的比对 P K D F C K A L VP100000000K10001000F0100000T00000K11000A100I00V1PKDFCKALVPK - FTKAIV字符的关系：字符的关系： 匹配匹配 删除或插入删除或插入1.替换替换点阵法的序列比对点阵法的序列比对 Sequence 1#1nSequence 2#1m“-” Insertion“-” Insertion寻找两条序列的最佳比对，实际寻找两条序列的最佳比对，实际上就是寻找在矩阵标记图中找非上就是寻找在矩阵标记图中找非重叠平行斜线最长的组合。重叠平行

9、斜线最长的组合。 序列序列1 1 序列序列2 2 实实 例例对于较长的对于较长的序列，有很序列，有很多匹配的字多匹配的字符，点阵图符，点阵图变得非常复变得非常复杂和模糊。杂和模糊。点阵法的滑动窗口技术：点阵法的滑动窗口技术： 使用滑动窗口代替一次一个位点的比较使用滑动窗口代替一次一个位点的比较, ,是是解决这个问题的有效方法。解决这个问题的有效方法。 假设窗口大小为假设窗口大小为1010，相似度阈值为，相似度阈值为8 8，则每，则每次比较取次比较取1010个连续的字符，如相同的字符超个连续的字符，如相同的字符超过过8 8个，则标记，个，则标记， 基于滑动窗口的点矩阵方法可以明显地降低基于滑动窗

10、口的点矩阵方法可以明显地降低点阵图的噪声，并且明确无误的指示出了两点阵图的噪声，并且明确无误的指示出了两条序列间具有显著相似性的区域。条序列间具有显著相似性的区域。 （a）对人类（）对人类（Homo sapiens）与黑猩猩（）与黑猩猩（Pongo pygmaeus）的）的球蛋白基因序列进行比较的完整点阵图。（球蛋白基因序列进行比较的完整点阵图。（b）利用滑动窗口对）利用滑动窗口对以上的两种球蛋白基因序列进行比较的点阵图，其中窗口大小为以上的两种球蛋白基因序列进行比较的点阵图，其中窗口大小为10个核苷酸，相似度阈值为个核苷酸，相似度阈值为8。 (a) (b) 滑动窗口和阈值的选择过于经验化滑动

11、窗口和阈值的选择过于经验化, , 信噪比较低信噪比较低 , , 不适合进行高通量的数据分析不适合进行高通量的数据分析, , 对长序列，计算时间长。对长序列，计算时间长。点阵序列比对的缺点点阵序列比对的缺点作为双序列比对的第一步，点阵图提供了作为双序列比对的第一步，点阵图提供了一个大尺度的轮廓信息一个大尺度的轮廓信息 http:/www.ebi.ac.uk/Tools/emboss/align/序列相似性打分矩阵序列相似性打分矩阵简单的得分函数：简单的得分函数： p(a, a)=1 p(a, b)=0 p(a, -)=p(-, b)=-1例例5：不同序列的比对：不同序列的比对 P K D F C

12、 K A L VP100000000K10001000F0100000T00000K11000A100I00V1PKDFCKALVPK - FTKAIV字符的关系：字符的关系： 匹配匹配 删除或插入删除或插入1.替换替换PKDFCKALVPK FTKA I V对于例五中的情况：对于例五中的情况：1 1 -11 0 1 1 0 1Score=1+1+(-1)+1+0+1+1+0+1 =5由于序列长度不同，因此相对长度的得分更有意义：由于序列长度不同，因此相对长度的得分更有意义： Sim(s,t)=2Score/(m+n)=25/(9+8)=0.588注意：注意： 不同类型的字符替换，其不同类型的

13、字符替换，其代价是不同的。代价是不同的。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 （A

14、TG） （TGG） 替换值为：替换值为：2蛋白质矩阵二蛋白质矩阵二疏水矩阵：疏水矩阵：利用氨基酸的疏水性利用氨基酸的疏水性蛋白质矩阵一蛋白质矩阵一PAM矩阵：矩阵：通过统计氨基酸的相互替换通过统计氨基酸的相互替换率得到的矩阵，最早是由率得到的矩阵，最早是由Dayhoff等研究了等研究了71个相关蛋白家族的个相关蛋白家族的1572个突变。个突变。BLOSUM矩阵：矩阵：也是通过统计氨基酸也是通过统计氨基酸的相互替换率。的相互替换率。蛋白质矩阵三蛋白质矩阵三二者区别：二者区别：前者前者Dayhoff模型，假设蛋白质序列各部位进化的模型，假设蛋白质序列各部位进化的速率是均等的。速率是均等的。但事实上

15、并非如此，因为保守区的进化速率但事实上并非如此，因为保守区的进化速率 显然显然低于非保守区。低于非保守区。后者由后者由Henikoff算法得到，对不同家族蛋白质序算法得到，对不同家族蛋白质序列片段进行比对，不加入列片段进行比对，不加入gaps，这些序列区间对，这些序列区间对应于高度保守的区域。应于高度保守的区域。 氨基酸匹配率可通过各区氨基酸匹配率可通过各区间可能的匹配率得到。再将这间可能的匹配率得到。再将这 些匹配率计入匹配些匹配率计入匹配率表。率表。PAM矩阵（矩阵（Point Accepted Mutation）序列相似度序列相似度 = 14% - 27% 40% 50% 60% | |

16、 | |打分矩阵打分矩阵 = PAM250 PAM120 PAM80 PAM60BLOSUM：BLOSUM60：序列相似度为：序列相似度为 60左右的序列使用；左右的序列使用；BLOSUM80：序列相似度为：序列相似度为 80左右的序列使用。左右的序列使用。 PAM1矩阵，除以矩阵，除以10000AA9867R19913Replacement amino acidOriginal amino acid2. 动态规划算法动态规划算法 动态规划往往被用于一个复杂的空间中寻找动态规划往往被用于一个复杂的空间中寻找一条最优路径。一条最优路径。 全局优化比对：全局优化比对：Needleman-Wunsc

17、h 局部优化比对：局部优化比对：Smith-Waterman为什么要设计动态规划算法为什么要设计动态规划算法? 直接的序列比对，要分别计算所有直接的序列比对，要分别计算所有比对情况的分值，以求得最大分值（或比对情况的分值，以求得最大分值（或最小分值）。最小分值）。但两序列比对数是序列长但两序列比对数是序列长度的指数函数，计算量很大，因此必须度的指数函数，计算量很大，因此必须设计高效的算法设计高效的算法。BLOSUM62替代矩阵替代矩阵以两序列以两序列VDSCY和和VESLCY为例为例空位罚分空位罚分d=-11动态规划算法：全局比对动态规划算法：全局比对(1)GapVDSCYGap01gap2g

18、apV1gapE2gapSLCYgdgr)(本例：线性罚分本例：线性罚分全局比对全局比对 (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)

19、ijGapVDSCYGap0-11-22-33-44-55V-11SijE-22S-33L-44C-55Y-66Needleman-Wunsch算法；算法； Si-1, j-1 + (xi, yj)Sij = max of Si-1, j +d (从上到下从上到下) Si, j-1 +d (从左到右从左到右)BLOSUM62替代矩阵替代矩阵全局比对全局比对 (4)(4)GapVDSCYGap0-11-22-33-44-55V-114E-22S-33L-44C-55Y-664-11-11Needleman-Wunsch算法；算法； Si-1, j-1 + (xi, yj)Sij = max of

20、 Si-1, j +d (从上到下从上到下) Si, j-1 +d (从左到右从左到右)BLOSUM62替代矩阵替代矩阵以两序列以两序列VDSCY和和VESLCY为例为例空位罚分空位罚分d=-11全局比对全局比对 (5)(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-6

21、6-51-38-23-315424-1-297比对结果：比对结果：V D S C YV E S L C YGapVDSCYGap0-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课后作业课后作业结论：结论： 比对结果与构建的打分矩阵和罚分比对结果与构建的打分矩阵和罚分函数有关函数有关 构建

22、合理的矩阵和罚分函数才能进构建合理的矩阵和罚分函数才能进行最佳比对行最佳比对1.目前较为公认的是目前较为公认的是BLOSUM和和PAM矩阵矩阵局部优化比对局部优化比对下例：局部优化打分下例：局部优化打分两条序列如下：两条序列如下：L D S C HG E S L C K目标：使用局部优化算法寻找最佳比对的目标：使用局部优化算法寻找最佳比对的结果结果对全局比对策略稍作修改可得到局部最对全局比对策略稍作修改可得到局部最优比对算法。优比对算法。比对的路径不需要到达搜索图的尽头比对的路径不需要到达搜索图的尽头 ,如如果某种比对的分值不会因为增加比对的果某种比对的分值不会因为增加比对的数量而增加时，这种

23、比对就是最佳的。数量而增加时，这种比对就是最佳的。依赖于记分系统的性质：因为某种路径依赖于记分系统的性质：因为某种路径的记分会在不匹配的序列段减少的记分会在不匹配的序列段减少 ,当分值当分值降为零时，路径的延展将会终止，一个降为零时，路径的延展将会终止，一个新的路径就会产生。新的路径就会产生。 局部优化比对局部优化比对 (1)GapLDSCHGap000000G0SijE0S0L0C0K0Smith-Waterman算法；算法； Si-1, j-1 + (xi, yj)Sij = max of Si-1, j +d (从上到下从上到下) Si, j-1 +d (从左到右从左到右) 0gap:

24、-11，线性罚分模型。，线性罚分模型。BLOSUM62替代矩阵替代矩阵局部优化比对局部优化比对 (2)GapLDSCHGap000000G00E0S0L0C0K0-11-11-3局部优化比对局部优化比对 (3)GapLDSCHGap000000G000E0S0L0C0K0-11-11-4比对结果：比对结果：GapLDSCHGap000000G000000E002210S002610L040052C001092K000008L D S C HG E S L C K序列比对的分值序列比对的分值1. Smith-waterman算法打分：算法打分：9分分2. 直接打分：直接打分：2+4-11+9=4

25、L D S C HG E S L C K比对的统计检验比对的统计检验用来判断两条序列的比对分数是否足够高，用来判断两条序列的比对分数是否足够高，是否有统计意义；是否有统计意义；没有数学理论描述全序列比对的期望分布；没有数学理论描述全序列比对的期望分布；通常做法：两序列分别打乱进行比对，获通常做法：两序列分别打乱进行比对，获得正态分布进行检验；得正态分布进行检验；Karlin-Altschul公式公式在众多序列比对软件中，均能计算显著性。在众多序列比对软件中，均能计算显著性。Karlin-Altschul公式公式在一定的序列长度在一定的序列长度m和和n限定下，高比值片段对的限定下，高比值片段对的

26、统计值可由统计值可由2个参数个参数(k和和)确定。最简单的形式，即确定。最简单的形式，即不小于比较值为不小于比较值为S的高比值片段对个数，可由下列的高比值片段对个数，可由下列公式算得其期望值：公式算得其期望值： 参数参数K和和可分别被简单地视为搜索步长和计分系统可分别被简单地视为搜索步长和计分系统的特征数的特征数第二节、多序列比对第二节、多序列比对1. 不同物种中，许多基因的功能保守，序列不同物种中，许多基因的功能保守，序列相似性较高，通过多条序列的比较，发现相似性较高，通过多条序列的比较，发现保守与变异的部分；保守与变异的部分；2. 构建进化的树的必须步骤；构建进化的树的必须步骤；3. 比较

27、基因组学研究需要。比较基因组学研究需要。动态规划算法：动态规划算法：hyperlattice注意注意最优的多序列比对，其两两序列之间的比对不最优的多序列比对，其两两序列之间的比对不一定最优。一定最优。 最优的多序列比对最优的多序列比对非最优的双序列比对非最优的双序列比对多序列比对软件的性能比较多序列比对软件的性能比较1. ProbCons：目前综合性能最好；：目前综合性能最好；2. T-Coffee：序列相似性高时最准确；：序列相似性高时最准确；3. DIALIGN: 序列相似性低时最准确；序列相似性低时最准确；4. POA：性能接近：性能接近T-Coffee和和DIALIGN，速度，速度最快；最快；5. ClustalW/X: 最经典、最被广泛接受最经典、最被广泛接受的工具；的工具；6. MUSCLE: 目前最流行的多序列比对工具。目前最流行的多序列比对工具。ClustalW/X发展历史：发展历史：1. Clustal: 1988年开发；年开发；2. ClustalW: 1994年，年，Julie D. Thompson等人等人改进、开发；改进、开发；3. ClustalX: 1997年，图形化软件年，图形化