基因组数学建模DOC.docx

1、基因组数学建模DOC封一答卷编号（参赛学校填写）：答卷编号（竞赛组委会填写）：论文题目： （标明A、B、C、D之一） B 组 别：(填写研究生、本科生、专科生或中学生) 本科生 参赛队员信息(必填)： 姓 名学 号联系电话参赛队员1参赛队员2参赛队员3 参赛学校： 沈阳理工大学 封二答卷编号（参赛学校填写）： 答卷编号（竞赛组委会填写）：评阅情况（学校评阅专家填写）：学校评阅1.学校评阅2.学校评阅3. 评阅情况（联赛评阅专家填写）：联赛评阅1.联赛评阅2.联赛评阅3.目录摘要： 1一、问题重述 2二、基本假设 4三、符号说明 4四、模型建立与求解 4五、模型的检验 8六、模型的评价 9参考文

2、献 10附录 10摘要： 在世界生物界的飞速发展的今天，我国生物界对于基因组组装也有了一定的突破。尤其是在测序技术上，我国从第一代技术发展到第二代技术，现在正步入第三代技术，但是目前能直接读取的碱基对序列长度远小于基因组序列长度，为了让组装效果更完整，连续和准确，我们利用概率论与数理统计，c语言，和多目标规划，以及c+、matlab等软件，参考了国内外生物界、医学界有关基因组组装技术，以olc技术为核心，建立了集完整、准确、连续为一体的优化模型，来解决测序中可能会出现的问题。具体如下：对于问题一： 首先利用olc技术，将附件中给出的测序之后的碱基对，建立优化模型，该模型能够利用多条约束条件来求

3、质量的最小值。利用c+软件来进行两两比较，选择出重复的部分作为结点，其他的部分作为有向线段，利用有向图来选择一条质量最重的一条链，但是可能会存在测序中个别碱基对的错误，所以，将所测出的的最重的这条链所有的结点全部去掉，再将剩下的碱基对重新利用olc法来组装，重复该种实验3-4次，见附录的程序，一定会有质量相同的两条链，则这两条链的组成就是最完整准确的。这个程序的试验次数少，完成组装的时间少。 针对组装后的基因组中可能会出现重复片段的问题，我们利用c+软件，建立了优化模型，采取的方法是在olc技术的基础上，检测有向图，判断结点的重复次数，进行标记，直到检测出同一结点前后两条有向线段并不相同，才是

4、无重复基因组的基因组装。 对于问题二：现有一个全长约为120,000个碱基对的细菌人工染色体（BAC）， 采用Hiseq2000测序仪进行测序，测序深度（sequencing depth）约为70，即基因组每个位置平均被测到约70次。利用c+软件，建立优化模型，在olc技术的基础上，基于问题一所做的程序（见附录），将附件一和附件二的碱基对输入到程序中，组装成完整、准确、连续的基因组，所得的结果就是最完整、准确、连续的基因组。关键字：olc技术 c+ 碱基对 优化模型 matlab 多目标规划一、问题重述 快速和准确地获取生物体的遗传信息对于生命科学研究具有重要的意义。对每个生物体来说，基因组包

5、含了整个生物体的遗传信息，这些信息通常由组成基因组的DNA或RNA分子中碱基对的排列顺序所决定。获得目标生物基因组的序列信息，进而比较全面地揭示基因组的复杂性和多样性，成为生命科学领域的重要研究内容。 测序技术始于20世纪70年代，伴随着人类基因组计划的实施而突飞猛进。从第一代到现在普遍应用的第二代，以及近年来正在兴起的第三代，测序技术正向着高通量、低成本的方向发展。尽管如此，目前能直接读取的碱基对序列长度远小于基因组序列长度，因此需要利用一定的方法将测序得到的短片段序列组装成更长的序列。通常的做法是，将基因组复制若干份，无规律地分断成短片段后进行测序，然后寻找测得的不同短片段序列之间的重合部

6、分，并利用这些信息进行组装。当然，由于技术的限制和实际情况的复杂性，最终组装得到的序列与真实基因组序列之间仍可能存在差异，甚至只能得到若干条无法进一步连接起来的序列。对组装效果的评价主要依据组装序列的连续性、完整性和准确性。连续性要求组装得到的（多条）序列长度尽可能长；完整性要求组装序列的总长度占基因组序列长度的比例尽可能大；准确性要求组装序列与真实序列尽可能符合。利用现有的测序技术，可按一定的测序策略获得长度约为50100个碱基对的序列，称为读长（reads）。基因组复制份数约为50100。基因组组装软件可根据得到的所有读长组装成基因组，这些软件的核心是某个组装算法。常用的组装算法主要基于O

7、LC（Overlap/Layout/Consensus）方法、贪婪图方法、de Bruijn图方法等。一个好的算法应具备组装效果好、时间短、内存小等特点。新一代测序技术在高通量、低成本的同时也带来了错误率略有增加、读长较短等缺点，现有算法的性能还有较大的改善空间。问题一：试建立数学模型，设计算法并编制程序，将读长序列组装成基因组。你的算法和程序应能较好地解决测序中可能出现的个别碱基对识别错误、基因组中存在重复片段等复杂情况。问题二：现有一个全长约为120,000个碱基对的细菌人工染色体（BAC）， 采用Hiseq2000测序仪进行测序，测序策略以及数据格式的简要说明见附录一和附录二，测得的读长

8、数据见附录三，测序深度（sequencing depth）约为70，即基因组每个位置平均被测到约70次。试利用你的算法和程序进行组装，并使之具有良好的组装效果。二、基本假设假设read1均为500bp片段的前88个，read2为后88个假设只有个别碱基对错误、基因组中存在重复片段的情况三、符号说明Ann表示DNA的第一条链；Bnn表示与之相对应的第二条链；Cnn表示结点；AAnn表示第一条链质量；BBnn 表示第二条链质量；CCnn 表示节点链质量；Pathnn表示路径节点；Distancenn表示链的质量；Read1表示测序后该片段的第一条链的前88个碱基对；Read2表示测序后该片段的第一

9、条链的后88个碱基对；四、模型建立与求解问题一：优化模型：4.1解决个别碱基对识别错误的问题（一）将所给的读长组成基因组方法一：OLC（Overlap/Layout/Consensus）方法1含义： 当前流行的很多序列拼接算法，包括PHRAP7,TIGR8等，都广泛的采用OLC（Overlap/Layout/Consensus）方法。先对所有shotgun小片段进行两两比对，获取重叠部分，并且对这些重叠进行衡量，再将重叠的片段相互连接，得到DNA目标序列的一个大致轮廓；最后构建所求的DNA目标序列。其实就是以每个shotgun片段为顶点，以重叠关系为边，构造一个有向图，并在图中找出找出一条经过

10、每个顶点一次且仅一次的最佳路径。2 OLC（Overlap/Layout/Consensus）方法对重叠序列的处理采用OLC（Overlap/Layout/Consensus）方式的拼接算法一直伴随着鸟枪测序 ，如PHRAP7,TIGR8等，这也是最传统和最经典的方式，利用计算机图论的知识，将shotgun集合中的所有片段fragments看做一个个结点如果两个片段之间重叠，那么就在这两个结点之间画上一边，构成一个有向图G，然后在图G中找每一个结点一次且仅一次的一条路径，就将DNA序列拼接问题转化成了一个图论中的hamilton路径问题。3 主要步骤 以附录一中的read1，read2，rea

11、d3，read4为例a)将read集合中所有片段fragments两两比对，获得存在的重叠部分信息;b)以片段为顶点，重叠的片段相互连接，形成有向图如图；在这一步中，phrap算法根据上一步得到的重叠信息将存在重叠的片段组合起来形成一个称为contig20的结构。c)在有向图找出一条路径，经过每个结点一次且仅一次，顺序抽取出改路径上的fragments，重叠拼接成“consensus”序列。OLC（Overlap/Layout/Consensus）方式构建的read示意图通过在有向图G中查找路径进行DNA序列拼接的时候，TIGR算法在shotgun集合中对repats采取如下方法进行处理1）含

12、repeats的fragments留到最后再进行拼接，以最大限度的利用不含repeats的fragments中的信息。2）碰上含repeats的fragments时，相应的提高序列匹配的标准，如果repeats的测序精度上和正常fragments相比有特别大的差异的话那么就可以排除掉含低测序精度repeats的fragments3）比fragments平均长度的一半还短的repeats通常不是问题因为这种repeats通常会被一个fragments所覆盖，不会影响到拼接4）对repeats的两端进行测序是很有帮助的这样能够帮助那些包含repeats的fragments依据其两端来定位。方法二：

13、 贪婪算法 1.含义： 贪婪算法（又称贪心算法）是指，在对问题求解时，总是做出在当前看来是最好的选择。也就是说，不从整体最优上加以考虑，他所做出的仅是在某种意义上的局部最优解。贪心算法不是对所有问题都能得到整体最优解，但对范围相当广泛的许多问题他能产生整体最优解或者是整体最优解的近似解2 贪婪算法的基本处理： 贪婪算法（Greedy algorithm）是一种对某些求最优解问题的更简单、更迅速的设计技术。用贪婪法设计算法的特点是一步一步地进行，常以当前情况为基础根据某个优化测度作最优选择，而不考虑各种可能的整体情况，它省去了为找最优解要穷尽所有可能而必须耗费的大量时间，它采用自顶向下，以迭代的

14、方法做出相继的贪心选择，每做一次贪心选择就将所求问题简化为一个规模更小的子问题，通过每一步贪心选择，可得到问题的一个最优解，虽然每一步上都要保证能获得局部最优解，但由此产生的全局解有时不一定是最优的，所以贪婪法不要回溯。贪婪算法是一种改进了的分级处理方法。其核心是根据题意选取一种量度标准。然后将这多个输入排成这种量度标准所要求的顺序，按这种顺序一次输入一个量。如果这个输入和当前已构成在这种量度意义下的部分最佳解加在一起不能产生一个可行解，则不把此输入加到这部分解中。这种能够得到某种量度意义下最优解的分级处理方法称为贪婪算法。对于一个给定的问题，往往可能有好几种量度标准。初看起来，这些量度标准似

15、乎都是可取的，但实际上，用其中的大多数量度标准作贪婪处理所得到该量度意义下的最优解并不是问题的最优解，而是次优解。因此，选择能产生问题最优解的最优量度标准是使用贪婪算法的核心。一般情况下，要选出最优量度标准并不是一件容易的事，但对某问题能选择出最优量度标准后，用贪婪算法求解则特别有效。最优解可以通过一系列局部最优的选择即贪婪选择来达到，根据当前状态做出在当前看来是最好的选择，即局部最优解选择，然后再去解做出这个选择后产生的相应的子问题。每做一次贪婪选择就将所求问题简化为一个规模更小的子问题，最终可得到问题的一个整体最优解。3.基本思路：A）建立数学模型来描述问题。B）把求解的问题分成若干个子问

16、题。C）对每一子问题求解，得到子问题的局部最优解。D）把子问题的解局部最优解合成原来解问题的一个解。实现该算法的过程：从问题的某一初始解出发；while 能朝给定总目标前进一步 do求出可行解的一个解元素；由所有解元素组合成问题的一个可行解。方法三：De Bruijn法：在 组合数学ak- aryde Bruijn序列B(k,n)秩序n以荷兰数学家命名 nicolaas Govert de Bruijn是a 循环序列 指定的字母表A以大小k为哪些每可能 subsequence 长度n在A一次确切地出现作为连贯字符序列。这样序列有以下物产：其中每一B(k,n)有长度k有k!/k分明De Brui

17、jn序列B(k,n).de Bruijn序列可以通过采取a修建 汉密尔顿的道路n-尺寸 de Bruijn图表k标志(或等效地， a Eulerian周期 a (n− 1) -尺寸de Bruijn图表)，或者通过 有限领域.例子 每个边缘在这张三维的de Bruijn图表对应于四个数字序列： 标记端点边缘离开的三个数字被标记边缘的那个跟随了。 如果你横断被标记的边缘1从000，你到达在001，从而表明subsequence 0001的出现在de Bruijn序列。 要横断每个边缘确切地一次是确切使用每一个16个四位数序列一次。例如，假设我们走以下Eulerian道路：000， 0

18、00， 001， 011， 111， 111， 110， 101， 011， 110， 100， 001， 010， 101， 010， 100， 000。这对应于以下de Bruijn序列：0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1八个端点接下来出现于序列：0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1

19、0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1 . 0 0 0 0 1

20、1 1 1 0 1 1 0 0 1 0 1 . 0 0 0 0 1 1 1 1 0 1 1 0 0 1 0 1我们然后回到出发点。 每一个八个3数字序列(对应于八个端点)两次确切地出现和四位数字的序列的每一十六(对应于16个边缘)一次确切地出现。序列可以用于缩短对a的强力攻击PIN-象没有Enter键并且不接受的代码锁持续n数字进入。 例如，一把数字式门锁以一个四位数字的代码将有B(10， 4)解答，以长度10000。 所以，只10000则+ 3则= 10003则(作为解答循环)新闻是需要的打开锁。 尝试所有代码将分开地要求4 10000 = 40000新闻。De Bruijn的标志程序化在通

21、报对象附近写(例如a轮子机器人)能使用辨认它角度通过审查n面对定点的连贯标志。 灰色代码 能使用作为相似的转台式位置内码机制。（二）将所测出的的最重的这条链所有的结点全部去掉，再将剩下的碱基对重新利用olc法来组装。（三）重复该种实验3-4次，一定会有质量相同的两条链，则这两条链的组成就是最完整准确的。42解决基因组中存在重复片段的问题（一）OLC（Overlap/Layout/Consensus）方法 以附录一中的read1，read2，read3，read4为例d)将read集合中所有片段fragments两两比对，获得存在的重叠部分信息;e)以片段为顶点，重叠的片段相互连接，形成有向图如

22、图；在这一步中，phrap算法根据上一步得到的重叠信息将存在重叠的片段组合起来形成一个称为contig20的结构。f)在有向图找出一条路径，经过每个结点一次且仅一次，顺序抽取出改路径上的fragments，重叠拼接成“consensus”序列。OLC（Overlap/Layout/Consensus）方式构建的read示意图通过在有向图G中查找路径进行DNA序列拼接的时候，TIGR算法在shotgun集合中对repats采取如下方法进行处理5）含repeats的fragments留到最后再进行拼接，以最大限度的利用不含repeats的fragments中的信息。6）碰上含repeats的fra

23、gments时，相应的提高序列匹配的标准，如果repeats的测序精度上和正常fragments相比有特别大的差异的话那么就可以排除掉含低测序精度repeats的fragments7）比fragments平均长度的一半还短的repeats通常不是问题因为这种repeats通常会被一个fragments所覆盖，不会影响到拼接8）对repeats的两端进行测序是很有帮助的这样能够帮助那些包含repeats的fragments依据其两端来定位。（二）判断结点的重复次数，进行标记，直到检测出同一结点前后两条有向线段并不相同，才是无重复基因组的基因组装。问题二： 基于问题一所做的程序，组装成完整、准确、

24、连续的基因组。五、模型的检验为了检验我们的模型是否具有较强的适用性，我们利用现有算法和程序，在olc技术的基础上，检测有向图并找出测序中可能出现的个别碱基对识别错误、基因组中存在重复片段等复杂情况。首先两两比较，找出重复片段如图一，然后将重复片段作为点，其他的作为线，连成有向图接着找出质量最重的一条线最后扔掉该线所有碱基如图二，剩下的重复实验3-4次，确定正确的序列。随后将一个全长约为120,000个碱基对的细菌人工染色体（BAC）， 采用Hiseq2000测序仪进行测序，测序策略以及数据格式的简要说明见附录一和附录二，测得的读长数据见附录三，测序深度（sequencing depth）约为7

25、0，即基因组每个位置平均被测到约70次。使之具有良好的组装效果。以此来检验算法的正确性经过与原序列进行比较，可以得出我们的模型算法对于普遍的情况具有一定的适用性。 图一 图二随后利用复原算法对其进行复原，以此来检验算法。其步骤如下：我们首先利用PDFXChangeViewer软件将下列所示PDF文档转换为bmp图像，然后利用Photoshop图像处理软件利用基准线和裁剪工具，将其划分为196*2145像素点的10个碎纸片，然后随机编码这些碎纸片，利用碎纸片复原的中文算法对这些随机编码的碎纸片进行复原。六、模型的评价（一）模型的优点本文采用优化模型。线性规划作为一种比较古老的优化和决策模型，它发

26、展的已经比较完善，他可以解决很多实际的问题，对于各种问题都是有用，而且解题的过程简单明了，并且对于不同的问题有多种方法进行求解，是一种比较实用且简单的方法。它有与问题领域无关切快速随机的搜索能力，搜索从群体出发，具有潜在的并行性，使用概率机制进行迭代，具有随机性，具有可扩展性，容易与其他算法结合。 试验次数较少，可以更快速地组成基因组。（二）模型的缺点但是在实际应用中如果出现大量的数据，他就必须借助其他的方法进行结合求解。计算量也比较大。编程实现比较复杂,首先需要对问题进行编码,找到最优解之后还需要对问题进行解码，另外三个算子的实现也有许多参数,如交叉率和变异率,并且这些参数的选择严重影响解的

27、品质,而目前这些参数的选择大部分是依靠经验。附件一和附件二数据庞大，一一输入程序中，会比较复杂。参考文献1、贪心算法在系统故障诊断策略生成中的应用-计算机系统应用-2011年 第1期2、基于分类树和贪心算法的测试数据自动生成方法-计算机工程与设计-2011年 第8期 (32)3、浅谈贪心算法在排课系统中的应用-电脑与电信-2011年 第7期4、许宁；基于OLC算法的软件可靠性预测分析D;西安理工大学;2007年5、张博峰；全基因组DNA测序中的片段拼接方法及其并行处理D；长沙；国防科学技术大学；2003附录#includeusing namespace std;const int W=88;

28、const int L=46847; class QYpublic: QY() int i,j; for(i=0;iL;i+) for(j=0;jw;j+) Aij=0；; OLC(); Show(); input(); private: int N; int n; int ALW; int BLW; int CLW;int AALW; int BBLW; int CCLW;int pathLL;int distanceLL; ;QY:input() int a,b,c,n,i,j,k=0,p,q; int MM=20000000; int allMM; int bllMM; cout输入附录

29、一allMM; cout输入附录二allMM; for(i=0;iMM;i+) if (alli=+) for(j=0;j88;j+) Akj=alli-88; AAkj=alli+1; i+; k+; for(i=0;iMM;i+) if (blli=+) for(j=0;j88;j+) Bkj=blli-88; BBkj=blli+1; i+;k+;QY:Show() coutpath is:endl; for(i=0;iL;i+) for(j=0;jL;j+) coutpathij ; coutendl; QY:OLC() int i,j,k,p,q,m=0,n=0,a=0,b=0;for(j=0;jL;j+) for(i=0;iL;i+) for(q=0;qW;q+) if(Bj88-q=Ai4) if(Bj87-q=Ai3) if(Bj86-q=Ai2) if(Bj85-q=Ai1) if(Bj84-q=Ai0) for(k=0;kq+5;k+) Cmn=Bj84-q; CCmn=BBj84-q; n+; q-;m+;for(a=0;am;a+) for(b=0;bm;b+)pathab=Aab;distanceab=CC