主讲人：孙 啸 制作人：刘志华

1. 第五章 DNA序列分析 主讲人：孙 啸 制作人：刘志华 东南大学 吴健雄实验室

2. 第五章 DNA序列分析 DNA序列分析 ——基因序列 ——基因表达调控信息 寻找基因牵涉到两个方面的工作 ： • 识别与基因相关的特殊序列信号 • 预测基因的编码区域 • 结合两个方面的结果确定基因的位置和结构 基因表达调控信息隐藏在基因的上游区域，在组成上具有一定的特征，可以通过序列分析识别这些特征。

3. 第一节 DNA序列分析步骤和分析结果评价 • 在DNA序列中，除了基因之外，还包含许多其它信息，这些信息大部分与核酸的结构特征相关联，通常决定了DNA与蛋白质或者DNA与RNA的相互作用。 • 存放这些信息的DNA片段称为功能位点 • 如启动子（Promoter）、基因终止序列（Terminator sequence）、剪切位点（Splice site）等。

4. 一个基本的DNA序列分析方案 发现重复元素 数据库搜索 分析功能位点 序列组成统计分析 综合分析

5. 功能序列分析的准确性来自于对“功能序列”和“非功能序列”的辨别能力。功能序列分析的准确性来自于对“功能序列”和“非功能序列”的辨别能力。 • 两个集合： 训练集（training set） • 用于建立完成识别任务的数学模型。 测试集或控制集（control set） • 用于检验所建模型的正确性。 • 用训练集中实例对预测模型进行训练，使之通过学习后具有正确处理和辨别能力。然后，用模型对测试集中的实例进行“功能”与“非功能”的判断，根据判断结果计算模识别的准确性。

6. 收集已知的功能序列和非功能序列实例 （这些序列之间是非相关的 ） 训练集 （training set） 测试集或控制集 （control set） 建立完成识别任务的模型 检验所建模型的正确性 对预测模型进行训练， 使之通过学习后具有 正确处理和辨别能力。 • 进行“功能”与“非功能”的 • 判断，根据判断结果计算 • 模识别的准确性。 识别“功能序列”和“非功能序列”的过程

7. Sn ——敏感性 • Sp——特异性 • Tp是正确识别的功能序列数， • Tn为正确识别的非功能序列数， • Fn是被错误识别为非功能序列的功能序列数， • Fp是被错误识别为功能序列的非功能序列数。

8. 敏感性和特异性的权衡 • 对于一个实用程序，既要求有较高的敏感性，也要求有较高的特异性。 • 如果敏感性很高，但特异性比较低，则在实际应用中会产生高比率的假阳性； • 相反，如果特异性很高，而敏感性比较低，则会产生高比率的假阴性。 • 对于敏感性和特异性需要进行权衡，给出综合评价指标。

9. 对于一个识别程序准确性可按下式进行综合评价：对于一个识别程序准确性可按下式进行综合评价： • 另一个综合评介指标为相关系数，其计算计算公式为：

10. 选择训练集和测试集 • 在检测算法的可行性时，需要从已知的数据中按照不同的方式选择训练集和测试集 • 测试集的构成非常关键 • 在不同的测试集上进行测试可能会得到不同的准确性结果，甚至准确性相差很大。 • 建立标准的功能序列测试集合。 • 如基因转录剪切位点的测试集合、编码区域的测试集合等。

11. 第二节 核苷酸关联分析 • 对于一个给定的基因组，最简单的计算就是统计DNA序列中各类核苷酸出现的频率。 • 对于随机分布的DNA序列，每种核苷酸的出现是均匀分布的 • 出现频率各为0.25。 • 而真实基因组的核苷酸分布则是非均匀的

12. 酵母基因组核苷酸出现频率

13. 在统计过程中，如果同时计算DNA的正反两条链，则根据碱基配对原则，A和T、C和G的出现频率相同。在统计过程中，如果同时计算DNA的正反两条链，则根据碱基配对原则，A和T、C和G的出现频率相同。 • 如果仅统计一条链，则虽然A和T、C和G的出现频率不同，但是非常接近。

14. 单链核苷酸出现频率

15. 基因和其它功能区域 在正反两条链上出现的 可能性通常一样 正反两条链在信息的 组织结构方面不应该有差别 核苷酸出现频率也不应该 有偏差 正反两条链碱基互补的原则 单链上A和T、C和G的 出现频率相近。 单链上A和T、C和G的出现频率相近的解释

16. 两联核苷酸频率 • 不同基因组中两个连续核苷酸出现的频率也是不相同的 • 4种核苷酸可以组合成16种两联核苷酸

17. 对酵母基因组两联核苷酸的统计结果 其中核苷酸对出现频率最高的达到0.119 而出现频率最低的只有0.028 酵母基因组两联核苷酸频率表

18. 关联性分析 令: Pij —— 代表两联核苷酸（i，j）的出现频率 Pi —— 代表核苷酸i的出现频率 则 ： Pij’= Pij/(PiPj) 的值反应核苷酸i和j的关联关系 如果Pij’=1，则在两个连续的位置上，核苷酸i和j的出现是相对独立的。

19. 关联性分析 • 对于酵母基因组 PA=0.3248 PAA=0.1193 PAA’ =0.1193/（0.3248*0.3248） =1.131 > 1 表明在两个连续位置上“A”的出现不是独立的，而是相关的。

20. 同样，对于相隔一定距离k（k代表核苷酸个数）的两个核苷酸，也可能具有一定的相关性。同样，对于相隔一定距离k（k代表核苷酸个数）的两个核苷酸，也可能具有一定的相关性。 • 假设Pij(k)代表核苷酸j出现在核苷酸i之后第k个位置的频率，则可定义一个反应统计相关性的互信息I(k) • I(k)值得大小实际上反应了距离为k的两个核苷酸之间的相关性的程度

21. 三联核苷酸——基因密码子 • 在进行编码区域识别时，常常需要对三联核苷酸进行统计分析，这实际上是分析密码子的使用偏性。 • 由于密码子的简并性（degeneracy），每个氨基酸至少对应1种密码子，最多有6种对应的密码子。 • 在基因中，同义密码子的使用并不是完全一致的。 • 不同物种、不同生物体的基因密码子使用存在着很大的差异

22. 基因密码子的使用与基因编码的蛋白的结构和功能有关，与基因表达的生理功能有着密切的联系 • 蛋白的三级结构与密码子使用概率有密切的关系 • 通过对密码子的聚类分析，可以很清晰地将具有不同三级结构蛋白质的编码基因分成不同的类，而具有相似三级结构蛋白的编码基因则大致聚在同一类中，从而证明基因密码子的使用偏性与蛋白质三级结构具有密切的相关性。 • 在不同物种中，类型相同的基因具有相近的同义密码子使用偏性 • 对于同一类型的基因由物种引起的同义密码子使用偏性的差异较小

23. 针对酵母第一染色体的分析结果

24. 第三节 功能位点分析 • 功能位点（functional site） • 与特定功能相关的位点，是生物分子序列上的一个功能单元，或者是生物分子序列上一个较短的片段。 • 功能位点又称为功能序列（functional sequence）、序列模式（motif）、信号（signal）等。 • 核酸序列中的功能位点包括转录因子结合位点、转录剪切位点、翻译起始位点等。 • 在蛋白质序列分析中，常使用序列模式这个名词，蛋白质的序列模式往往与蛋白质结构域或者作用部位有关。

25. 功能位点示意

26. 基因组序列中若干个相邻的功能位点组合形成功能区域（functional region）。 • 功能位点分析的任务 • 发现功能位点特征 • 识别功能位点

27. 1、利用共有序列搜索功能位点 • 共有序列（consensus）又称一致性片段 • 共有序列是关于功能位点特征的描述，它描述了功能位点每个位置上核苷酸进化的保守性 例如: NTATN • 利用共有序列进行功能位点分析牵涉到两个方面的问题， • 如何构造共有序列 • 如何利用共有序列在给定的核酸序列上搜索寻找功能位点，并计算所找到的功能位点的可靠性

28. 共有序列具有以下几个方面的特征： （1）共有序列中既有保守的位置，也有可变的位置； （2）任何位置上的核苷酸可以用15种类型之一来表示：

29. 核苷酸表示符号

30. 共有序列构造过程: (1)初始化共有序列为一系列可变位置，以“N”代表； (2)在可变位置寻找出现次数最多的核苷酸，并将该位置转化为保守位置； (3)对当前所得到的共有序列进行特异性检查，若通过检查，转（5），否则转（4）； (4)形成与当前共有序列一致的位点子集，转（2）； (5)从原位点集合中删除与当前共有序列一致的位点，若还有剩余位点，则转（1），构造另外的共有序列。

31. [1] [2] [3] [4] [2] [3] TTATG TTATG TNNNN tTATG tTATG TNNNC ATATA ATATA 非特异 tACGC tACGC 非特异 TACGC TACGC tTGTC tTGTC TTGTC TTGTC tCCAC tCCAC TCCAC TCCAC T N N N C N N N N N T N N N N [4] [2] [3] [5] [5] tACGc tACGc TNSNC Consensus1： Consensus2： tTGTc tTGTc 特异 TNSNC NTATN tCCAc tCCAc 剩余位点： TTATG S T N N C ATATA

32. 在给定的序列中搜索与共有序列一致的序列片段在给定的序列中搜索与共有序列一致的序列片段 • 数据库搜索 • 共有序列表示方法的缺点： • 是关于序列特征的一种定性描述，对于DNA序列，它能够说明序列每个位置可能出现的碱基类型，但是不能准确地说明各位置上不同类型碱基出现的可能性大小。

33. 2、用感知矩阵分析功能位点 • 用权系数描述功能位点各位置上每种核苷酸的相对重要性 • 感知矩阵（或加权矩阵） • 根据一系列功能位点的多重对比排列结果而建立的 • 其大小为4n 4代表碱基的种类数目，n代表功能位点的长度

34. 矩阵的每一个元素M(a,j)的值代表第a种核苷酸在功能位点第j个位置上出现的得分,a {A,T,G,C}。 感知矩阵示例

35. 对于一个序列s=a1a2…an，根据对应位置上核苷酸的类型，取感知矩阵中对应的权值，加和以后得到该序列的得分对于一个序列s=a1a2…an，根据对应位置上核苷酸的类型，取感知矩阵中对应的权值，加和以后得到该序列的得分 • 设S=ATTGCA，则 Ws = 1+6+14-5+8+19=43 • T——功能位点阈值 • T‘——非功能位点阈值 • 如果Ws T，则S是功能位点； • 如果WsT'，则S是非功能位点。

36. 感知矩阵M的构造算法 令A+代表功能位点集合 A-代表非功能位点集合 • 过程如下： （1）初始化M为零矩阵； （2）执行过程（3）-（6）的循环； （3）逐步取训练集合中的每个实例Si，如果Si A+，转 过程（4）；如果Si A-，转过程（5）； （4）如果W（Si）T，M不变，否则根据Si的核苷酸分布 将M中所有对应元素的值加1；转（6）； （5）如果W（Si）T‘，M不变，否则根据Si的核苷酸分 布将M中所有对应元素的值减1；转（6）； （6）若训练集合中的所有实例都处理过，则循环结束， 转（7），否则继续执行循环体，直到处理完所有实 例； （7）如果M稳定，则结束；否则转（2）。

37. 上述算法反复调整感知矩阵M的元素值，直到M矩阵能够正确识别训练集中的所有功能位点和非功能位点。上述算法反复调整感知矩阵M的元素值，直到M矩阵能够正确识别训练集中的所有功能位点和非功能位点。 • 对于最终得到的感知矩阵，要求其具有敏感性和特异性，每一列上的元素值应该尽可能地有明显的差别，以便反应功能位点各个位置上的特点。

38. 与感知矩阵类似，如果令矩阵每一个元素M(a,j)的值代表第a种核苷酸在功能位点第j个位置上出现的概率，则M是一个概率矩阵。与感知矩阵类似，如果令矩阵每一个元素M(a,j)的值代表第a种核苷酸在功能位点第j个位置上出现的概率，则M是一个概率矩阵。 • 假设各个位置上出现的碱基是相互独立的，即任何两个位置上的碱基是不相关的，那么对于给定一个序列s=a1a2…an，可以计算出功能位点序列为s的概率：

39. 如果分别统计功能位点和非功能位点，通过计算可以形成两个矩阵M和M’，进一步计算可以判断一个给定的序列究竟属于功能位点，还是属于非功能位点。给定一个序列s=a1a2…an，定义似然比LR(M,M’,s)：如果分别统计功能位点和非功能位点，通过计算可以形成两个矩阵M和M’，进一步计算可以判断一个给定的序列究竟属于功能位点，还是属于非功能位点。给定一个序列s=a1a2…an，定义似然比LR(M,M’,s)： • 在进行功能位点检测时，计算LR(M,M’,s)，并与给定的阈值L比较，如果LR(M,M’,s)>L，则序列s可能是一个功能位点。

40. 概率矩阵M和M’的每个元素是一个0和1之间的正数。概率矩阵M和M’的每个元素是一个0和1之间的正数。 • 如果令一个4n新矩阵U的元素(a,j)的值为 • log2(M(a,j)/M’(a,j)) • 则矩阵U的每个元素值可能是正值，也可能是负值。 实际上，矩阵U就是感知矩阵。

41. 第四节 隐马尔柯夫模型 1、马尔柯夫链（Markov chain） • 考虑一个具有多个状态的系统S，S={s1,s2,…,s|s|},令S0、S1、…、St为一系列在各个时刻系统状态的变量，即状态链。 • 对于每个1到|S|的整数，它们分别与状态链中的一个状态相联系，并且在任何时刻，这条链都处于一个特殊的状态。当且仅当对于任何t有 则St形成一条马尔柯夫链。

42. 简单地说，就是系统未来的状态仅依赖于当前状态。St称为在时刻t系统链的状态。一条马尔柯夫链完全决定于初始分布P(S0)和转换概率Pt=P(St+1|St)。 • 令状态转换矩阵为 F =(fij) fij代表从状态si移动到状态sj的概率。

43. 生物序列可以被描述为一个随机过程的输出，其中对于一个给定的核酸在位置p出现的概率依赖于已占据前面k个位置的核酸，这样一种表示称为k阶马尔柯夫模型。生物序列可以被描述为一个随机过程的输出，其中对于一个给定的核酸在位置p出现的概率依赖于已占据前面k个位置的核酸，这样一种表示称为k阶马尔柯夫模型。 ATCGTAGCAT…….

44. 一个序列具有不同的统计性质 （如二目频率或三目周期性） • 不同的功能区域（如编码区域、非编码区域）对应于不同的马尔柯夫模型。

45. 马尔柯夫链在识别CpG岛中的应用 • CpG岛是一类长度在几百bp的特殊DNA序列，其中CG核苷酸对出现的频率非常高。 ACGCGCGTACGCGAAT • CpG岛在基因组中有重要的生物学意义，而识别CpG岛有助于在基因组序列中确定我们感兴趣的区域。

46. CpG岛的识别问题表述为：给定一段DNA序列 X =(x1，x2，… ,xL)，确定X是否是一个CpG岛。 • 设字母表A={a,t,c,g}，对于字母表中的任何两个字符s、t，定义转换概率为fst=p(xi=t|xi-1=s)，即字符s后面出现字符t的概率。 • 假设{xi}是一个随机过程，随机变量xi的取值仅依赖于xi-1，即对于所有x1,x2,…,xiA，

47. 整个序列X的发生概率为 • 为了处理方便，添加两个特殊的字符‘B’（begin）和‘E’（end），使得x0=‘B’, xL+1=‘E’，则上述公式简化为：

48. 令fst+为CpG 岛内的字符转换概率 fst-为CpG 岛外的字符转换概率 则X的对数似然得分为 • 上述计算值越大，则X越可能是CpG岛。

49. CpG岛内部和外部的转换概率

50. 另外一个待解决的问题是： 给定DNA序列，确定CpG岛的位置。 直接的方法： 对窗口内的子序列计算得分Score(Xk)， 具有正值的Xk就是可能的CpG岛 子序列起始位置为k+1 ，长度为l