序列分析（一） 一一序列比对

生物信息学研究的三个层面 初级层面： 基于现有的生物信息数据库和资源，利用成熟的生物信息学工具（专业网站、软件）解决生物信息学问题 ——生物信息数据库（NCBI、EBI、DDBJ、UniProt等） ——基因组序列分析、序列比对软件（BLAST、CLUSTAL等） ——系统发育树构建软件的简单使用（PHYLIP、PALM等） ——搜集、整理有特色的生物信息学数据库

中级层面： 利用数理统计方法和相关的工具，研究生物信息学问题 ——现有的数理统计和科学计算工具（EXCEL、SPSS等） ——概率、数理统计基础 ——现有的数理统计和科学计算工具（EXCEL、SPSS等）

高级层面： 提出有重要意义的生物信息学问题；自主创新，发展新方法，开发新工具，引领生物信息学领域研究方向。 ——面向生物学领域，解决重要生物学问题 ——利用数学、物理、化学、计算科学等思想和方法 ——建立模型，发展算法 ——自行编程，开发软件

序列分析内容 与DNA和protein序列相关的研究都可称为序列分析。 主要包括： 1.序列比对 2.基因组序列分析 3.蛋白质序列分析 4.综合序列分析

为什么要进行序列比对？？？ 序列拼接 数据库搜索方面 进化方面 功能方面

与进化相关的几个概念 同源性与相似性：极易被混淆的两个概念！！！ 同源性(homology)：是指序列们是由共同祖先进化而来，两条序列的同源关系。只有两种情况：同源、不同源。 相似性：指序列间的差别，是一个度量。 同源与相似的关系：一般认为序列相似性达到一定程 度，即可认为是同源，但不绝对。

序列的相似性是由 随机因素 产生。 Ortholog（直系同源）： 两个基因通过 物种形成 的事件而产生，或源于不同物种的具有共同祖先的两个基因，或者两个物种中的同一基因，一般具有相同的功能。 Paralog（旁系同源）： 指相同的基因组内因 基因复制 形成的多个具有不同功能的基因。 Xenolog（异同源）： 由某一个 水平基因转移 事件而得到的同源序列。 Convergence（趋同）： 序列的相似性是由 随机因素 产生。

第一节 双序列比对 点阵分析 动态规划

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

例1：自身的比对 A K G F C D E 1

例2：重复序列 A K G F D E 1

例3：反向重复/回文 A U G C 1

例4：RNA stem/loop A U G C 1

例5：不同序列的比对 P K D F C A L V 1 T I PKDFCKALV PK - FTKAIV 字符的关系： 匹配 删除或插入 T I PKDFCKALV PK - FTKAIV 字符的关系： 匹配 删除或插入 替换

点阵法的序列比对 Sequence 1# 1 n 1 “-” Insertion Sequence 2# “-” Insertion m

寻找两条序列的最佳比对，实际上就是寻找矩阵标记图中非重叠平行斜线最长的组合。

→ 序列1 → 实 例 →序列2 → 对于较长的序列，有很多匹配的字符，点阵图变得非常复杂和模糊。

点阵法的滑动窗口技术： 使用滑动窗口代替一次一个位点的比较,是解决这个问题的有效方法。 假设窗口大小为10，相似度阈值为8，则每次比较取10个连续的字符，如相同的字符超过8个，则标记， 基于滑动窗口的点矩阵方法可以明显地降低点阵图的噪声，并且明确无误的指示出了两条序列间具有显著相似性的区域。

(a) (b) （a）对人类（Homo sapiens）与黑猩猩（Pongo pygmaeus）的β球蛋白基因序列进行比较的完整点阵图。（b）利用滑动窗口对以上的两种球蛋白基因序列进行比较的点阵图，其中窗口大小为10个核苷酸，相似度阈值为8。

点阵序列比对的缺点 滑动窗口和阈值的选择过于经验化, 信噪比较低 , 不适合进行高通量的数据分析, 对长序列，计算时间长。 作为双序列比对的第一步，点阵图提供了一个大尺度的轮廓信息

http://www.ebi.ac.uk/Tools/psa/

序列相似性打分矩阵 简单的得分函数： p(a, a)=1 p(a, b)=0 p(a, -)=p(-, b)=-1

例5：不同序列的比对 P K D F C A L V 1 T I PKDFCKALV PK - FTKAIV 字符的关系： 匹配 删除或插入 T I PKDFCKALV PK - FTKAIV 字符的关系： 匹配 删除或插入 替换

对于例五中的情况： PKDFCKALV PK – FTKA I V 1 1 -11 0 1 1 0 1 Score=1+1+(-1)+1+0+1+1+0+1 =5 由于序列长度不同，因此相对长度的得分更有意义： Sim(s,t)=2×Score/(m+n)=2×5/(9+8)=0.588

注意： 不同类型的字符替换，其代价是不同的。

核酸矩阵一 BLAST矩阵 A T C G 5 -4 目前最流行的序列比较程序BLAST使用的矩阵

核酸矩阵二 转移矩阵 A T C G 1 -5 -1 基于颠换、转换原理。

遗传密码矩阵：通过计算一个氨基酸变化为另一个氨基酸，所需密码子中碱基的变换数目而得到。 蛋白质矩阵一 遗传密码矩阵：通过计算一个氨基酸变化为另一个氨基酸，所需密码子中碱基的变换数目而得到。 例如：苯丙氨酸M=>色氨酸W （ATG） （TGG） 替换值为：2 蛋白质矩阵二 疏水矩阵：利用氨基酸的疏水性

蛋白质矩阵三 PAM矩阵：通过统计氨基酸的相互替换率得到的矩阵，最早是由Dayhoff等研究了71个相关蛋白家族的1572个突变。 BLOSUM矩阵：也是通过统计氨基酸的相互替换率。

二者区别： 前者Dayhoff模型，假设蛋白质序列各部位进化的速率是均等的。 但事实上并非如此，因为保守区的进化速率显然低于非保守区。 后者由Henikoff算法得到，对不同家族蛋白质序列片段进行比对，不加入gaps，这些序列区间对应于高度保守的区域。氨基酸匹配率可通过各区间可能的匹配率得到。再将这些匹配率计入匹配率表。

PAM矩阵（Point Accepted Mutation） 序列相似度 = 14% - 27% 40% 50%  60%            |  |  |   | 打分矩阵 = PAM250 PAM120 PAM80 PAM60 BLOSUM： BLOSUM60：序列相似度为 60％左右的序列使用； BLOSUM80：序列相似度为 80％左右的序列使用。

PAM1矩阵，除以10000 A R 9867 2 1 9913 Original amino acid Replacement amino acid PAM1矩阵，除以10000

2. 动态规划算法 动态规划往往被用于一个复杂的空间中寻找一条最优路径。 全局优化比对：Needleman-Wunsch 局部优化比对：Smith-Waterman

为什么要设计动态规划算法? 直接的序列比对，要分别计算所有比对情况的分值，以求得最大分值（或最小分值）。但两序列比对数是序列长度的指数函数，计算量很大，因此必须设计高效的算法。

BLOSUM62替代矩阵 以两序列VDSCY和VESLCY为例 空位罚分d=-11

动态规划算法：全局比对(1) Gap V D S C Y 1gap 2gap … E L 本例：线性罚分

全局比对 (2) Gap V D S C Y -11 -22 -33 -44 -55 Sij E L -66 -11 -22 -33 -44 -55 Sij E L -66 要求解Sij的分数，我们必须先知道 Si-1, j-1, Si-1, j, Si, j-1的分数，这种方法叫做递归算法； 采用这种方法，可以把大的问题分割成小的问题逐一解决，即动态规划算法；需要存储如何得到Sij分数的过程。

全局比对 (3) j i Gap V D S C Y -11 -22 -33 -44 -55 Sij E L -66 -11 -22 -33 -44 -55 Sij E L -66 Needleman-Wunsch算法； Si-1, j-1 + σ(xi, yj) Sij = max of Si-1, j +d (从上到下) Si, j-1 +d (从左到右)

BLOSUM62替代矩阵 以两序列VDSCY和VESLCY为例 空位罚分d=-11

全局比对 (4) Gap V D S C Y -11 -22 -33 -44 -55 4 E L -66 4 -11 -11 -11 -22 -33 -44 -55 4 E L -66 4 -11 -11 Needleman-Wunsch算法； Si-1, j-1 + σ(xi, yj) Sij = max of Si-1, j +d (从上到下) Si, j-1 +d (从左到右)

全局比对 (5) VD: -3 Gap V D S C Y -11 -22 -33 -44 -55 4 -7 E L -66 -3 -11 -11 -22 -33 -44 -55 4 -7 E L -66 -3 -11 -11 VD: -3

全局比对 (6) Gap V D S C Y -11 -22 -33 -44 -55 4 -7 -18 -29 -40 E 6 -5 -16 -11 -22 -33 -44 -55 4 -7 -18 -29 -40 E 6 -5 -16 -27 10 -1 -12 L 9 -2 8 7 -66 -51 -38 -23 -3 15 4 2 4 -1 -2 9 7

比对结果： V D S – C Y V E S L C Y Gap V D S C Y -11 -22 -33 -44 -55 4 -7 -11 -22 -33 -44 -55 4 -7 -18 -29 -40 E 6 -5 -16 -27 10 -1 -12 L 9 -3 8 7 -66 -51 -38 -23 15

基于蛋白质疏水矩阵的全局比对 课后作业 以两序列VDSCY和VESLCY为例 空位线性罚分d=-20

结论： 比对结果与构建的打分矩阵和罚分函数有关 构建合理的矩阵和罚分函数才能进行最佳比对 目前较为公认的是BLOSUM和PAM矩阵

局部优化比对 下例：局部优化打分 两条序列如下： L D S C H G E S L C K 目标：使用局部优化算法寻找最佳比对的结果

对全局比对策略稍作修改可得到局部最优比对算法。 比对的路径不需要到达搜索图的尽头，如果某种比对的分值不会因为增加比对的数量而增加时，这种比对就是最佳的。 依赖于记分系统的性质：因为某种路径的记分会在不匹配的序列段减少，当分值降为零时，路径的延展将会终止，一个新的路径就会产生。

局部优化比对 (1) Gap L D S C H G Sij E K Smith-Waterman算法； G Sij E K Smith-Waterman算法； Si-1, j-1 + σ(xi, yj) Sij = max of Si-1, j +d (从上到下) Si, j-1 +d (从左到右) gap: -11，线性罚分模型。

BLOSUM62替代矩阵

局部优化比对 (2) Gap L D S C H G E K -3 -11 -11

局部优化比对 (3) Gap L D S C H G E K -4 -11 -11

比对结果： L D S – C H G E S L C K Gap L D S C H G E 2 1 6 4 5 9 K 8

序列比对的分值 L D S – C H G E S L C K 1. Smith-waterman算法打分：9分 2. 直接打分：2+4-11+9=4

比对的统计检验 用来判断两条序列的比对分数是否足够高，是否有统计意义； 没有数学理论描述全序列比对的期望分布； 通常做法：两序列分别打乱进行比对，获得正态分布进行检验； Karlin-Altschul公式 在众多序列比对软件中，均能计算显著性。

Karlin-Altschul公式 在一定的序列长度m和n限定下，高比值片段对的统计值可由2个参数(k和λ)确定。最简单的形式，即不小于比较值为S的高比值片段对个数，可由下列公式算得其期望值： 参数k和λ可分别被简单地视为搜索步长和计分系统的特征数

第二节 多序列比对 1. 不同物种中，许多基因的功能保守，序列相似性较高，通过多条序列的比较，发现保守与变异的部分； 第二节 多序列比对 1. 不同物种中，许多基因的功能保守，序列相似性较高，通过多条序列的比较，发现保守与变异的部分； 2. 构建进化树的必需步骤； 3. 比较基因组学研究需要。

动态规划算法：hyperlattice

注 意 最优的多序列比对，其两两序列之间的比对不一定最优。 最优的多序列比对 非最优的双序列比对

多序列比对软件的性能比较 1. ProbCons：目前综合性能最好； 5. ClustalW/X: 最经典、最被广泛接受的工具； 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年，图形化软件

ClustalW/X ClustalW/X对输入序列的格式比较灵活，可以是FASTA格式，还可以是其它格式；用户可以根据自己的需要选择合适的输出格式。 用ClustalW/X得到的多序列比对结果中，所有序列排列在一起，并以特定的符号代表各个位点上残基的保守性，“*”号表示保守性极高的残基位点；“.”号代表保守性略低的残基位点。 ClustalX在windows下运行，而ClustalW可在dos下运行

ClustalX: 使用指南 EBI的Clustal网址是：http://www.ebi.ac.uk/Tools/msa/clustalw2/ 下载ClustalW/X的网址是：ftp://ftp.ebi.ac.uk/pub/software/clustalw2 需要安装 详细指南： http://akira.ruc.dk/~olesk/sekvens/Treedraw.htm/ http://www.dbbm.fiocruz.br/james/ClustalX_tutorial.html

导入序列文件 Load Sequences 注意：序列文件不可放在中文名命的文件夹下

以conotoxin O-superfamily Delta-family 为例 调节字符的大小

执行比对 Do Complete Alignment

文件导出 可被其它软件(bioedit)用来进一步分析

多序列比对结果 保守位点

Parameter selection

可利用自己定义的矩阵进行连配

结果输出

Muscle使用 下载地址： http://www.drive5.com/muscle/ 解压即可使用 自带使用说明

命令界面

基本命令： muscle –in <inputfile> -out <outputfile>

muscle –in <inputfile> -html -out <outputfile>