第八届华中地区大学生数学建模邀请赛承诺书我们仔细阅读了第八届华中地区大学生数学建模邀请赛的竞赛细则。我们完全明白,在竞赛开始后参赛队员不能以任何方式(包括电话、电子邮件、网上咨询等)与队外的任何人(包括指导教师)研究、讨论与赛题有关的问题。我们知道,抄袭别人的成果是违反竞赛规则的,如果引用别人的成果或其他公开的资料(包括网上查到的资料),必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。我们郑重承诺,严格遵守竞赛规则,以保证竞赛的公正、公平性。如有违反竞赛规则的行为,我们将受到严肃处理。我们的参赛报名号为:参赛队员(签名):队员1:队员2:队员3:武汉工业与应用数学学会第八届华中地区大学生数学建模邀请赛组委会第八届华中地区大学生数学建模邀请赛编号专用页选择的题号:B参赛的编号:(以下内容参赛队伍不需要填写)竞赛评阅编号:1第八届华中地区大学生数学建模邀请赛题目:基因调控网络的重构及病毒感染的致病机制【摘要】一个基因的表达受其他基因的影响,而这个基因又影响其他基因的表达,这种相互影响相互制约的关系构成了复杂的基因调控网络。基因调控网络的研究是从基因之间相互作用的角度揭示复杂的生命现象,是当前生物信息学研究的前沿。疾病的发病因素和原理,对于医疗领域有着十分重要的作用。这不仅仅能够让更多的患者免受病痛的困扰,还能促进人类医学史的进步。所以根据基因数据谱来重构基因调控网络,以及某个疾病症状产生的原因的研究具有很大的意义。本文对基因调控网络的重构以及导致严重临床症状的蛋白质进行了研究和推测。由于所给的基因数据谱(附录一)十分庞大,所以首先要对数据进行降维处理。本题基于时间序列给出了272组基因数据,为了减小噪声以及缺失值对实验精度的干扰,在实验前对四组噪声较大或有缺失的数据进行剔除。具体的降维方式采用了多元统计法中的主成分分析和聚类分析:先对这一万多个数据做主成分分析,从这一万多个数据中,通过线性变化选出了1000个左右的重要变量来组成新的样本。既降低了数据的处理难度,又尽量保持了新数据和原数据相比,尽可能保持原数据的信息。然后用spss两阶聚类法粗略地对要聚类的数目进行一个估计,根据此估计用K-means算法对数据进行处理,得到相应的30组数据。对这30组数据建立模型,来重构基因调控网络。本文中采用的模型是线性回归模型,并对它的合理性,以及相较贝叶斯网络作了对比。最后依据所得到的系数矩阵进行基因网络图的绘制与呈现。问题二在第一问的基础上,寻求导致产生严重临床症状的蛋白质。根据附录二给出的个体出现感染症状时间节点示意图,1代表此志愿者在该时间节点表现出了临床症状,0则表示没有,这是一个二分类。本题采用逻辑回归模型,利用LR分类器模型去寻找该重要蛋白质。用268组数据,其中每一个基因视为该组数据的一个属性,对这些基因进行LR分类,并得到相应的系数矩阵。然后对系数矩阵进行分析,取出影响比较大的几个基因,然后对照基因表对基因作用的描述去寻求该重要蛋白。本题最终找出四个导致志愿者产生严重的临床症状的蛋白质。所有代码实现,以及每次得到的系数矩阵均在附录中给出。关键词:线性回归模型,基因调控网络重构,多元统计法,主成分分析,聚类分析,逻辑回归(LR分类器模型)21.问题重述通过基因之间的相互调控,生物体可以实现细胞的生长,器官的发育、以及免疫等各种生物机能。随着测序技术的发展,产生了越来越多的高通量实验数据。基于这些实验数据重建基因调控网络(Generegulatorynetworks,GRNs),对于深入了解生物机能的实现过程具有重要作用。生物实验中,在17个健康志愿者鼻内接种流感病毒H3N2/Wisconsin,其中9个人出现了严重的感染症状,另外的8个人没有出现症状。接种后,每隔大约8h从血液中采集样本测量基因表达谱数据,实验数据一共有16个时间点(单位:h),包括baseline(-24),0,5,12,21,29,36,45,53,60,69,77,84,93,101,108,共268个样本。基因表达谱数据见附件1,其中前8个为未出现严重感染症状的数据,后9个为出现严重感染症状的数据。(其中行代表探针号,对应着不同的基因;列为各个个体血液样本在各个时间节点的数据)个体出现感染症状的时间节点示意图见附件2。问题:1)根据实验数据重构基因调控网络;2)通过比较出现感染症状的志愿者和健康志愿者的样本数据,试确定病毒感染人体后导致志愿者是否会出现严重临床症状的重要蛋白。2.问题分析一个基因的表达受其他基因的影响,而这个基因又影响其他基因的表达,这种相互影响相互制约的关系构成了复杂的基因调控网络。更一般些,几乎所有的细胞活动都被基因网络所控制。生命是存储并加工信息的复杂系统,孤立地研究单个基因及其表达往往不能确切地反映生命现象本身的内在规律。因此,需要从复杂系统的角度研究基因网络。对于问题一,考察我们如何根据已有的基因表达谱(附录一)去重构基因调控网络,从而推断调控网络各节点之间潜在的调控关系。考虑“反向分析法”来重构基因调控网络,常见的基因调控网络模型有布尔网络模型、线性组合模型和贝叶斯网络模型等等。然而题目所给的数据集十分庞大,如果直接将这一万个基因全部带入模型,那么计算量是惊人的。所以需要用到多元统计方法中的主成分分析和聚类分析去实现降维的操作。对于问题二,在已经重构好的基因网络的基础上寻找导致病毒感染人体以后导致志援者是否产生严重临床症状的蛋白质。首先我们要对数据进行分析,寻找与染病相关系数大的基因,然后依据附录一的sheet2中对于基因的描述去进一步确定关键蛋白质。33.模型假设针对本问题,建立如下合理假设:(1)题目所给数据准确可靠;(2)假设不考虑个体差异性;(3)基因表达呈高斯分布;4.符号说明mn,X表示第n个基因基于时间序列的第m组数据;KA表示一个基因;i为回归系数;Xt,a代表基因X在时间点t具有的表达值;321,,bbb为常数;21,为误差项。5.问题一的建模与算法实现求解5.1数据的分析问题一需要根据所给的基因表达谱数据来重构基因调控网络,附录一中的sheet1中给出了17个志愿者体内的10000种基因,随着注入病毒后的时间变化而出现的数值变化。由于数据集过大,所以第一步要做的就是对这一万种基因进行筛选降维操作。只选取部分具有代表性的数据代入模型,从而减少计算量。对于数据的处理部分,采用多元统计中的常用方法,主成分分析和聚类分析。5.2数据预处理5.2.1数据处理方法选择由于这道题目的数据量庞大,所以,如何筛选数据就成了很重要的一步。我们这里采取先对10000组数据做主成分分析,形成1000组新变量,再对这些新变量进行聚类分析,进一步降维。5.2.2主成分分析(1)主成分分析的基本思想:主成分分析的基本思想是通过构造10000个基因初始数据的适当的线性组合,以产生一系列互不相关的新变量,从中选出少数几个新变量并使它们尽可能多地包含原先所有基因的信息(降维),从而使得用这几个新变量替代原变量分4析问题成为可能。即在尽可能少丢失信息的前提下从所研究的m个变量中求出几个新变量,它们能综合原有变量的信息,相互之间又尽可能不含重复信息。(2)主成分分析的实现:设有n个样品,m个变量(指标)的数据矩阵。本题中n=10000,表示10000种基因;m=268,表示基于时间序列的基因数据变化指标。(1)11121(2)21222()12mmnmnnnnmxxxxxxxxXxxxx寻找k个新变量12,,,()kyyykm,使得1、1122,(1,2,,)llllmmyaxaxaxlk2、12,,kyyy彼此不相关主成分的系数向量12(,,,)llllmaaaa的分量lja刻划出第j个变量关于第l个主成分的重要性。可以证明,若12(,,,)Tmxxxx为m维随机向量,它的协方差矩阵V的m个特征值为120m,相应的标准正交化的特征向量为12,,,muuu,则12(,,,)Tmxxxx的第i主成分为(1,2,,)Tiiyuxim。称1/mijj为主成分(1,2,,)Tiiyuxim的贡献率,11/kmjjjj为主成分12,,kyyy的累计贡献率,它表达了前k个主成分中包含原变量12,,,mxxx的信息量大小,通常取k使累计贡献率在85%以上即可。当然这不是一个绝对不变的标准,可以根据实际效果作取舍,例如当后面几个主成分的贡献率较接近时,只选取其中一个就不公平了,若都选入又达不到简化变量的目的,那时常常将它们一同割舍。计算步骤如下:1、由已知的原始数据矩阵nmX计算样本均值向量12ˆ(,,,)Tmxxxx;其中11(1,2,,)niijjxximn2、计算样本协方差矩阵1ˆ()()ˆ1ijijVsn其中1()()(,1,2,,)nijliiljjlsxxxxijm53、把原始数据标准化,即ijjijjjxxx,记()nmijXx。形成样本相关矩阵ˆTRXX;4、求ˆR的特征根120m及相应的标准正交化的特征向量12,,,muuu,可得主成分为(1,2,,)Tiiyuxim。(3)主成分分析降维结果用Matlab实现以上算法(代码见附录),实现结果如下:图5.1主成分分析结果如图可见是一个1000组新的变量,由于数据集比较大,在这里只截出一部分。下面再对这1000组新变量做聚类分析处理。5.2.3聚类分析(1)聚类分析的基本思想:聚类(clustering),简单的讲就是将一个给定的数据集分成若干个不同簇的过程。聚类算法中的簇指的是数据对象的集合,且这种数据对象集合必须满足条件:同一簇中的数据对象间具有较大的相似性,而不同簇中的数据对象间具有较小的相似性。聚类的主要指导思想就是尽可能使同一簇内对象相似度达到最大,且不同簇间对象相异度达到最大。(2)K-means算法:6首先从含有n个数据对象的数据集中随机选择K个数据对象作为初始中心然后计算每个数据对象到各中心的距离,根据最近邻原则,所有数据对象将会被划分到离它最近的那个中心所代表的簇中,接着分别计算新生成的各簇中数据对象的均值作为各簇新的中心,比较新的中心和上一次得到的中心,如果新的中心没有发生变化,则算法收敛,输出结果,如果新的中心和上一次的中心相比发生变化,则要根据新的中心对所有数据对象重新进行划分。直到满足算法的收敛条件为止。(3)K-means算法的实现:从含有1000个数据对象的基因表达谱(附录一)中随机选择100个数据对象作为初始的聚类中心;两基因数据的相似度可通过计算两个基因数据的欧式距离来得到,再根据最近邻原则将数据对象逐个划分到离其最近的聚类中心所代表的簇中,计算误差平方和准则函数E的值;更新聚类中心,即分别计算各个簇中所有数据对象的均值作为各个簇的新的中心,以新的聚类中心来计算误差平方和准则函数E的值;○4将步骤计算得到的E值和前一次计算得到的E值进行比较,若两者差值的绝对⑤值小于等于预先设定的阈值,即聚类准则函数收敛,则转步骤(5),否则转步骤;⑤输出K个聚类。5.2.4数据降维后的结果显示将10000个基因数据通过主成分分析和聚类分析缩小到30组特征数据。图5.2各类中基因数目分布75.3模型一的分析目前重构基因调控网络的模型主要有贝叶斯网络,线性组合模型等,本题我们采用的是线性回归模型,下面先对基因调控网络与线性回归模型的对应关系作分析,以证明模型的可行性:假设基因kA受基因,2K1A,Ak.....kmA调控,则认为相应的基因表达值,有以下线性关系:kkmtmktktktabababa,12,121,11,...式中kta,代表基因KA在t时刻的表达水平,1b....mb为常数。若基因A,B,C之间的真实调控关系如图所示,其中A,B,C代表基因,而边代表调控关系,比如从A到B有一条有向边,代表了基因A对基因B有调控关系。A