数学建模论文-癌症基因与正常基因的区别.docVIP

癌症基因与正常基因的区别摘要随着人类社会的发展，各种癌症的发病率不断增加。癌症成为影响人类健康的第一大杀手。由科学研究可知其主要源头为基因，为此我们针对癌症样本与正常样本在基因表达上的区别，以及癌症样本具体分类的问题，在合理科学的假设下对题目中给出的40个已知基因表达信息样本的数据用进行处理建立神经网络数学模型，从而有效的实现对未知样本的预测以及对预测癌症样本进行合理的分类，并最终以可视化的方法展现出来。针对问题一：利用附件前40个已知的基因表达信息样本，对其利用神经网络，主成分分析法，遗传算法对数据进行分析与对比并建立相应的数学模型，从而得到其癌症人基因与正常人基因的特征值，利用其与未知样本比较即可得到未知样本的基因信息。 针对问题二:针对癌症基因中可能出现的子类问题，我们主要采用动态模糊聚类分析的方法对问题求解。首先对20个癌症基因样本利用相关系数法求得其相关性。而后抓住主要矛盾以最低程度降低同组样本相似性为代价求得对其分类的最佳阈值，从而将前20个癌症样本最佳分为4个子类。最后用神经网络对已分类的癌症样本进行训练，通过可视化的方法呈现出我们的分类结果。关键词：模糊聚类分析 神经网络 可视化 主成分分析法 1. 问题重述癌症作为当今人类的第一大杀手，其主要特点是破坏性强，难预防。基因作为遗传控制因子对癌症的发病率有着最直接的联系。为此我们需利用题目附件中给出的60个人114个基因表达水平的样本解决以下两个问题。在附件中前20个是癌症病人的基因表达水平的样本(其中还可能有子类), 其后的是20个正常人的基因表达信息样本, 其余的20个是待检测的样本(未知它们是否正常)。问题一：试设法找出描述癌症与正常样本在基因表达水平上的区别, 建立数学模型,及识别方法,去预测待检测样本是癌症还是正常样本.问题二：设计图示 (可视化) 方法,使得在你的数学模型下, 尽量清楚地表现癌症与正常样本在基因表达水平上的区别, 以及癌症样本中是否有子类.2. 问题分析对于问题一：根据附件中前40个已知基因信息样本使用作图工具箱进行分析，将所有样本的基因信息反映在图表上，由此我们发现每个基因都有着其一定的正常表达范围。由于基因种类太多，生物网络过于复杂，为了得到精确的标准，我们将生物技术和计算机技术有效的结合起来对已知信息利用基于遗传算法的神经网络数理化进行分析处理从而建立相对应的数学模型，首先用-1,1对40个已知样本量化使数据简化。由于每个样本中包括都由114个基因，所包含的信息太多。所以我们需要用到主成分分析法对所有已知数据处理找出其特征值致癌因子，而后可通过神经网络随机挑选出的34个训练样本对所建立神经网络进行训练，用其余6个已知样本对网络测试检验。对于问题二：基于问题一所建立的数学模型，将待测的样本分为正常和癌症基因样本。但对于将癌症基因的进一步分类，则需要运用动态模糊矩阵，同时，还要找出阈值的合理求解方法，从而得到最佳的答案。最后，利用可视化的方法将样本的分组具体描绘出来。3.模型的假设与符号说明3.1 模型的假设假设一：基因样本的信息准确无误，真实可信，不含特殊样本。假设二：每个样本间相互独立，每个基因间相互独立。假设三：样本的基因表达方式一致，由统一计算方法所得并可直接参与计算。假设四：样本中只有正常人与癌症病人的基因区别，忽略其他遗传病症。假设五：癌症病人的主要致癌基因均在这114个中。假设六：基因样本的数据中纵向表达一一对应。假设七：癌症患者由致癌因子直接导致，不考虑其他客观因素的影响。假设八：在研究过程中，114个基因不再发生变化影响研究结果。3.2 符号说明:癌症基因样本:第个样本:样本与相似程度:模糊等价矩阵:特征值:相邻阈值之比:第组样本相似度:正常基因样本:第个样本所包含的第个基因信息:模糊相似矩阵：阈值：样本所分组数：样本组数之比：神经网络目标值4.模型的准备本文针对癌症样本与正常样本在基因表达上的不同，利用作图工具箱（程序见附录一）使基因表达信息反映在图表上。观察以下图像可以发现每个正常基因都是有一定的正常表达范围的，超过这个范围的基因可以确定为癌症基因。但是通过作图模型做出的图表只能对样本进行粗略的判断，不能得到精准的结果，所以我们接下来需要利用合理的数理分析对数据进行再次处理建立有效的数学模型来区分癌症与正常的基因。如图一所示：图15.模型的建立与求解5.1 问题一的解决:第一步：针对问题的特殊性，我们主要选择bp神经网络对其求解。根据题中所给60个人的基因表达，对前40个人已知基因进行量化:将前40个正常人与癌症患者具体量化为：=则目标值为：t=-1,-1，-1，-1，-1，-1，-1，-1，-1，-1，-1，-1，-1，-1，-1，-1，-1，-1，-1，-1，1，1，1，1，1，1，1，1，1，1,1，1，1，1，1，1，1，1，1，1第二步：题中所给数据每个人均有114个基因来决定是属于癌症还是正常，变量太多，无疑会增加分析问题的难度与复杂性，并且在许多实际问题中，多个变量之间是具有一定的相关关系的。因此我们选取60个人的主要致癌基因，也即对60行114列数据进行处理，选取主要特征。利用主成分分析法减少输入变量，减少训练时间，提高训练效率，达到提高神经网络泛化能力的目的。得到特征值的贡献率如图2所示（程序见附录）： 图2 特征值的贡献率经计算得到：第一，第二，第三，第四，第五，第六，第七，第八主成分的累积贡献率已达到86.82%，所以我们选取八个主成分，。如图3所示：图3结果得到8行60列为主要成分，利用中的函数，再通过编程得到新的主成分矩阵（程序见附录）。由此可得出影响一个人是正常还是癌症主要由此8个主要成分决定。在选择训练样本和检验样本时，我们取癌症人群和正常人群各17个样本（占样本总数56.7%），这两类人群中的剩余的各3个样本作为检验样本（占样本总数10%），剩余20个作为测试样本（占样本总数33.3%）。首先，对训练样本归一化，设置网络隐含层的神经元数为5个，之所以选5个，是因为当神经元数为3,4,5时，其输出精度都相仿，一般的讲，网络神经元的选择原则是：在能够解决问题的前提下，再加上一个到两个神经元以加快误差的下降速度即可。而当神经元数过大时，会产生其它的问题。然后，运用网络中的函数来创建网络，其中隐含层选择作为传递函数，用作为网络训练函数，网络的全局误差为0.00001，训练速率为0.01，学习次数为5000000次，得出以下结果为（程序见附录）：训练样本为：=-1.0001 -0.9998 -0.9999 -0.9997 -1.0001 -1.0002 -1.0000 -0.9999 -0.9999 -1.0000 -1.0001 -0.9997 -1.0000 -0.9999 -1.0001 -1.0001 1.0006 0.9978 0.9962 1.0075 0.9972 0.9959 1.0136 1.0000 0.9993 0.9939 1.0001 1.0001 0.9999 0.9972 1.0024 1.0000 0.9979误差性能之和：0.00034均方差： 0.00001网络训练误差如图4：检验样本为：-1.0042 -1.0298 -0.7234 0.9820 0.8615 1.0345误差为：0.42% 2.98% 37.66% -2.8% -13.85% 3.45%预测的结果为：-0.2336 0.4908 -0.3489 0.7683 2.0489 -0.0097 0.4940 0.8415 -0.8534 0.4664 -0.7463 -1.2671 0.1883 -0.7680 -0.4632 -1.0935 1.1057 -0.8225 -2.3802 -1.4352图4故在20个测试样本中有8个是正常人，12个癌症患者为41,43，46,49,51,52,54,55,56,58,59,60。5.2 问题二的解决: 对于问题二，我们选用基于动态模糊矩阵的模糊聚类分析方法，建立模型。为了能够对已知癌症基因样本，进行细化，使得可以对癌症基因样本再分类，得到子类。绘制动态聚类图,如图5所示： 17.13.20.8.19.14.10.18.15.16.4.11.1.5.7.6.9.2.3.120.970.960.950.940.930.920.910.890.800.78图5成组的准则：阈值1,0.97,0.96,0.94,0.93,0.92,0.91,0.89,0.80,0.78,组数20,19,16,14,11,9,7,6,4,2,1。随着阈值的减小，组数也将减小，同组中样本数增加，而样本的相似程度降低，即扩大各组样本数量是以降低样本相似性为代价的。在分组的过程之中，我们发现，最优解为两个癌症子类（程序见附录），但是考虑实际情况，有一组中只有一个样本，所以这样的分组并不具有代表性，所以我们选取次优解，即：0.89时，组数为4。即有4个癌症子类。6.模型结果的分析与检验判别致癌因子是一个十分复杂的非线性问题，针对这一问题的特殊性，我们主要利用bp神经网络兼用主成分分析法的方法对问题求解。由于问题中所给数据基因种类太多，我们采用主成分分析法以抓住问题的实质德目的综合指标作为主成分简化数据，在一定程度上没有将可能影响问题结果的全部因素考虑进去。这使问题的求解结果在存在一定范围的误差。但是在问题一中，由所得的结果我们可以清晰明了的看出利用主成分分析法简化同时在求解问题二时，利用微分方程的思想求解，用对数据处理后得到了黑匣子在海底运动轨迹的三维立体模型，由图可知其运动为变加速的曲线运动大致符合黑匣子在海底的实际情况，同时问题二及问题三的结果均在题目所给范围内。验证了所得求解结果的正确性。7.模型的推广与改进方向7.1 模型的推广在如何正确区别正常基因与致癌基因从而判别致癌因子的问题中，我们根据问题的特殊性采用了神经网络算法对问题进行求解。神经网络算法由众多的神经元可调的连接权值连接而成，具有大规模并行处理、分析式信息存储、良好的自组织自学习能力等特点，同时具有很强的非线性映射能力。所以我们认为在解决如何判别致癌因子的问题中所建立的数学模型可推广至一些其他非线性问题上（如公路的运量预测、艾滋病治疗最佳停药时间的确定等）。7.2 模型的改进在解决问题的过程中我们首先利用的是模糊聚类分析,将数据归一化处理后由预测的结果看出其检验样本误差较大.为此我们认为对于正常基因与癌症基因这两类差距较大样本的判别可使用bp神经网络模型,它对于解决非线性问题有着独特的优势.同时对于基因基数大,种类多的问题.我们在解决过程中加入了主成分分析法对数据简化.从而改进了所建立的模型是结果在一定得误差范围更快速的得出。8.模型的优缺点8.1 模型的优点: 1. 对于可能出现影响的次要因素进行假设，抓住主要影响因素，对问题进行分析。2. 利用软件对数据进行处理,使解题效率提高。3. 使用作图使所求解问题的结果简单明了。4. 致癌因子的判别是一个十分复杂的非线性问题，针对这一问题的特殊性我们主要采用神经网络模型对问题求解。5. 在解题过程中兼用主成分分析法求主要特征值简化问题。使所建立的数学模型可普遍适用。8.2 模型的缺点：在问题一中，在采用主成分分析法时没有将其他可能影响的因素考虑进去，使结果存在一定的误差。 在问题二中，动态模糊聚类分析具有一定的误差，相识度归类并不能达到100%。参考文献1姜启源，谢金星，叶俊，数学建模（第三版）m.高等教育出版社，2003.8 2黄洪钟.模糊设计m.北京：机械工业出版社，1999 3鲁开讲，张宁，杨峰，基于零件分类编码的模糊聚类分析j.文理学院20024朱雪龙.应用信息论基础m.北京：清华大学出版社，20025汤服成，一种聚类分析方法j.数学的实践与认识，19886李增，迟道才.基于遗传算法的改进elman神经网络模型的降雨量预测j.沈阳农 业大学学报，2010附录：1. 基因图表分析的matlab程序：a=load(附件2.txt);t=1:114;subplot(3,1,1)for i=1:20 b=a(i,:); plot(t,b,r.-) hold onendgrid xlabel(基因序列)ylabel(基因表达水平)title(癌症病人基因样本)subplot(3,1,2)for i=21:40 b=a(i,:); plot(t,b,g.-) hold onendgrid xlabel(基因序列)ylabel(基因表达水平)title(正常人基因样本)subplot(3,1,3)for i=41:60 b=a(i,:); plot(t,b,b.-) hold onend grid xlabel(基因序列)ylabel(基因表达水平)title(待检测基因样本)2. 求截距阵的matlab程序：clccleara=load(附件2.txt)x=a(1:40,2:115);xi=mean(x);xi=xi;for i=1:114 xx(:,i)=xi;endxxfor i=1:40 for j=1:40 temp=sum(abs(x(i,:)-xx(i,:).*(x(j,:)-xx(j,:); flag=sqrt(sum(x(i,:)-xx(i,:).2)*sqrt(sum(x(j,:)-xx(j,:).2); p(i,j)=temp/flag; endendp=pn=input(参数n=);for i=1:40 for j=1:40 if p(i,j)=n p(i,j)=1; else p(i,j)=0; end endendfuzequmat(p)3. 用相关系数法求模糊相似矩阵的matlab程序：clccleara=load(附件2.txt)x=a(1:20,2:115);xi=mean(x);xi=xi;for i=1:114 xx(:,i)=xi;endxxfor i=1:20 for j=1:20 temp=sum(abs(x(i,:)-xx(i,:).*(x(j,:)-xx(j,:); flag=sqrt(sum(x(i,:)-xx(i,:).2)*sqrt(sum(x(j,:)-xx(j,:).2); p(i,j)=temp/flag; endendp4. 求模糊等价矩阵的matlab程序：function rp=rpower(dati,p)%模糊关系乘法%rp=rpower(r,p)%r为模糊关系,p为r自乘的次数%error(nargchk(2,2,nargin)%if(p=floor(p),error(p must be an integer.);endn m=size(datin);%if(n=m),error(matrix length of each demensoin must agree.);endrp=datin;while p1 for i=1:n for j=1:n rp(i,j)=max(min(rp(i,:);datin(:,j); end end p=p-1;endfunction rk=fuzequmat(r)%fuzzy equivalent matrix 求模糊等价矩阵%rk k=fuzequmat(r)%rk即求得的等价矩阵，k为自乘的次数，r为相似矩阵%error(nargchk(1,1,nargin);k=0;rk=r;while all(all(rk=rpower(rk,2) rk=rpower(rk,2); k=k+1;endrkk=2k%n=input(参数n=);%在r*上对u分类，直接利用floor函数%r*n=floor(r+(1-n)5. bp神经网络模型的matlab程序：clcclose allcleara=load(附件2.txt)%bp建模x=a(1:60,2:115); coef score latent t=princomp(x);z=coef*x;latenta=sum(latent);m=size(latent);temp=0;for i=1:m temp=temp+latent(i); if temp/a=0.85 n=i break; endendz=z(1:n,:);p=z(1:n,1:17,21:37);%原始数据归一化t=-ones(1,17),ones(1,17);p minp maxp=premnmx(p);%设置网络隐单元的神经元数（530验证后5个最好）n=5;%建立相应的bp网络net=newff(minmax(p),n,1,tansig,purelin,traingdm);inputweights=net.iw1,1;inputbias=net.b1;%bp网络参数的设置net.trainparam.show=50;net.trainparam.lr=0.01;net.trainparam.epochs=2000000;net.trainparam.goal=1e-5;%训练网络net tr=train(net,p,t);%对bp网络进行仿真a=sim(net,p)e=a-t;%误差性能函数之和m=sse(e);%求均方差n=mse(e);test1=z(1:n,18:20,38:40);t1=-ones(1,3),ones(1,3);%对新样本进行量化test1=tramnmx(test1,minp,maxp);%仿真a1=sim(net,test1)e1=a1-t1;%误差性能函数之和m1=sse(e1);n1=mse(e1);m,n,m1,n1%预测test=z(1:n,41:60);%对新样本进行量化test=tramnmx(test,minp,maxp);%仿真a=sim(net,test)6. 改进后的bp神经网络模型的matlab程序：clcclose allclear%下载数据load adata.txt%bp建模x=adata(1:60,2:115);coef score latent t=princomp(x);z=coef*x;latenta=sum(latent);m=size(latent);temp=0;for i=1:m temp=temp+latent(i); if temp/a=0.85 n=i break; endendz=z(1:n,:);%p=z(1:n,:);p=z(1:n,2:17,20,24:40);%原始数据归一化t=-ones(1,17),ones(1,17);p minp maxp=premnmx(p);%设置网络隐单元的神经元数n=10;m=1;%遗传算法nvar ncase=size(p);cl=nvar*n+n*m+n+m;popu=100;varb=ones(cl,1)*-1 1;initpop=initializega(popu,varb,gabpeval);%遗传算法g=200;%调用matlab工具箱中的遗传函数x,endpop,bpop,trace=ga(gabpeval,varb,initpop,1e-6 10,maxgenterm,g,normgeomselect,0.09,arithxover,2,nonunifmutation,2 g3);iw,b1,1w,b2,val=gadecode(x);%调用traingdm算法训练bp网络net=newff(minmax(p),n,m,tansig,purelin,traingdm);%设置初始权值和偏置参数net.iw1,1=iw;net.b1=b1;net.lw2,1=1w;net.b2=b2;%bp网络参数的设置net.trainparam.show=50;net.trainparam.lr=0.05;net.trainparam.mc=0.9;net.trainparam.epochs=5000000;net.trainparam.goal=1e-5;%训练网络net tr=train(net,p,t);%对bp网络进行仿真a=sim(net,p)e=a-t;%误差性能函数之和m=sse(e);%求均方差n=mse(e);r=corroef(t,st);r2=r(1,2)2;fprintf(决定系数bata:%.4gn,r2);%自变量对因变量贡献值分析p=ones(1,ncase);p;w=net.b1 net.iw1;v=net.lw2;h=w*p;h=2./(1+exp(-2*h)-1;y=net.lw2*h+net.b2;c=zeros(s1,nvar);c=zeros(1,s1);for j=1:nvar for i=1:s1 pm1=cov(h(i,:),p(j+1,:),1); pm1=pm1(1,2); pm2=var(h(i,:),1); c(i,j)=p