数学建模_DNA序列分类模型_(终稿).pdf

DNADNA 序列分类模型序列分类模型 摘要摘要 本文分析了已知类别的人工 DNA 序列的特征 建立了聚类分析延拓模型和 马尔可夫模型 分别对未知类别的人工 DNA 序列和自然序列进行分类 根据分 类效果选出了较优模型 首先对数据进行预处理 得到人工 DNA 序列的单个碱基丰度和不同碱基丰度 之比等特征量 进而分析 A B 两类的差异 得到合适的特征判定条件对未知类 别的 DNA 序列进行分类 计算人工 DNA 序列的特征量 给出各序列的统计数据 其次用聚类分析延拓模型进行分类 用 A B 两类具有明显差异的特征作为 样品特征变量 得到欧式空间中表征编号 1 20 人工 DNA 序列的特征向量 计算 两两之间的 Lance 和 Williams 距离进行相似性度量 逐步选择相似性较大的归为 一类 同时不断更新类内的标准比较特征向量 对聚类方法进行延拓 最终得到 类内差异小 类间差异大的 A B 两类 建立了聚类分析延拓模型 再对选取的 特征变量进行改进 提高模型的分类效果 最后 借助均值 方差和相关系数等 参数对改进模型的分类效果进行分析 再次用马尔可夫模型进行分类 将 DNA 序列看成是马尔可夫链 求出编号 1 10 和 11 20 人工 DNA 序列在已知当前碱基种类的条件下 下一个碱基出现任 一种的概率 结果存入概率转移矩阵 1 和 2 再利用矩阵 1 和 2 分别求出编号 1 20 中任一条 DNA 序列出现的概率 选择较大的一个作为该 DNA 序列的分类 建立 马尔可夫模型 再进行与聚类分析延拓模型类似的改进和检验工作 然后对编号 21 40 人工 DNA 序列和 182 条自然序列进行分类 得到最终结果 最后 用层次分析法综合评价模型一与模型二 选择聚类分析延拓模型作为 最终模型 其分类结果作为最终结果 具体如下 编号 21 40 人工 DNA 序列中属于 A 类的样品编号为 22 23 25 27 29 30 34 35 36 37 39 属于 B 类的样品编号为 21 24 26 28 31 32 33 38 40 182 条自然序列中 属于 B 类的样品编号为 7 10 12 22 23 24 26 28 30 34 43 48 50 54 57 65 75 76 80 84 85 86 92 98 103 107 110 114 116 119 121 122 123 127 128 129 130 131 137 138 140 142 143 144 146 151 156 159 161 162 163 166 168 170 173 174 175 179 180 181 182 其余为 A 类 关键词关键词 DNA 序列分类 聚类分析延拓法 Lance 和 Williams 距离 马尔可夫 法 一 问题重述一 问题重述 1 11 1 题目背景题目背景 1 2000 年 6 月 人类基因组计划中 DNA 全序列草图完成 预计 2001 年可 以完成精确的全序列图 此后人类将拥有一本记录着自身生老病死及遗传进化的 全部信息的 天书 2 这本 天书 是由 4 个字符 A T C G 按一定顺序排成的无间隔的长 约 30 亿的序列 除了这 4 个字符表示 4 种碱基以外 人们对它包含的 内容 知之甚少 因此 破译这部世界上最巨量信息的 天书 是二十一世纪最重要的 任务之一 3 为解读这部 天书 首先要研究 DNA 全序列具有什么结构 以及由这 4 个字符排成的看似随机的序列中隐藏着什么规律 这也是生物信息学最重要的 课题 1 21 2 题目信息题目信息 1 DNA 序列分为编码区与非编码区 编码区是用于编码蛋白质的序列片段 即由这 4 个字符组成的 64 种不同的 3 字符串 其中大多数用于编码构成蛋白质 的 20 种氨基酸 2 在不用于编码蛋白质的序列片段中 A 和 T 的含量特别多些 于是以某 些碱基特别丰富作为特征去研究 DNA 序列的结构也取得了一些结果 3 利用统计的方法还发现序列的某些片段之间具有相关性 这些发现说明 DNA 序列中存在着局部的和全局性的结构 充分发掘序列的结 构对理解 DNA 全序列有十分重要的意义 目前在这项研究中最普通的思想是省略 序列的某些细节 突出特征 然后将其表示成适当的数学对象 1 31 3 题目要求题目要求 1 有 20 个已知类别的人工制造的 DNA 序列 见附件 1 其中序列标号 1 10 为 A 类 11 20 为 B 类 从中提取特征 构造分类方法 并用这些已知类 别的序列 衡量所选分类方法是否足够好 2 用 1 中的分类方法对另外 20 个未标明类别的人工序列 见附件 1 标号 21 40 进行分类 根据分类效果对方法不断完善 将得到的最终结果用 序号 按从小到大的顺序 标明它们的类别 A 类或 B 类 无法分类的不写入 要求详细描述所选的分类方法 给出计算程序 若论文中部分地使用了现成 的分类方法 应将方法名称准确注明 3 已知 182 个自然 DNA 序列 见附件 2 它们都较长 同样用以上所选 的分类方法对它们进行分类 并根据分类效果对方法不断完善 像 2 中一样 给出最终的分类结果 二 二 名词解释名词解释 1 编码区与非编码区 编码区是指 DNA 上编码蛋白质的序列片段 而非编码区不 用于编码蛋白质 2 聚类分析 由已知数据 计算各个观察个体或变量之间亲疏关系的统计量 再 根据某种准则 最短距离法 最长距离法 中间距离法 重心法等 使同一 类内的差别较小 而类与类之间的差别较大 最终将观察个体或变量分为若干 类的分类方法 其中 对样品所作的分类为Q 型聚类 对变量所作的分类为R 型聚类 3 相似性度量 对数值型数据而言 两个个体的相似度是指它们在欧氏空间中互 相邻近的程度 而对分类型数据而言 两个个体的相似度与它们取值相同的属 性的个数有关 4 样品 每个观察个体即每条DNA序列为一个样品 5 样品变量 每个样品所具有的不同特征用不同的变量来表示 变量数等于特征 数 6 碱基丰度 每条 DNA 序列中碱基 A G C 或 T 出现的频率 三三 问题分析问题分析 DNA 序列分类问题要求在对 DNA 序列的一些规律和结构有所了解的基础上 从 20 个已知类别的人工制造的 DNA 序列中提取特征 构造分类方法 并用所选 择的分类方法对其余未知类别的 20 个人工制造的 DNA 序列以及 182 个自然 DNA 序列进行分类 3 13 1 建模目标的分析建模目标的分析 DNA 序列分类是一个复杂的统计分析问题 数据量大 影响因素多 无法直 接从 20 条已知类别的人工制造的 DNA 序列中提取出所有的有效特征 因此有必 要对这 20 条 DNA 序列进行预处理 观察并分析数据预处理结果 归纳总结出 A 类和 B 类的有效特征 将其表示 成适当的数学对象 并选择适当的分类方法 建立普遍意义下数学模型 再用得 到的模型对其余未知类别的 20 个人工制造的 DNA 序列以及 182 个自然 DNA 序列 进行分类 由题意 建立的数学模型应该保证分类结果具有以下特点 1 类别间差异尽量大 2 类别内差异尽量小 3 样品能够尽可能的落入 A B 范围 且只能落入其中的一个 3 23 2 建模及求解方向建模及求解方向 1 分析已知类别的 DNA 序列 1 20 的结构 提取出相应的特征 主要的特征有 碱基的丰度 碱基或碱基序列的重复出现情况 碱基或碱基 序列之间的相邻情况 不同碱基的丰度之比 如碱基 A 与碱基 T 的丰度之比 等 2 根据提取出的特征 选用合适的分类方法 对数据进行预处理后 尝试以下方法建立模型 1 根据聚类分析法 建立模型一 由题意 DNA 序列分类属于对样品所做的分类 为 Q 型聚类 首先引入样品 变量 例如可选择碱基 T 的丰度 碱基 G 的丰度 碱基 T 与碱基 G 的丰度之比 碱基 A 与碱基 T 的丰度之比等 由已知数据 计算出每条已知类别的人工制造的 DNA 序列的各个样品变量值 存入向量中 根据相似性度量原理 计算 20 个样品两两之间的 Lance 和 Williams 距离 选择相距最远的两个样品 假设为样品 3 和样品 16 分别作为 A 类和 B 类 再 分别以样品 3 和样品 16 为标准点 通过分别计算样品 3 和样品 16 与其余 18 个 样品之间的 Lance 和 Williams 距离 找出与其相距最近的一个样品 假设为样 品 1 和样品 18 归为一类 此时 新的标准点变为样品 1 与样品 3 的中点 样 品 16 与样品 18 的中点 然后再以新的标准点为基准 分别找出与其相距最近的 一个样品归为一类 逐步进行下去 直至 20 个样品被明显分成 A B 两类 2 根据马尔可夫法 建立模型二 以单个碱基为单位 分别统计编号 1 10 和编号 11 20 人工制造的已知类别 的 DNA 序列中 4 种碱基出现的次数 再以相邻的两个碱基为单位 共 16 种组合 情况 分别统计编号 1 10 和编号 11 20 的 DNA 序列中 16 种碱基对出现的次数 为满足大样本需求 将 A 类和 B 类中的 10 条 DNA 序列组合起来看作两个大样品 单个碱基或碱基对出现 不包括上一条链的末尾碱基与下一条链的初始碱基组合 的情况 的次数为 10 条序列之和 由条件概率的思想 分别求出 A 类和 B 类大样品中在已知当前碱基种类 可 以为 A G C T 中任何一个 的条件下 下一个碱基分别为 A G C T 的概率 存入两个矩阵 1 和 2 中 对于任何一条给定的 DNA 序列 可将其看作一个已经发 生的事件 说明该事件发生的概率比较大 用矩阵 1 和矩阵 2 分别求出这一事件 发生 即形成当前 DNA 序列 的概率 若用矩阵 1 算出该编号的 DNA 序列出现的 概率较大 则该编号的 DNA 序列属于 A 类 否则属于 B 类 3 模型的初步检验与改进 用编号 1 20 已知类别的序列 分别衡量模型一与模型二中所选方法是否足 够好 不断改进 尽可能使 1 20 号 DNA 序列在所选分类方法下 所得结果与已 知分类完全一致 改进时 对于聚类分析法 可以尝试改变样品变量的个数或者 改变样品变量的组合方式 对于马尔可夫法 可以尝试引进中间变量 运用隐马 尔可夫法求解 4 模型的进一步检验与完善 1 用以上的得到的两种分类方法对编号 20 40 未知类别的人工序列 182 个自然序列进行分类 2 通过计算样品方差 均值等比较两种分类方法得到的分类结果与建模目 标 类别间差异尽量大 类别内差异尽量小 样品能够尽可能的落入 A B 范 围 且只能落入其中的一个 的接近程度 3 选择更接近建模目标的一种分类方法作为最终的分类方法 其分类结果 即为最终结果 四 基本假设四 基本假设 1 假设所给的 DNA 序列片段中没有断句和标点符号 2 假设具有特殊碱基的 DNA 序列中 特殊碱基可以剔除 其影响可以忽略 3 较长的 182 个自然序列与已知类别的 20 个样本序列具有共同的特征 4 假设给定的DNA序列均是从全序列中随机截取出来的 无法确定序列的起始位 无法从序列中辨认出氨基酸 所以 在对DNA 序列分类时 从碱基层次上进行 分类 而不是从氨基酸层次上分类 五 定义与符号说明五 定义与符号说明 i n 各个样品中碱基i出现的数量 i 为 A T C 或 G i N 第 i 个样品的总碱基数目 i F 各个样品中碱基i的丰度 i 为 A T C 或 G i x 各个样品的第 i 个特征变量 ij f 各个样品中碱基 i 和碱基 j 的比值 i j 为 A T C 或 G i Y 第 i 个样品的特征向量 ij d 向量 i Y和向量 j Y间的 Lance 和 Williams 距离 p 特征向量的分量个数 即向量的维数 k 特征向量的第 k 个分量 n 样品的个数 ik x 特征向量 i 的第 k 个分量 max d 不同向量代表的p维空间中任意两点间 Lance 和 Williams 距离的最大值 min d 不同向量代表的p维空间中任意两点间 Lance 和 Williams 距离的最小值 i y 聚类分析中 i 类的标准向量 i 为 A 或 B 六 六 数据预处理数据预处理 1 A 类和 B 类样品单个碱基丰度的计算 用 maTlab 编写程序 见附件 3 分别求出 20 条已知类别的人工制造的 DNA 序列中 4 种碱基的丰度 绘出散点图如下 图 6 1 1 单个碱基丰度比较图 分析上图可得 A 类和 B 类 DNA 序列中碱基 T 和碱基 G 的丰度有明显差异 而碱基 A 和碱基 C 的丰度则比较接近 2 A 类和 B 类样品不同碱基丰度之比的计算 用 matlab 编写程序 见附件 4 分别求出 20 条已知类别的人工制造的 DNA 序列中 不同碱基的丰度之比 包括 TA f CA f GA f CT f GT f GC f 绘出散点 图如下 图 6 1 2 不同碱基丰度之比的比较图 分析上图可得 A 类和 B 类 DNA 序列中 碱基 T 与碱基 A 的丰度之比 碱 基 G 与碱基 A 的丰度之比 碱基 C 与碱基 T 的丰度之比 碱基 G 与碱基 T 的丰度 之比有明显差异 而碱基 C 与碱基 A 的丰度之比 碱基 G 与碱基 C 的丰度之比则 比较接近 3 将编号 1 40 人工制造的 DNA 序列的中 碱基 T 的丰度 碱基 G 的丰度 碱基 T 与碱基 A 的丰度之比 碱基 G 与碱基 A 的丰度之比 碱基 C 与碱基 T 的丰度之 比 碱基 G 与碱基 T 的丰度之比 用表格的形式加以表达 见附件 5 表 1 4 统计所有 DNA 序列中碱基 A T C G 的比例 发现在未知类别的人工制造的 DNA 序列以及自然序列中并非只存在 A T C G 四种碱基 还存在 n s w y 等特殊碱基 这可能和生物自身需要完成的特定功能有关 具体列表如下 表 2 特殊的 DNA 序列及特殊碱基种类 DNA 序列 特殊碱基 DNA 序列 特殊碱基 人工 37 号 s 自然 131 n 自然 71 n 自然 147 n 自然 101 n s 自然 169 n 自然 105 r s w y 由上表可知 编号 1 20 的人工制造的 DNA 序列中并未出现特殊碱基 所以 在提取特征时不需要考虑特殊碱基的影响 同样 在处理编号 21 40 的人工制造 的 DNA 序列以及 182 条自然序列时 也不必考虑特殊碱基的影响 使用数据时 可将特殊碱基直接剔除 七 模型的建立与求解七 模型的建立与求解 7 1 模型一 聚类分析模型一 聚类分析延拓延拓模型模型 要使 DNA 序列的分类能够尽量科学合理 集中要解决的问题是让分类后的 样品满足 同类样品间的差异性尽可能小 不同类样品间的差异性尽可能大 为达到上述目的 引入聚类分析模型对不同的 DNA 序列进行分类 7 1 1 模型模型一一的的建立建立 聚类分析方法根据分类对象的不同可以分为两类 1 对样品所作的分类 即 Q 型聚类 2 对变量所作的分类 即 R 型聚类 此问题将给出的不同 DNA 序列 看成是不同的样品 选用 Q 型聚类进行具体求解 1 样品特征变量的引入 为了刻画不同样品的性质 需要对样品引入统一的特征作为样品特征变量 特征变量的确定来源于聚类分析前对数据进行预处理得到的分析结果 1 样品中 A C T G 的碱基丰度 样品 i 中 A 碱基丰度的计算 AAi FnN 1 其他碱基丰度的计算方法同上 绘出编号 1 20 的人工制造的已知类别的 DNA 序列中 4 种碱基丰度的离散统 计图 图 6 1 1 观察该散点图 进行数据分析可得 DNA 序列中碱基 A 和碱基 C 在分类 A 和 B 中的区分不大 均大致在相同的频率区间内波动 故不选用碱基 A 和碱基 C 的 丰度作为特征区分 而 DNA 序列中碱基 T 和碱基 G 在分类 A 和 B 中的区分较大 A 类和 B 类相应的碱基丰度分别集中在不同的频率区间范围内 故选用碱基 T 和 碱基 G 的丰度作为特征区分 将 T 的碱基丰度作为样品的第 1 个特征变量 记为 1 x 将 G 的碱基丰度作为样品的第 2 个特征变量 记为 2 x 2 样品不同碱基间的比例 样品 i 中碱基 T 和碱基 A 的比值计算 GTGT fnn 2 其他碱基比例的计算方法同上 绘出编号 1 20 的人工制造的已知类别的 DNA 序列中不同碱基的丰度之比的 离散统计图 图 6 1 2 观察该散点图 进行数据分析可得 DNA 序列中碱基 T 和碱基 A 的丰度之比 以及碱基 G 和碱基 T 的丰度之比在分类 A 和 B 中的区分较大 A 类和 B 类相应 的碱基丰度之比分别集中在不同的频率区间范围内 故选用碱基 T 和碱基 A 的丰 度之比以及碱基 G 和碱基 T 的丰度之比作为特征区分 将碱基 T 和碱基 A 的比值作为样品的第 3 个特征变量 记为 3 x 将碱基 G 和碱基 T 的比值作为样品的第 4 个特征变量 记为 4 x 2 样品特征数据的向量转化 把上述得到的 4 种特征变量分别作为一个向量的四个分量 用该向量作为样 品特征向量来描述不同样品 由附件 5 表 1 编号 1 40 样品的 1 x 2 x 3 x 和 4 x 的值分别为表中的第 1 2 3 6 列 于是得到编号 1 20 的样品的 20 个特征向量如下 1 0 1351 0 3964 0 4545 2 9333 Y 2 0 1532 0 4144 0 5667 2 7059 Y 3 0 0631 0 4505 0 2333 7 1429 Y 4 0 2883 0 1802 0 6809 0 6250 Y 5 0 1081 0 4234 0 4615 3 9167 Y 6 0 1261 0 3964 0 3590 3 1429 Y 7 0 1892 0 3604 0 5385 1 9048 Y 8 0 1892 0 3694 0 6774 1 9524 Y 9 0 1532 0 4324 0 7391 2 8235 Y 10 0 1364 0 4091 0 7500 3 0000 Y 11 0 5000 0 1000 1 4103 0 2000 Y 12 0 5000 0 1455 1 5278 0 2909 Y 13 0 5182 0 1273 2 0357 0 2456 Y 14 0 5000 0 1182 1 6667 0 2364 Y 15 0 6455 0 0636 2 2188 0 0986 Y 16 0 4636 0 0909 1 2750 0 1961 Y 17 0 2636 0 1364 0 7436 0 5172 Y 18 0 5000 0 0909 1 7188 0 1818 Y 19 0 5636 0 0727 2 5833 0 1290 Y 20 0 5636 0 0636 2 8182 0 1129 Y 3 不同样品的相似性度量 分析编号 1 20 的样品 因为 20 个已知类别的 DNA 序列的样品变量均属于数值型数据 所以两个个 体的相似度是指它们在欧氏空间中互相邻近的程度 据此 引用距离测度来描述 不同样品的相似性 距离测度小的两个样品 相似性较高 反之 距离测度大的 两个样品 相似性较低 为了排除不同变量之间的相互影响 以及减弱较大数据出现时对结果的不良 影响 即减弱较大值 包括异常值 的敏感度 选用 Lance 和 Williams 距离来描 述距离测度 进而衡量不同样品间的相似性 此外 Lance 和 Williams 距离还与 样品变量的单位无关 使结果无量纲化 向量 i Y和向量 j Y间的 Lance 和 Williams 距离为 1 p ikjk ij k ikjk xx d xx 3 用公式 3 计算所有向量所代表的p维空间中所有样品点之间的两两距离 由排列组合知识 所有向量 n 个 进行两两组合的个数为 2 n C 分别计算 出每个组合的 Lance 和 Williams 距离 本次聚类中选用的向量个数为 n 20 一共有 2 20 190C 种组合 用 matlab 编 程 见附件 6 求解出所有组合的 Lance 和 Williams 距离 并对数据进行比较 得出 max 3 3671d 4 根据距离测度进行分类 1 样品数据分成两类 由上述得到的 max d 查找 max d所对应的向量组合 假定该向量组合是向量 i Y 和向量 j Y 则将第 i 个样品和第 j 个样品分为 A B 两类 可以令 i 样品为 A 类 令 j 样品为 B 类 分别将 i Y和 j Y作为 A B 两类的标准向量 A y B y 对剩余样品进 行分类 2 剩余样品分类 样品 i 和样品 j 分完类后 还剩余 n 2 个样品未进行分类 将这 n 2 个 样品数据分别和 A 类的标准向量 A y 进行组合 计算出每个组合的 Lance 和 Williams 距离 将所得的距离进行比较 得出最小的 min d 查找 min d所对应的向 量 假定该向量是 a Y 则将该向量和样品 i 分为一类 同属于 A 类 用同样的方 法把这 n 2 个样品数据分别和 B 类的标准向量 B y 进行组合 得出最小的 min d 假定该组合所对应的向量是 b Y 则将该向量和样品 j 分为一类 同属于 B 类 此时得到 A 组为 a Y i Y B 组为 b Y j Y A B 两类标准的重新计算 将此时 A B 组中的所有向量分别求出平均值得 到 A B 类的新的标准向量 A 类的标准向量 2 Aai yYY 4 B 类的标准向量 2 Bbj yYY 5 3 上述步骤后还剩余 n 4 个样品未进行分类 依照 2 剩余样品分类给 出的方法不断重复进行计算 对所有的剩余样品均实现分类 7 1 2 模型模型一一的的求解求解 按照上述方法首先计算得到这些样品中向量 3 Y 和向量 20 Y间的 Lance 和 Williams 距离最大 则将第 3 个样品和第 20 个样品分为 A B 两类 令第 3 个样 品为 A 类 第 20 个样品为 B 类 按照 7 1 1 中的步骤依次进行分类 用 matlab 编程 见附件 7 求解得到分类结果如下 A 类的样品编号为 1 2 3 5 6 7 8 9 10 17 B 类的样品编号为 4 11 12 13 14 15 16 18 19 20 7 1 3 模型模型一一的的检验检验与与改进改进 1 模型一的改进与可行性分析 由以上分类结果可知 用聚类分析延拓法对编号 1 20 人工制造的 DNA 序列 进行分类的结果与已知分类结果并非完全一致 在此分类方法下 第 4 条 DNA 序列不再属于 A 类 而属于 B 类 第 17 条 DNA 序列不再属于 B 类 而属于 A 类 因此 有必要对模型进行改进 可以改变样品变量的组合方式 选择碱基 T 的丰度 碱基 T 与碱基 A 的丰度 之比 碱基 C 与碱基 T 的丰度之比 碱基 G 与碱基 T 的丰度之比作为四个样品变 量 分别设为 1 x 2 x 3 x 和 4 x 由附件 5 表 1 编号 1 40 样品的 1 x 2 x 3 x 和 4 x 的值分别为表中的第 1 3 5 6 列 得到编号 1 20 的样品的 20 个特征向量如下 1 0 1351 0 4545 1 2667 2 9333 Y 2 0 1532 0 5667 1 0588 2 7059 Y 3 0 0631 0 2333 3 4286 7 1429 Y 4 0 2883 0 6809 0 3750 0 6250 Y 17 0 2636 0 7436 0 9310 0 5172 Y 18 0 5000 1 7188 0 2364 0 1818 Y 19 0 5636 2 5833 0 2581 0 1290 Y 20 0 5636 2 8182 0 3065 0 1129 Y 用公式 3 计 算 20 个向量所代表的 4 维空间中所有样品点两两之间的 Lance 和 Williams 距离 并按照 7 1 1 中的距离测度法对编号 1 20 人工制造的 DNA 序 列进行分类得到的分类结果如下 A 类的样品编号为 1 2 3 4 5 6 7 8 9 10 B 类的样品编号为 11 12 13 14 15 16 17 18 19 20 由以上分类结果可知 改变样品变量的组合方式 选择碱基 T 的丰度 碱基 T 与碱基 A 的丰度之比 碱基 C 与碱基 T 的丰度之比 碱基 G 与碱基 T 的丰度之 比作为四个样品变量后 用聚类分析延拓法对编号 1 20 人工制造的 DNA 序列进 行分类的结果与已知分类结果完全一致 所以 该分类方法可行 2 模型一的进一步检验与实践 1 用模型一中改进后的聚类分析延拓法 对编号 21 40 人工制造的 DNA 序列进行分类 对附件 7 中的程序稍作修改 求解得到分类结果如下 A 类的样品编号为 22 23 25 27 29 30 34 35 36 37 39 B 类的样品编号为 21 24 26 28 31 32 33 38 40 2 用模型一中改进后的聚类分析延拓法 对 182 个自然 DNA 序列进行分 类 同样对附件 7 中的程序稍作修改 求解得到分类结果如下 B 类的样品编号为 7 10 12 22 23 24 26 28 30 34 43 48 50 54 57 65 75 76 80 84 85 86 92 98 103 107 110 114 116 119 121 122 123 127 128 129 130 131 137 138 140 142 143 144 146 151 156 159 161 162 163 166 168 170 173 174 175 179 180 181 182 其余的自然 DNA 序列为 A 类 7 1 4 模型模型一一改进后改进后分类效果的评价分类效果的评价 1 求出 A 类中 10 条 DNA 序列 4 个样品变量 碱基 T 的丰度 碱基 T 与碱基 A 的丰度之比 碱基 C 与碱基 T 的丰度之比 碱基 G 与碱基 T 的丰度之比 的平 均值 作为 A 类的标准点 a 求 出 B 类中 10 条 DNA 序列 4 个样品变量的平均值 作为 B 类的标准点 b a 0 1393 0 5311 1 5172 3 2803 b 0 5020 1 8300 0 2618 0 2131 2 计算 A 类中 10 个样品点与标准点 a 之间的 Lance 和 Williams 距离 并求出 距离的平均值和标准差 平均值 1 0 7224a 标准差 1 0 5298 a 计算 B 类中 10 个样品点与标准点 a 之间的 Lance 和 Williams 距离 并求出距离的平均值和 标准差 平均值 2 2 6923a 标准差 2 0 4646 a 3 计 算 A 类中 10 个样品点与标准点 b 之间的 Lance 和 Williams 距离 并求出 距离的平均值和标准差 平均值 1 2 5116b 标准差 1 0 5193 b 计算 B 类中 10 个样品点与标准点 b 之间的 Lance 和 Williams 距离 并求出距离的平均值和 标准差 平均值 2 0 7227b 标准差 2 0 5323 b 4 对以上数据进行分析 若分类方法合理 那么不同类别之间的差别应尽可能大 即 1 a 与 2 a的差别 1 b与 2 b 的差别应尽可能大 同类之间的差别应尽可能小 即 1a 2a 1b 和 2b 应尽可能小 此外 定义相关系数 12 12 XX X r XX X 为 a时 表示选择标准点a进行评价时的相关系数 12 12 0 5048 aa a r aa X 为b时 表示选择标准点b进行评价时的相关系数 12 12 0 5879 bb b r bb 由均值和标准差的含义 为使 A 类与 B 类之间的差别尽可能大 那么相关系 数 r 应该尽可能小 由以上结果 a r和 b r的大小均为 0 5 左右 可知该分类方法合 理 且能够达到较好的分类效果 7 2 模型二模型二 马尔可夫马尔可夫模型模型 7 2 1 模型的建立模型的建立与与求解求解 1 DNA 序列的马尔可夫链转换 把 DNA 的一个样品序列看成是一个系统 组成该 DNA 序列的不同位置的碱 基看成是这个系统中的相应的不同状态 DNA 的长度为 N 则该系统有 N 个状 态 分别记为 1 2 N S SS 每个状态对应一个碱基 这样给定的一条长度为 N 的 DNA 序列转化成有 N 个状态组成的系统 即为 1 2 N S SS 随着时间的推移 系 统从某一状态转移到另一状态 设 t q 为时间 t 的状态 tj qS 系统在 t 时间 的状态 1t q 只与其在时间 t 1 的状态相关 1ti qS 其概率为 1 tt P qq 这样 将该系统转换成一个离散的一阶马尔可夫链 2 不同碱基的组合情况 将 4 个碱基进行两两组合 用表格的形式进行考虑 表 3 两个碱基的组合情况 碱基 A 碱基 T 碱基 C 碱基 G 碱基 A AA AT AC AG 碱基 T TA TT TC TG 碱基 C CA CT CC CG 碱基 G GA GT GC GG 两个碱基组合排列一共有 16 种情况 P AT表示在前一状态为碱基 A 的 情况下后一状态出现碱基 T 的概率 其他字母表示意义和上述同 3 中间状态发生的概率 中间状态 系统中除第一个状态和最后一个状态外均称为中间状态 以 A 类情况下 AT 的情况为例进行计算 给定的 A 类样品为编号 1 10 的 DNA 序列 将这 10 条 DNA 序列组合作为一个 大样品 A n表示在该样品中出现碱基 A 的个数 AT n表示一条 DNA 链中碱基 A 后出现 碱基 T 的情况组合起来看成一个新的碱基 AT 该样品中碱基 AT 的个数 P AT的概率 P AT ATA P T Ann 6 根据上述的计算公式 算出 A 类其余 15 种组合的概率 B 类的计算情况和 A 类的情况相同 根据 A 类 B 类给出的 DNA 序列具体进行计算 得到以下表格 表 4 A 类不同中间状态发生的概率表 碱基 A 碱基 T 碱基 C 碱基 G 碱基 A P AA 0 3662 P AT 0 1911 P AC 0 1911 P AG 0 2516 碱基 T P TA 0 2456 P TT 0 3275 P TC 0 1696 P TG 0 2573 碱基 C P CA 0 2359 P CT 0 1385 P CC 0 0974 P CG 0 5282 碱基 G P GA 0 2601 P GT 0 0644 P GC 0 2029 P GG 0 4726 表 5 B 类不同中间状态发生的概率表 碱基 A 碱基 T 碱基 C 碱基 G 碱基 A P AA 0 3707 P AT 0 4081 P AC 0 1184 P AG 0 1028 碱基 T P TA 0 2701 P TT 0 5931 P TC 0 0858 P TG 0 0511 碱基 C P CA 0 2182 P CT 0 4727 P CC 0 1636 P CG 0 1455 碱基 G P GA 0 3063 P GT 0 3964 P GC 0 0811 P GG 0 2162 4 待判定分类 DNA 序列的概率计算 给定一条长度为 N 的 DNA 将其转换为系统状态序列 1 2 N S SS 每一个系 统状态对应同一位置 DNA 序列给出的一个碱基 计算该 DNA 序列产生的概率 该 DNA 序列系统产生的概率计算公式 123121321 P nnn S S S S P S P S S P S S P S S 7 第一个状态的出现概率均设为 1 即 1 P S 1 分别根据 A 类 B 类给出的中间状态出现的概率 得到该 DNA 序列产生概 率 5 DNA 分类的判定 将上面得到的两个 DNA 序列产生概率经行比较 如果通过 A 类中间状态的 概率计算值远远大于 B 类中间状态的概率计算值 则将该状态归为 A 类 同样 若通过 B 类中间状态的概率计算值远远大于 A 类中间状态的概率计算值 则将该 状态归为 B 类 6 实际数据的代入计算 对已知类别的 20 个样品依照上述方法进行分类 表 6 编号为 1 10 个样本产生概率统计表 样品 1 样品 2 样品 3 样品 4 样品 5 样品 6 样品 7 样品 8 样品 9 样品 10 A类 数 据计算 0 10 e 60 0 14 e 59 0 14 e 56 0 21 e 65 0 14 e 59 0 35 e 56 0 63 e 59 0 31 e 61 0 74 e 61 0 79 e 62 B类数据 计算 0 14 e 76 0 13 e 75 0 86 e 80 0 89 e 65 0 60 e 80 0 12 e 73 0 22 e 71 0 24 e 74 0 74 e 82 0 52 e 82 所属类别 A A A B A A A A A A 表 7 编号为 11 20 个样本产生概率统计表 样品 11 样品 12 样品 13 样品 14 样品 15 样品 16 样品 17 样品 18 样品 19 样品 20 A 类 数据 计算 0 25 e 63 0 44 e 64 0 28 e 65 0 64 e 67 0 28 e 60 0 27 e 65 0 51 e 71 0 28 e 65 0 20 e 68 0 71 e 70 B 类 数据 计算 0 20 e 49 0 89 e 51 0 97 e 54 0 97 e 55 0 26 e 41 0 61 e 54 0 29 e 68 0 11 e 52 0 19 e 53 0 12 e 54 所属 类别 B B B B B B B B B B 根据上面的判断结果 只有 4 号样品的类别出现了偏差 说明马尔可夫模型 进行判定具有一定的合理性 可以进一步推广 对其他的 DNA 序列进行判别 7 2 2 分类结果统计分类结果统计 用 matlab 编写程序 见附件 8 对编号 21 40 人工 DNA 序列以及 182 个自 然序列进行分类 结果如下 1 编号 21 40 人工 DNA 序列分类结果 A类 22 23 25 26 27 29 30 32 33 34 35 36 37 39 B类 21 24 28 31 38 40 2 182 个自然序列分类结果 A类 B类 1 2 3 4 5 6 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 24 25 26 27 28 29 31 32 33 34 35 36 37 38 39 40 41 42 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 66 67 68 69 70 71 72 73 74 77 78 79 80 81 82 83 84 86 87 88 89 90 91 93 94 95 96 97 98 99 100 101 102 104 105 106 108 109 111 112 113 115 116 117 118 120 121 123 124 125 126 127 129 130 132 133 134 135 136 137 139 140 141 142 143 145 146 147 148 149 150 152 153 154 155 157 158 160 164 165 167 168 169 171 172 173 174 175 176 177 178 179 180 181 总计 182 个自然序列分入 A 类的个数 154 个 182 个自然序列分入 B 类 的个数 28 个 7 23 30 43 65 75 76 85 92 103 107 110 114 119 122 128 131 138 144 151 156 159 161 162 163 166 170 182 7 2 3 模型的评价模型的评价 对模型的分类效果进行评价 评价标准 同类样品间的差异性较小 不同类 样品间的差异性较大 则这样的分类效果较好 反之同类样品间的差异性较大 而不同类样品间的差异性较小 这样的分类效果就不够理想 1 检验样本的数据处理 选取前 20 号 DNA 样品作为检验样本 依照上述马尔可夫模型的分类结果 A 类为 1 2 3 5 6 7 8 9 10 B 类为 4 11 12 13 14 15 16 17 18 19 20 将每一个样品的最大产生概率 每个样品有两个生成概率 分别用 A 类和 B 类的相应数据进行求解 其中较大的一个即为该样品的最大产生概率 作为该样 品的特征变量 由于最大产生概率的数量级很小 为了后期数据处理的方便和准确性 将各 样品最大产生概率通过分别取对数的方法进行数据处理 A 类样品对应的计算结果为 61 0000 59 8539 56 8539 59 8539 56 4559 59 2007 61 5086 61 1308 62 1024 B 类样品对应的计算结果为 65 0506 49 6990 51 0506 54 0132 55 0132 41 5850 54 2147 68 5376 52 9586 53 7212 54 9208 2 同类样品间的差异性判定 对上面的数据 按照类为单位分别计算相应的平均值 得到 A 类的平均值为 60 3638 A B 类的平均值为 53 5714 B 用公式 2 2 1 1 1 n i i S n 分别求 A 类和 B 类的无偏方差 2 S为方差 i 为 第 i 个最大概率取对数的值 为平均概率取对数的值 为标准差 这样进一步得到 A 类和 B 类的标准差分别为2 6464 A 6 6158 B 该 标准差的值较小说明同类样品间的差异性较小 3 不同类样品间的差异性判定 将A类的9个样品分别用A类的中间状态概率和B类的中间状态概率进行计 算得到两组平均值 1 60 3638 A 2 66 2520 A 1A 表示 A 类的 9 个样 品用 A 类的中间状态概率得到平均值 2A 表示 A 类的 9 个样品用 B 类的中间状 态概率得到平均值 同理将 B 类的 11 个样品分别用 A 类的中间状态概率和 B 类的中间状态概率 进行计算得到两组平均值 1 76 2691 B 2 53 5714 B 1B 表示 B 类的 11 个样品用 A 类的中间状态概率得到平均值 2B 表示 B 类的 11 个样品用 B 类的中 间状态概率得到平均值 A 类样品的两个标准差 2 2 11 1 1 1 n AiA i S n 1 2 6464A 2 2 22 1 1 1 n AiA i S n 2 3 1973A 1A S表示 A 类样品通过 A 类的中间状态概率计算得到的方差 1A 表示 A 类样品 通过 A 类的中间状态概率计算得到的标准差 2A S表示 A 类样品通过 B 类的中间 状态概率计算得到的方差 2A 表示 A 类样品通过 B 类的中间状态概率计算得到 的标准差 A 类样品中分别用 A 类的中间状态概率和 B 类的中间状态概率进行计算 所 得结果的差异性 1A2 A 1A2 0 3273 A A r A r表示 A 类样品中两种中间状态概率 计算结果的相关系数 把 B 类的 11 个样品用上述的方法做同样的处理 得到 B 类样品中两种中间 状态概率计算结果的相关系数 B1B2 B B1B2 0 8809 r 上述相关系数间接刻画了两类别间的差异性大小 值越小 则说明不同类样 品间的差异性越大 反之 值越大 则说明不同类样品间的差异性越小 选取两个系数中的较差值 B r 作为最后的相关性系数 依照最差值选取原则 若最差值也能取得较好的结果 则得到整个结果的合理性 该值较小 从而得 到不同类样品间的差异性较大 依照上述给定的判断标准 得到马尔可夫模型的分类效果较好的结论 7 37 3 模型一模型一与与模型二选择模型二选择的层次分析模型的层次分析模型 本文给出了两种不同的模型对 DNA 序列进行了分类 现用层次分析法综合 对两者进行比较 选出其中更优的一个模型 该层次分析中 目标层为选取两个模型中更优的一个 准则层为模型好坏的 评价标准 确定为两点 1 模型分类情况和给定前 20 个样品的分类情况的相似 程度 2 模型自身的分类效果 同类样品间的差异性较小 不同类样品间的差异 性较大 则自身的分类效果较好 方案层为聚类分析延拓模型和马尔科夫模型 根据上述的层次 给出层次分析的结构图 图 7 3 1 层次分析法的结构图 1 构造准则层对于目标层的成对比较矩阵 1 1 1 a A a 其中 a 是表示模型与给定分类相似程度与模型自身分类效果重要性之比 求得权向量为 1 11 T a W aa 参数确定 a 取 3 表示模型与给定分类相似程度比自身分类效果稍微重要 2 构造方案层对准则层各因素的成对比较矩阵并求出权向量相对于指标层模型 与给定分类相似程度 1 1 1 1 1 1 b B b 其中 b1 表示聚类分析延拓模型与马尔科夫模型对于指标层模型与给定分 类相似程度 求得权向量为 1 1 11 1 11 T b W bb 参数确定 由于聚类分析延拓模型与模型给定人工序列的分类结果完全吻 合 而马尔科夫模型有一条链出现偏差 取 b1 5 求得相对于指标层模型自身分类效果矩阵 2 2 2 1 1 1 b B b 其中 b2 表示聚类分析延拓模型与马尔科夫模型对于指标层模型自身分类 效果 求得权向量为 2 2 22 1 11 T b W bb 参数确定 由于聚类分析延拓模型和马尔科夫模型相比 模型自身分类效果 较佳 取 b2 3 由此求得方案最终的权重矩阵为 1 12 1212 1111 11111111 T bbaa W abababab 代入数据解得 1 0 8125 0 1875 T W 由结果可得聚类分析延拓模型明显优于马尔科夫模型 所以采用聚类分析延 拓模型作为最终的 DNA 序列分类模型 八八 模型的评价与改进 模型的评价与改进方向方向 8 8 1 1 模型的优点模型的优点 1 对整个问题的解决 基于聚类分析法和马尔可夫法分别建立了模型一和模 型二 从中选择较优的模型 避免单一模型解决问题 引起较大误差 2 建立模型一时 对传统的聚类分析法进行改进 建立了聚类分析延拓模型 3 进行模型二的求解时 将结果存入矩阵 简化了运算 增加了效率 4 在进行模型一和模型二的选择时 采用层次分析法 综合考虑评判分类效 果的多个因素 8 8 2 2 模型的缺点模型的缺点 1 在 用 马尔可夫模型进行分类时 编号 1 20 人工 DNA 序列的分类结果没有达 到与已知分类结