《DNA序列研究数学建模问题.doc》由会员分享，可在线阅读，更多相关《DNA序列研究数学建模问题.doc（23页珍藏版）》请在三一文库上搜索。
1、题目 序列摘要本文主要研究序列的结构问题，通过建立相应的数学模型，对序列中所隐藏的规律进行研究和分析，给出了解决问题的最优方案，并且对模型进行了评价和推广。对于问题一，为了挖掘序列的特征将其分为类和类，以20种基本氨基酸为目标，利用软件编程得出每一行每一种氨基酸出现的概率；再运用主成分分析法进行降维，利用软件进行数据处理得到矩阵；然后再将模糊聚类问题转化为如下优化问题： 用模糊聚类分析方法来获取样本与聚类中心的加权距离最小的最佳分类，使其分类结果与问题所给的答案一致，如果如下表：类1235678910类411121314151617181920由此可见，误差较小可忽略不计，证明此种方法可行。再
2、以此种方法将21到40行序列分为两类，其中类，14681011131820； 类，235791214151617对于问题二：在问题一的基础上，对182个氨基酸序列进行分类，采取与问题一相同的方法进行分类，分类结果见问题二的求解。 总的来说，本模型在未知数据特征的情况下很好的将数据进行分类，成功地解决了此次数学建模的序列问题，是聚类分析问题的一个有效而且具有较强实用性的方法。关键词：主成分分析 模糊聚类分析 软件 软件一、问题重述1.1背景分析 随着测序时代的到来，越来越多生物的全基因组序列正逐渐展现于人们的眼前。如何从中挖掘有用的信息成为对当今生物学乃至整个科学领域的一个挑战。本文主要致力于对
3、序列结构以及序列中所隐藏规律的研究。1.2问题重述2000年6月，人类基因组计划中DNA全序列草图完成，预计2001年可以完成精确的全序列图，此后人类将拥有一本记录着自身生老病死及遗传进化的全部信息的“天书”。这本大自然写成的“天书”是由4个字符A，T，C，G按一定顺序排成的长约30亿的序列，其中没有“断句”也没有标点符号，除了这4个字符表示4种碱基以外，人们对它包含的“内容”知之甚少，难以读懂。破译这部世界上最巨量信息的“天书”是二十一世纪最重要的任务之一。在这个目标中，研究DNA全序列具有什么结构，由这4个字符排成的看似随机的序列中隐藏着什么规律，又是解读这部天书的基础，是生物信息学（Bi
4、oinformatics）最重要的课题之一。虽然人类对这部“天书”知之甚少，但也发现了DNA序列中的一些规律性和结构。例如，在全序列中有一些是用于编码蛋白质的序列片段，即由这4个字符组成的64种不同的3字符串，其中大多数用于编码构成蛋白质的20种氨基酸。又例如，在不用于编码蛋白质的序列片段中，A和T的含量特别多些，于是以某些碱基特别丰富作为特征去研究DNA序列的结构也取得了一些结果。此外，利用统计的方法还发现序列的某些片段之间具有相关性，等等。这些发现让人们相信，DNA序列中存在着局部的和全局性的结构，充分发掘序列的结构对理解DNA全序列是十分有意义的。目前在这项研究中最普通的思想是省略序列的
5、某些细节，突出特征，然后将其表示成适当的数学对象。这种被称为粗粒化和模型化的方法往往有助于研究规律性和结构。作为研究DNA序列的结构的尝试，提出以下对序列集合进行分类的问题：问题一：下面有20个已知类别的人工制造的序列（见附录），其中序列标号110 为A类，11-20为B类。请从中提取特征，构造分类方法，并用这些已知类别的序列，衡量你的方法是否足够好。然后用你认为满意的方法，对另外20个未标明类别的人工序列（标号2140）进行分类，把结果用序号（按从小到大的顺序）标明它们的类别（无法分类的不写入）：A类 ；B类 。 请详细描述你的方法，给出计算程序。如果你部分地使用了现成的分类方法，也要将方法
6、名称准确注明。问题二：在同样网址的数据文件Nat-model-data 中给出了182个自然DNA序列，它们都较长。用你的分类方法对它们进行分类，像1）一样地给出分类结果。二、 模型假设结合本题实际，为确保模型求解的准确性和合理性，我们排除了一些因素的干扰，提出以下几点假设：1、序列具有连续有序性；2、较长的182个自然序列与已知类别的20个样本序列具有共同的特征；3、64种3字符串压缩为20组后不影响分类的结果；4、各序列中碱基三联组（即3字符串）的起始位置和基因表达不影响分类的结果。5、分类的结果在一定误差范围内认为合理；三、符号说明3.1为了便于问题的求解，我们给出以下符号说明：表示第种
7、氨基酸表示20个变量构成的20维随机变量原始数据的平均值随机变量的均值协方差矩阵矩阵的特征值主成分的贡献率主成分与原始指标的相关系数氨基酸综合评分函数样本隶属于第类的隶属度模糊划分矩阵样本与第类聚类中心之间的距离分类标准原始数据标准差3.2名词解释:1、氨基酸：含有氨基和羧基的一类有机化合物的通称。生物功能大分子蛋白质的基本组成单位，是构成动物营养所需蛋白质的基本物质。2、模糊聚类分析法：聚类分析是对事物按一定要求进行分类的数学方法。实际的分类问题常伴有模糊性,因此,聚类问题用模糊数学的方法解决更确切。3、隶属度：若对论域（研究的范围）中的任一元素，都有一个数与之对应，则称为上的模糊集，称为对
8、的隶属度。四、 问题分析4.1对问题一的分析问题一中需要20行序列的特征并对其进行分类，使其分为两类，并使前十行属于类后十行属于类中，因为序列较长，且为字符，不能直接对其进行分类，这里我们将其转化为三字符的氨基酸形式，以查询资料得到的20种基本氨基酸为目标，利用软件编程求出每一行这20种氨基酸分别出现的概率，以此作为实验数据。然而20种氨基酸种类仍然较多，需要对其进行主成分分析将其降维以方便下面的运算。数据处理完之后，建立相应的20*6的矩阵，利用模糊聚类分析方法从事先给出的个划分出发，通过不断地反复修改样本的类别、聚类中心以及样本隶属于各类别的隶属度，来获得样本与聚类中心的加权距离最小的最佳
9、分类。再将分类结果与已分好的类别作比较，判断误差是否较小结果是否正确。在以相同的模糊聚类分析方法将后20行分为、两类，分别填入题中横线。4.2对问题二的分析因为前面的模糊聚类分析方法经第一次分类的检验，证明其误差较小，结果比较准确，可以继续用这种方法对接下来的182个序列进行分类。因为序列较长且多，所以可以用编程直接导入数据，算出每行20种氨基酸出现的概率，然而数据太过庞大，根据第一题可以看出序列中各种氨基酸出现的概率之间方差较大，可以先对其进行主成分分析，然后重新建立目标矩阵直接编程对其进行分类。五、模型的建立与求解5.1问题一的模型建立与求解5.11根据对问题一第一小问的分析可知，附件里的
10、序列较长较多，都为字符，而且总体上只有四个字符数，不能直接对这些字符进行分类，并且这些字符所代表的在生物上会三三结合形成氨基酸以创造出蛋白质，所以其有一定的连续性，必须将这些单个字符转化成三个字符相结合而成的氨基酸形式来进行研究以找出其相应特征进行分类。根据资料可知，这四个字符组成了64种不同的3字符串，这64种3字符串构成生物蛋白质的20种氨基酸；通过对这20种氨基酸出现概率的特征研究结合软件编程（程序见附录），得出前20行序列的20*20概率矩阵，见附录。然而数据还是太大，氨基酸种类太多，不利于特征的研究；所以根据这概率矩阵利用降维的思想将多个变量转化为少数几个综合变量对其进行主成分分析。
11、主成分分析模型建立如下：假设原来研究对象是20种氨基酸，分别用来表示，这20个变量构成的20维随机向量为 。设随机向量的均值为，为反应标准化的数据之间相关关系密切程度的统计指标建立协方差矩阵为，值越大，说明有必要对数据进行主成分分析，其中，为原始变量，的相关系数。为对称矩阵（即），只需计算其上三角元素或下三角元素即可，其计算公式为：然后对进行线性变化，考虑原始变量的线性组合：再根据协方差矩阵求出特征值、主成分贡献率和累计方差贡献率，确定主成分个数。特征值是各主成分的方差，它的大小反映了各个主成分的影响力。解特征方程求出。因为是正定矩阵，所以其特征值都为正数，将其按大小顺序排列，即。特征值是各主
12、成分的方差，它的大小反映了各个主成分的影响力。主成分的贡献率，累计贡献率为。根据选取主成分个数的原则，特征值要求大于1且累计贡献率达80%-95%的特征值所对应的1,2，其中整数即为主成分的个数。接下来建立初始因子载荷矩阵，解释主成分。因子载荷矩阵是主成分与原始指标的相关系数，揭示了主成分与氨基酸比率之间的相关程度。过程详细数据见附件软件运行结果图。最后根据特征值及累计贡献率的数值，分析得到3个主成分，根据筛选出来的3个主成分重新设立成分矩阵（详情见附录图1）。利用模糊聚类分析将前二十行序列按照各个样本点之间、样本点与样本点子集合之间的关系体系确定氨基酸之间的“亲疏关系”从而对其进行正确而合理
13、的分类。模糊聚类分析模型建立过程如下：设待分类的个样本为其中每个样本均有个特征值，即，根据题意可知需要将样本划分为2类。为了获得一个最佳的模糊分类，需要一个分类准则，定义目标函数为样本与聚类中心的加权距离，加权系数为样本的隶属度函数的次方。于是，模糊聚类问题可以转化为如下优化问题如下：其中为样本隶属于第类的隶属度，模糊划分矩阵；为两个聚类中心集合；为加权指数，当时，模糊聚类就退化为硬2均值聚类；通常的最佳选择范围为，一般是比较理想的取值。表示样本与第类聚类中心之间的距离，定义为式中为的正定矩阵，当时，即为欧氏距离。样本隶属于类的程度，即隶属度，可在0到1之间取值。这样，样本不在明确的属于某一类
14、，而是对于每一个类别都有一个隶属度，隶属度的数值越大表明样本隶属于该类别的程度越大，反之则说明隶属于该样本的程度越小。模糊聚类的这种模糊划分描述了样本聚类过程中的模糊现象，从而可以获得更为合理的聚类结果。为了使目标函数到达极小点,其应该满足：利用上述公式，运用迭代方法通过软件（编程见附录）可以获得聚类分析结果，如下：1235678910411121314151617181920根据模糊聚类分析得出的结果，1到10行除了第4行均属于一类，11到20行以及第4行属于另一类。与题中给出的序列标号110 为A类，11-20为B类的分类情况相比，只有第4行出现误差。与数据庞大的题目相比，误差相对较小，可
15、以忽略不计；所以用模糊聚类分析结合主成分分析法得出的分类结果较为准确，可以推广到更多此类问题的研究上。5.12由第一题的第一问可知此种分类方法较为准确，所以根据序列中4个字符组成的64种不同的3字符串，其中大多数用于编码构成蛋白质的20种氨基酸这种生物特性，按照前一小问序列分类方法将其转化为3字符串的氨基酸来进行分析；用软件编程对其进行统计，统计出每一行各种氨基酸出现的概率；再运用软件对其进行主成分分析，分析得到四个主成分，经过主成分分析后得到的成分矩阵（见附录图2）。以此达到简化计算过程的目的后，再结合软件利用模糊聚类分析法得出分类，如下：14681011131820235791214151
16、61719下图为聚类结果图，通过此图可以清晰地看出20行氨基酸被分成各具特征的两类，两类各自的数据特别紧凑，证明用此种方法分类准确，误差较小。5.2问题二的模型建立与求解由第一题可知上面使用的分类方法较为准确，所以按照前二十行序列分类方法对其先进行统计，统计出每一行各种氨基酸出现的概率，再对其进行主成分分析，得到8个主成分；再结合软件利用模糊聚类分析法得出分类。因为数据太多太杂乱，可以先将其标准化，使之便于分析和比较。令其中为原始数据；为原始数据的平均值，为原始数据的标准差。再根据实际情况按一个准则（如绝对值减数法）给需要进行研究的元素两两之间都赋以区间内的一个相似系数，得到一个模糊矩阵；再依
17、次取的截矩阵，以此种等价关系矩阵对元素进行分类，最后确定最佳分类。运用此种方法结合软件进行编程。分类结果如下：类：236791012172223242526283034363738404346474850515354565760636566727475767780838485868788929394969798991021031071101131141161191211221231251261271281291301311331381401421431441461471491501511521531561591601611621631641651661671681691701731741751
18、77178179180181182类：1458111314151618192021272931323335394142444549525558596162646768697071737879818289909195100101104105106108109111112115117118120124132134135136137139141145148154155157158171172176根据下面的聚类结果图，可以清晰直观地看出182行氨基酸根据自身数据间的特征被分成互相独立的两类，两组内的数据都比较紧凑，证明用此种方法分类可行，较为准确。六、模型的评价与改进6.1模型的优点1、本题研究的对
19、象是尚未分类的群体，具有模糊性，所以用模糊聚类分析法比较适宜。2、主成分分析以降维的思想将多个变量中综合为少数几个变量，既能够代表原始变量的绝大多数信息，又互不相关，并且在新的变量基础上，可以更方便进一步统计分析。3、聚类分析是将特征空间中的样本集合按照各个样本点之间、样本点与样本点子集合之间以及样本点子集合之间的相似性测度进行聚类，得出的分类组能充分把握住数据间的特征；并且分类组内的数据具有凝聚性，能较好的反应分类组的特征。4、模糊聚类分析中样本的隶属度可在0到1之间取值，而不是如经典聚类算法（如型聚类法和型聚类）中隶属度只有0和1两种选择，样本聚类的这种模糊划分能够是聚类更加合理。5、软件
20、提供了一个更为直观的图形工具来帮助选择主成分。6、软件可以高效的求解复杂的问题，并可以对系统进行动态的仿真，得到聚类结果图，用强大的图形功能对数值计算结果进行显示。6.2模型的缺点1、未充分考虑分子序列的实际特性，主观认为它的序列具有连续有序性。2、主观的以构成蛋白质的20种氨基酸这种生物特性作为基准处理数据，未充分考虑其他生物特性，具有一定的局限性。2、只考虑了样本序列中3字符串出现的频率作为特征，不一定与实际情况相符。3、主成分分析是从空间生成的角度看问题，具有一定的局限性。4、循环运算效率低，封装性不好。5、在图形工具上不大全面，不容易掌握而且统计结果使用不方便。七、模型的推广在实际问题
21、中，我们经常会遇到研究多个变量的问题，而且在多数情况下，多个变量之间常常存在一定的相关性，比较复杂；这时就需要进行主成分分析。主成分分析是采取一种数学降维的方法，找出几个综合变量来代替原来众多的变量，使这些综合变量能尽可能地代表原来变量的信息量，而且彼此之间互不相关。主成分分析可以与回归分析结合，进行主成分回归分析，以及利用主成分分析进行挑选变量，选择变量子集合的研究。另外，还可以用预报经济统计研究中，用于对经济效益的综合评价研究，对不同地区经济发展水平的评价研究，不同地区经济发展竞争力的评价研究，人民生活水平、生活质量的评价研究等。模糊聚类分析是将事物根据一定的特征矩阵，并按某种特定要求或规
22、律分类的方法，其对象是尚未分类的群体，即存在模糊性。根据这种模糊聚类分析方法，在实际中更能很好的反映出数据的特征，不同种类别中数据更加具有凝聚性、更能反映各个分类的特征。对研究DNA序列的规律性和结构提供了一种有效的分类模型。对人类基因组的研究有现实意义，有利于加快科研步伐。因为现实中大部分的样本具有模糊性，所以此种方法在各种领域被广泛应用。八、参考文献1数学模型，姜启源 主编，高等数学出版，19932谭永基，数学模型，上海复旦大学出版社，19973The Invention of the Genetic Code，Brain Hayes(美)，American ScientistComput
23、ing Science，Jan.-Feb., 19984叶海军，模糊聚类分析技术及其应用研究D，合肥，合肥工业大学，20065张秀梅，王涛，模糊聚类分析方法在学生成绩评价中的应用J，渤海大学学报（自然科学版），2007,28（2）：169-1726主成分分析论文，http： 2011.7.26九、附录Art-model-data1.aggcacggaaaaacgggaataacggaggaggacttggcacggcattacacggaggacgaggtaaaggaggcttgtctacggccggaagtgaagggggatatgaccgcttgg2.cggaggacaaacgggatgg
24、cggtattggaggtggcggactgttcggggaattattcggtttaaacgggacaaggaaggcggctggaacaaccggacggtggcagcaaagga3.gggacggatacggattctggccacggacggaaaggaggacacggcggacatacacggcggcaacggacggaacggaggaaggagggcggcaatcggtacggaggcggcgga4.atggataacggaaacaaaccagacaaacttcggtagaaatacagaagcttagatgcatatgttttttaaataaaatttgtattattatggtat
25、cataaaaaaaggttgcga5.cggctggcggacaacggactggcggattccaaaaacggaggaggcggacggaggctacaccaccgtttcggcggaaaggcggagggctggcaggaggctcattacggggag6.atggaaaattttcggaaaggcggcaggcaggaggcaaaggcggaaaggaaggaaacggcggatatttcggaagtggatattaggagggcggaataaaggaacggcggcaca7.atgggattattgaatggcggaggaagatccggaataaaatatggcggaaagaac
26、ttgttttcggaaatggaaaaaggactaggaatcggcggcaggaaggatatggaggcg8.atggccgatcggcttaggctggaaggaacaaataggcggaattaaggaaggcgttctcgcttttcgacaaggaggcggaccataggaggcggattaggaacggttatgagg9.atggcggaaaaaggaaatgtttggcatcggcgggctccggcaactggaggttcggccatggaggcgaaaatcgtgggcggcggcagcgctggccggagtttgaggagcgcg10.tggccgcggagggg
27、cccgtcgggcgcggatttctacaagggcttcctgttaaggaggtggcatccaggcgtcgcacgctcggcgcggcaggaggcacgcgggaaaaaacg11.gttagatttaacgttttttatggaatttatggaattataaatttaaaaatttatattttttaggtaagtaatccaacgtttttattactttttaaaattaaatatttatt12.gtttaattactttatcatttaatttaggttttaattttaaatttaatttaggtaagatgaatttggttttttttaaggtagttattta
28、attatcgttaaggaaagttaaa13.gtattacaggcagaccttatttaggttattattattatttggattttttttttttttttttttaagttaaccgaattattttctttaaagacgttacttaatgtcaatgc14.gttagtcttttttagattaaattattagattatgcagtttttttacataagaaaatttttttttcggagttcatattctaatctgtctttattaaatcttagagatatta15.gtattatatttttttatttttattattttagaatataatttgaggtatg
29、tgtttaaaaaaaatttttttttttttttttttttttttttttttaaaatttataaatttaa16.gttatttttaaatttaattttaattttaaaatacaaaatttttactttctaaaattggtctctggatcgataatgtaaacttattgaatctatagaattacattattgat17.gtatgtctatttcacggaagaatgcaccactatatgatttgaaattatctatggctaaaaaccctcagtaaaatcaatccctaaacccttaaaaaacggcggcctatccc18.gttaattatt
30、tattccttacgggcaattaattatttattacggttttatttacaattttttttttttgtcctatagagaaattacttacaaaacgttattttacatactt19.gttacattatttattattatccgttatcgataattttttacctcttttttcgctgagtttttattcttactttttttcttctttatataggatctcatttaatatcttaa20.gtatttaactctctttactttttttttcactctctacattttcatcttctaaaactgtttgatttaaacttttgtttctttaag
31、gattttttttacttatcctctgttat21.tttagctcagtccagctagctagtttacaatttcgacaccagtttcgcaccatcttaaatttcgatccgtaccgtaatttagcttagatttggatttaaaggatttagattga22.tttagtacagtagctcagtccaagaacgatgtttaccgtaacgtqacgtaccgtacgctaccgttaccggattccggaaagccgattaaggaccgatcgaaaggg 23.cgggcggatttaggccgacggggacccgggattcgggacccgagga
32、aattcccggattaaggtttagcttcccgggatttagggcccggatggctgggaccc24.tttagctagctactttagctatttttagtagctagccagcctttaaggctagctttagctagcattgttctttattgggacccaagttcgacttttacgatttagttttgaccgt25.gaccaaaggtgggctttagggacccgatgctttagtcgcagctggaccagttccccagggtattaggcaaaagctgacgggcaattgcaatttaggcttaggcca26.gatttactttagcat
33、ttttagctgacgttagcaagcattagctttagccaatttcgcatttgccagtttcgcagctcagttttaacgcgggatctttagcttcaagctttttac 27.ggattcggatttacccggggattggcggaacgggacctttaggtcgggacccattaggagtaaatgccaaaggacgctggtttagccagtccgttaaggcttag28.tccttagatttcagttactatatttgacttacagtctttgagatttcccttacgattttgacttaaaatttagacgttagggcttatcag
34、ttatggattaatttagcttattttcga29.ggccaattccggtaggaaggtgatggcccgggggttcccgggaggatttaggctgacgggccggccatttcggtttagggagggccgggacgcgttagggc30.cgctaagcagctcaagctcagtcagtcacgtttgccaagtcagtaatttgccaaagttaaccgttagctgacgctgaacgctaaacagtattagctgatgactcgta31.ttaaggacttaggctttagcagttactttagtttagttccaagctacgtttacggg
35、accagatgctagctagcaatttattatccgtattaggcttaccgtaggtttagcgt32.gctaccgggcagtctttaacgtagctaccgtttagtttgggcccagccttgcggtgtttcggattaaattcgttgtcagtcgctctrtgggtttagtcattcccaaaagg33.cagttagctgaatcgtttagccatttgacgtaaacatgattttacgtacgtaaattttagccctgacgtttagctaggaatttatgctgacgtagcgatcgactttagcac34.cggttagggcaaag
36、gttggatttcgacccagggggaaagcccgggacccgaacccagggctttagcgtaggctgacgctaggcttaggttggaacccggaaa35.gcggaagggcgtaggtttgggatgcttagccgtaggctagctttcgacacgatcgattcgcaccacaggataaaagttaagggaccggtaagtcgcggtagcc36.ctagctacgaacgctttaggcgcccccgggagtagtcgttaccgttagtatagcagtcgcagtcgcaattcgcaaaagtccccagctttagccccagagtcga
37、cg37.gggatgctgacgctggttagctttaggcttagcgtagctttagggccccagtctgcaggaaatgcccaaaggaggcccaccgggtagatgccasagtgcaccgt38.aacttttagggcatttccagttttacgggttattttcccagttaaactttgcaccattttacgtgttacgatttacgtataatttgaccttattttggacactttagtttgggttac39.ttagggccaagtcccgaggcaaggaattctgatccaagtccaatcacgtacagtccaagtcaccgttt
38、gcagctaccgtttaccgtacgttgcaagtcaaatccat40.ccattagggtttatttacctgtttattttttcccgagaccttaggtttaccgtactttttaacggtttacctttgaaatttttggactagcttaccctggatttaacggccagttt氨基酸概率矩阵编程clc clear allclose allfid1=fopen(data.txt,r); i=1;while (feof(fid1) data=fgetl(fid1); b1=length(findstr(data,aaa)+length(findstr(data,
39、ata); b2=length(findstr(data,aca)+length(findstr(data,aga); b3=length(findstr(data,cac)+length(findstr(data,ctc); b4=length(findstr(data,ccc)+length(findstr(data,cgc); b5=length(findstr(data,gag)+length(findstr(data,gtg); b6=length(findstr(data,gcg)+length(findstr(data,ggg); b7=length(findstr(data,tat)+length(findstr(data,ttt); b8=length(findstr(data,tct)+length(findstr(data,tgt); b9=length(findst