确定肿瘤的重要基因信息

全国第七届研究生数学建模竞赛 题 目 确定肿瘤的重要基因信息 摘 要： 随着生物分子学的发展，人们已经发现癌症与基因之间存在密切关系。 本文通过对比基因表达谱中不同样本的表达水平的差别，区分出无关基因和信息基因，在信息基因的空间中搜索分类能力强的特征子集,然后将所有样本划分为训练集和测试集，使用支持向量机SVM和人工神经网络检查特征子集的分类能力。其次，再考虑噪声的定义，以及如何去除噪音，并分析噪声对特征子集分类能力的影响。最后，生成特征子集的时候需要考虑已有的医学发现，提出基于知识库的基因图谱分析模型KFS模型，有效利用了医学成果。 第一问，利用巴氏距离模型和理想基因模型区分无关基因和相关基因，剔除基因表达谱中无关基因，达到降维的效果。 第二问，使用FSSM算法在信息基因空间里寻找出候选特征子集，然后将样本划分为训练集和测试集，使用支持向量机SVM和人工神经网络，对FSSM搜索出来的特征子集的分类能力进行判定。本文得到由5个信息基因组成的特征子集，样本分类正确率达到95.79%。 第三问，考虑了阈值滤波和主成分分析两种去噪模型，并阐述噪音模型在高斯过程分类器的构建中的作用，最后论述噪音能够在学习算法中防止过渡拟合从而可以孵化出泛化能力更强的分类器用于确定基因标签。 第四问，提出基于知识库的基因图谱分析模型KFS模型，该模型在引入信息基因知识库的基础上，对基因图谱进行去噪处理、样本评价函数增益、剔除无关基因，并采用基于知识库的KFSSM算法获得特征子集，最后分别采用SVM及ANN方法获取信息基因集合。本文最后对已知临床经验的结肠癌数据进行处理得到一组信息基因组合，样本分类正确率达到94.52%。 关键词：基因表达谱，信息基因，巴氏距离，FSSM，噪声，KFS模型 一、问题重述 癌症起源于正常组织在物理或化学致癌物的诱导下，基因组发生的突变，即基因在结构上发生碱基对的组成或排列顺序的改变，因而改变了基因原来的正常分布（即所包含基因的种类和各类基因以该基因转录的mRNA的多少来衡量的表达水平）。所以探讨基因分布的改变与癌症发生之间的关系具有深远的意义。 DNA微阵列（DNA microarray）,也叫基因芯片，是最近数年发展起来的一种能快速、高效检测DNA片段序列、基因表达水平的新技术。它将数目从几百个到上百万个不等的称之为探针的核苷酸序列固定在小的（约1 ）玻璃或硅片等固体基片或膜上，该固定有探针的基片就称之为DNA微阵列。根据核苷酸分子在形成双链时遵循碱基互补原则，就可以检测出样本中与探针阵列中互补的核苷酸片段，从而得到样本中关于基因表达的信息，这就是基因表达谱，因此基因表达谱可以用一个矩阵或一个向量来表示，矩阵或向量元素的数值大小即该基因的表达水平（见附件）。 随着大规模基因表达谱 ( Gene expression profile ,或称为基因表达分布图)  技术的发展，人类各种组织的正常的基因表达已经获得，各类病人的基因表达分布图都有了参考的基准，因此基因表达数据的分析与建模已经成为生物信息学研究领域中的重要课题。如果可以在分子水平上利用基因表达分布图准确地进行肿瘤亚型的识别，对诊断和治疗肿瘤具有重要意义。因为每一种肿瘤都有其基因的特征表达谱（见附图）。从DNA 芯片所测量的成千上万个基因中，找出决定样本类别的一组基因“标签”,即“信息基因” （informative genes ）是正确识别肿瘤类型、给出可靠诊断和简化实验分析的关键所在，同时也为抗癌药物的研制提供了捷径。 通常由于基因数目很大，在判断肿瘤基因标签的过程中，需要剔除掉大量“无关基因”，从而大大缩小需要搜索的致癌基因范围。事实上，在基因表达谱中,一些基因的表达水平在所有样本中都非常接近。例如，不少基因在急性白血病亚型（ALL,AML）两个类别中的分布无论其均值还是方差均无明显差别,可以认为这些基因与样本类别无关,没有对样本类型的判别提供有用信息,反而增加信息基因搜索的计算复杂度。因此,必须对这些“无关基因”进行剔除。1999 年《Science》发表了Golub 等针对上述急性白血病亚型识别与信息基因选取问题的研究结果[1] 。Golub 等以“信噪比”(Signal to noise ratio) 指标作为衡量基因对样本分类贡献大小的量度,采用加权投票的方法进行亚型的识别,仅根据72个样本就从7 129 个基因中选出了50 个可能与亚型分类相关的信息基因。Golub 的工作大大缩小了决定急性白血病亚型差异的基因范围，给出了亚型识别的基因依据，富有创造性。Guyon 等则利用支持向量机的方法再从中选出了8个可能的信息基因[2]。 但信噪比肯定不是衡量基因对样本分类贡献大小的唯一 标准 excel标准偏差excel标准偏差函数exl标准差函数国标检验抽样标准表免费下载红头文件格式标准下载 ，肿瘤是致癌基因、抑癌基因、促癌基因和蛋白质通过多种方式作用的结果，在确定某种肿瘤的基因标签时，应该设法充分利用其他有价值的信息。有专家认为[3]在基因分类研究中忽略基因低水平表达、差异不大的表达的倾向应该被纠正，与临床问题相关的主要生理学信息（见问题4）应该融合到基因分类研究中。 面对提取基因表达谱信息这样前沿性课题，命题人根据自己科学研究的经历和思考，猜测以下几点是解决前沿性课题的有价值的工作。这种猜测是科学研究中的重要环节，当然猜测不会总是可行的，更不一定总是正确的。但不探索就不能前进，如果能够通过数学建模，得到的部分结果可以佐证你们的猜测或为新探索提供若干依据，就很有价值。我们的目的只是给研究生以启发，鼓励研究生培养这样的创造性发现的能力。所以研究生完全可以独立设计自己的技术路线，只要能够有效提取附件的基因表达谱信息就行。 1、由于基因表示之间存在着很强的相关性，所以对于某种特定的肿瘤，似乎会有大量的基因都与该肿瘤类型识别相关，但一般认为与一种肿瘤直接相关的突变基因数目很少。对于给定的数据（见附件），如何从上述观点出发，选择最好的分类因素？ 2、相对于基因数目，样本往往很小，如果直接用于分类会造成小样本的学习问题，如何减少用于分类识别的基因特征是分类问题的核心，事实上只有当这种特征较少时，分类的效果才更好些。对于给定的结肠癌数据如何从分类的角度确定相应的基因“标签”？ 3、基因表达谱中不可避免地含有噪声（见1999 年Golub在《Science》发表的文章），有的噪声强度甚至较大，对含有噪声的基因表达谱提取信息时会产生偏差。通过建立噪声模型，分析给定数据中的噪声能否对确定基因标签产生有利的影响？ 4、在肿瘤研究领域通常会已知若干个信息基因与某种癌症的关系密切，建立融入了这些有助于诊断肿瘤信息的确定基因“标签”的数学模型。比如临床有下面的生理学信息：大约90%结肠癌在早期有5号染色体长臂APC基因的失活，而只有40%~50%的ras相关基因突变。 二、基本的模型假设 1、基因表达谱中的样本类别没有错误。 2、基因表达谱中的癌症病人样本都是结肠癌病人。 3、基因表达谱中有重复的基因标签，我们假定重复的基因标签的样本分类能力类似，所以只处理其中一个基因标签 4、所有的癌症病人都处于同一时期，不分早期和晚期。 5、基因表达谱中有一些重复的基因标签，比如HSAC07、UMGAP和i都出现了4次，而Has.13491、Has.44472等基因出现了两次。我们假定每一个标签都具有代表性，所以我们只处理一个基因标签，而不考虑其他重复的基因标签。这样的话，基因表达谱信息表中总共有1911个不同类型的基因标签。本文剩下所有的数据处理都是针对这1911个基因标签。 三、名词解释 基因表达谱：关于基因表达的信息，可以看成是一个矩阵或者一个向量，矩阵或者向量元素的数值大小就是该基因的表达水平。 基因表达水平：可以理解为样本中某种基因的数量或者密度。 理想基因：一种可以完全分辨出样本类别的基因，它在不同类型的样本中的表达水平相差很大。 信息基因：可以决定样本类别的一组基因。 无关基因：有一些基因的表达水平在所有的样本中非常接近，对样本分类没有帮助的基因。 特征子集：每一种信息基因的组合。 信噪比：作为衡量基因对样本分类贡献大小的量度。 训练集：用来训练分类器学习能力的样本集，包含正常人的样本和结肠癌病人的样本。 测试集：用来测试分类信息基因分类能力的样本集，包含正常人的样本和结肠癌病人的样本。 分类器：通过学习训练集中样本之后，可以自动的对给定的测试集中样本进行分类的一类程序。 四、符号化 Bi：基因标签i的巴氏距离（1≤i≤1911）。 SN：无关基因集合。 SI ：信息基因集合。 e：理想基因，与肿瘤有很强的关联性。 Num(SN)：无关基因集合的大小。 Num(SI)：信息基因集合的大小。 正常样本nSi：第i个正常人样本（1≤i≤22）。 病人样本cSj：第j个结肠癌病人样本(1≤j≤40)。 基因表达谱信息：A[m×n]，用一个m×n的矩阵来表示基因，其中m=62，n=1911 vij：样本i在基因j上的表达水平。 基因矢量Vg：基因g在各个样本上的表达水平的一个矢量，Vg={vg1, vg2,…, vgn}。 五、模型的建立与求解 5.1 问题一的分析、建模与求解 5.1.1 问题一的分析 目前人们通过生物芯片的技术可以快速检测样本的基因表达水平，人类各类组织的正常基因表达谱已经获得，但是还没有有效的方法能够定位与癌症直接有关的基因。 问题一的出发点在于，人类基因表达谱中包含有太多与癌症无关的基因，这大大地增加了人们从基因表达谱中搜索与癌症相关的信息基因的难度，而一般情况下，人们认为直接与特定类型癌症相关的突变基因数目很少，所以我们应该考虑首先从基因图谱中删除大量的无关基因，缩小搜索信息基因的范围。这个 步骤 新产品开发流程的步骤课题研究的五个步骤成本核算步骤微型课题研究步骤数控铣床操作步骤 可以称为基因表达谱去除无关信息的降维处理。本题我们从下面两个角度对基因表达谱初步降维： 1、 信息基因在不同类型样本表达水平的差异。 2、 基因与理想基因的相似度。 我们首先从癌症的“基本致病机理”角度分析肿瘤与基因之间的可能存在关系。癌症会导致信息基因在不同类型样本的表达水平上产生一些差异。我们应当用一种合理的指标将这种差异量化。目前比较通用的方法是比较样本的平均值和方差，我们考虑一种模型，可以综合考虑这两方面的因素。 其次，我们假设有一种理想基因，这种基因在不同类型上的样本上的表达水平差异非常大。我们通过样本在理想基因上的表达水平就可以直接判断样本的类型。利用这种理想基因，我们通过比较基因表达谱中各个基因与理想基因的相似度。相似度高的基因可以认为是信息基因，相反，如果某个基因标签与理想基因的相似度很小，我们基本可以认为它是无关基因。 本文在处理第一题时，综合利用了巴氏距离模型和理想基因模型，以巴氏距离模型为主，但是由于基因表达谱中的噪声会影响巴氏距离模型选出来的信息基因的效果，所以再以理想基因模型为辅，选取一些与理想基因相似度高的基因，防止将一些信息基因剔除。最后选取大小为250的信息基因集合，大约占题目所给基因总数的20%，作为第二问的特征子集的搜索空间。 下面详细介绍我们建立的降维模型。 5.1.2 问题一的模型建立 我们首先考虑下如何量化信息基因在不同类型样本中表达水平的差异以及如何利用这种差异将基因分类，区分出无关基因和信息基因。 Golub等人以“信噪比”(Signal to noise ratio)[1]指标作为衡量基因对样本分类贡献大小的度量，信噪比的定义如下： （1） 其中：d是基因的信噪比，μ1和μ2分别是该基因在两种样本中表达水平的均值，σ1和σ2是该基因在两种样本中表达水平的标准差。 但是使用这种方法来区分信息基因和无关基因存在应的问题。如果d=0，该基因就会被当做无用基因删除，而实际上，如果该基因在两种样本中表达水平方差有很大差异，那么很有可能这个基因与癌症有很密切的关系。所以，我们需要选用一个模型，综合考虑平均值和方差的差异。 5.1.2.1巴氏距离模型 巴氏距离既考虑到基因在样本中均值，也考虑到基因在样本中的方差分布，是一个很好的信息度量指标。它综合考虑了均值和方差差异对样本分类的作用。它的定义如下： （2） 其中B为基因的巴氏距离。由式(2) 知，巴氏距离由两部分构成：第一项体现了基因在两个类别中分布均值的差异对样本分类的贡献；第二项体现了分布方差的不同对分类的贡献。依据该距离公式,即使基因在两类不同样本中分布的均值相同，只要分布的方差出现大的差异，仍然可以获得较大的距离值 [3]。 从模式分类的角度来看，基因的巴氏距离越大，说明该基因的分类能力越强，基因的分类信息越多。设SN是无关基因集合，SI是信息基因集合，我们设置一个阈值θ，巴氏距离大于θ可以认为是无关基因，巴氏距离小于θ的可以认为是信息基因。 （3） 其中，g是基因，B(g)为基因g的巴氏距离，θ是选取的巴氏距离的阈值。利用公式3，选取好的阈值，我们就可以区分出无关基因集合SN以及信息基因集合SI。 5.1.2.2 理想基因模型 当然我们不能只从巴氏距离这一个标准来衡量基因分类信息的大小，同时由于生物基因芯片不可避免地存在一些噪声（噪声的处理方法会在第三问中提出解决 方案 气瓶 现场处置方案 .pdf气瓶 现场处置方案 .doc见习基地管理方案.doc关于群访事件的化解方案建筑工地扬尘治理专项方案下载 ），噪声会对样本的均值和方差产生较大的影响，而均值和方差是计算基因巴氏距离的两个重要因素。所以为了使基因的分类更为合理，除了以巴氏距离，我们还提出理想基因的概念。 题目中给定的基因表达谱数据可以用矩阵A[m×n]表示，其中m表示基因的数目，去除重复的基因之后，只有1911个，n表示样本的数量，总共有62个。矩阵中元素vij表示第j个样本对基因i的表达水平。我们首先对基因表达谱的数据作归一化处理，使得矩阵A中每个元素的值都在[-1, 1]之间： （4） 其中vmax是指矩阵A中的最大值，vmin是指矩阵A中的最小值， 是vij归一化之后的数值。 本题中的样本总数有62个，分为两类，一类是正常人样本，另一类是结肠癌患者样本。正常人样本总共有22个，样本编号从N1—N2，结肠癌患者样本有40个，样本编号从C1—C40。基因g在每个样本中表达水平Vg={vg1, vg2,…, vgn}可以看成是一维向量。我们将理想基因e定义为： （5） 如果基因g是结肠癌的信息基因，它本身携带有分类信息越多，分类能力越强，那么它越接近于理想基因e。我们从两个方面来考虑基因g向量与理想基因e的接近程度，一个是基因向量之间的夹角的余弦值，如果余弦值靠近-1或者1，那么这两个基因向量夹角越小。此外，这两个向量之间的欧拉距离也可以作为一个衡量的标准，两个基因向量之间的欧拉距离越小，说明这两个基因向量越靠近。最后，我们用相似度这个值来量化信息基因与理想基因之间的接近程度[8]。 （1）基因g与基因e的夹角定义为 : （6） 其中 （2）基因g与基因e的欧拉距离为： （7） （3）基因g与基因e的相似度为： （8） 从相似度的定义中，我们可以看出，基因g与理想基因e向量的欧拉距离不变，夹角越小，余弦值的绝对值越大，相似度越大。基因g和e的夹角不变，欧拉距离越小，相似度越大。相似度定义很好地量化了基因g与理想基因e之间的相关性，相似度越高，说明基因g的分类能力越强。同样的，我们可以通过设置适当的阈值θ，将基因表达谱中信息基因和无关基因区分开，达到降维的效果。 （9） 5.1.2.3 综合模型 本文在处理第一题时，综合利用了巴氏距离模型和理想基因模型，以巴氏距离模型为主，但是由于基因表达谱中的噪声会影响巴氏距离模型选出来的信息基因的效果，所以再以理想基因模型为辅，选取一些与理想基因相似度高的基因，防止将一些信息基因剔除。 我们首先计算基因的巴氏距离前200的基因，然后再计算与理想基因的相似度，选取相似度值前50，并且不与前面重复的基因，组成大小为250的信息基因集合。这个基因集合大约占题目所给基因总数的20%，大幅压缩了冗余基因。这250个基因集合作为第二问的特征子集的搜索空间。 图5.1.1 综合模型的分类基因 流程 快递问题件怎么处理流程河南自建厂房流程下载关于规范招聘需求审批流程制作流程表下载邮件下载流程设计 图 下图是综合模型得到的信息基因集合的构成。 图5.1.2 综合模型的分类基因组成 5.1.3 问题一的模型求解 1、计算所有基因标签的巴氏距离。 根据巴氏距离计算公式，我们得到了1911个基因的巴氏距离分布情况，见表2。 表5.1.2 基因标签的巴氏距离分布 巴氏距离 基因个数 百分比 0~0.05 1492 78.04% 0.05~0.1 290 15.18% 0.1~0.15 78 4.08% 0.15~0.2 31 1.62% 0.2~0.25 13 0.68% 0.25~0.4 7 0.36% 图5.1.3 基因标签的巴氏距离分布直方图 2、计算剩余基因与理想基因的相似度 接下来，我们基因表达谱中所有基因与理想基因的相似度，然后取相似度前50，且不与巴氏模型的结果重复的基因。 表5.1.3 基因的相似度分布表 相似度 基因数目 百分比 0.0-0.02 638 33.38% 0.02-0.04 934 48.87% 0.04-0.06 291 15.23% 0.06-0.08 43 2.25% 0.08-0.10 5 0.26% 图5.1.4 基因标签的相似度分布直方图 5.1.4 问题一的结果及分析 问题一的主要目的在于剔除与肿瘤无关的基因，通过巴氏距离模型和理想基因模型处理之后，基因的数量从原来的1911个大幅度地减少到了250个，降低了第二问FSSM算法的搜索特征子集的空间。 5.2 问题二的分析、建模与求解 5.2.1 问题二的分析 问题一和问题二其实都是对基因降维，问题一从单个基因分类能力的角度出发，剔除癌症无关基因，缩小了人们搜索与癌症相关的信息基因的范围。而问题二则是建立在问题一的基础上，从剩下的250个信息基因中搜索特征子集，而这250个基因可以组成2250个不同的特征子集，这就需要一个很高效的搜索算法，同时也需要有一个合理的评价函数，能够评价不同的特征子集的分类能力，从而筛选出分类能力强的特征子集。最后将样本分成训练集和测试集，再选择合适的具有学习能力的分类器，查看候选特征子集的分类能力。特征子集大小和分类准确率为评价指标可以作为衡量指标。 5.2.2基于FSSM算法的特征子集的生成 本文采用FSSM(Floating Sequential Search Method)搜索算法[2][3]，对特征子集所构成的子空间进行搜索，从中选取30个具有不同维数的待选分类特征子集，然后使用SVM和人工神经网络检验这30个特征子集的分类能力。 FSSM搜索算法中采用关键函数J作为动态搜索特征子集过程中的评价函数，评价函数的描述如下[3]： （10） 其中，Fi表示含有i个信息基因的特征子集，μ1 、μ2表示特征子集Fi中的信息基因在正常样本和结肠癌样本中的均值向量，∑1、∑2表示Fi中的信息基因在正常样本和结肠癌样本中数据的的协方差矩阵，本文中∑1为i×22矩阵，∑2为i×40矩阵。 J函数也是基于巴氏距离的，不过计算的是信息基因集合的巴氏距离。第一问我们利用巴氏距离模型，计算出单个基因的巴氏距离，很好地区分出信息基因和无关基因。FSSM中J函数的作用在于，它从同样大小的信息集合中选择J值较大的，具有最强的分类能力子集代表。 FSSM算法中，令数组Fmax[i]表示含有i个信息基因的具有最大评价函数值的特征子集，本文采用matlab语言实现了FSSM算法，数组SelectMax [i]表示计算过程中，计算出含有i+1个基因的最大特征子集时选择i个的特征子集的最大J值，参考文献[2]中的算法思想给出matlab算法的处理过程如图2.1所示。 图 5.2.1 FSSM算法matlab实现的算法流程 通过运行FSSM算法，最终生成30个信息基因的特征子集，本文的下一节分别通过支持向量机方法(SVM)和人工神经网络方法考察选出的特征子集的分类能力。 5.2.3.1 基于支持向量机（SVM）的分类特征子集选择 本文的上一节采用FSSM算法生成了34个具有不同维数的特征子集，该部分以支持向量机为分类器对34个特征子集进行样本识别，最终获得具有最大分类正确率的基因组合。 支持向量机是一种基于统计学习理论，采用结构风险最小化原理的机器学习算法。机器学习的目的是根据给定的训练样本对输入输出之间的依赖关系的估计，使得可以对未知的输出尽可能准确的预测。支持向量机的核心思想就是调整评价函数使得最好地利用边界样本点的分类信息，从而构造出最佳分类超平面,因此支持向量机可以获得很好的泛化能力并且与样本的具体分布无关[7]。 结合基因图谱数据特点发现，支持向量机可以有效地处理高维样本的分类问题，计算复杂度受样本维数的影响较小，适合处理小样本、高维数的基因表达谱数据的样本分类问题。 本文使用34组具有不同维数的特征子集依次作为特征属性使用支持向量机学习出分类模型后验证其分类能力。由于样本实例的数目有限，我们采用10-fold交叉验证的方法来评估学习得到的模型的分类能力，进而评估选择的特征基因组合的识别能力。具体方法是将数据集分成10份，轮流将其中9份做训练1份做测试，10次的结果的均值作为对算法精度的估计。 5.2.3.2 基于人工神经网络（ANN）的分类特征子集选择 人工神经网络是一种应用类似于大脑神经突触联接的结构进行信息处理的数学模型. 在这一模型中, 大量的节点之间相互联接构成神经网络，以达到处理信息的目的。人工神经网络对矢量进行样本识别前需要进行训练，训练的过程就是应用一系列输入矢量，通过某种算法逐步调整权值和阈值的过程，通过训练人工神经网络对一组输入矢量产生希望的输出。训练后的人工神经网络即可以用于对正常样本和结肠癌样本的分类。 本文也采用了基于人工神经网络的分类器，使用34组具有不同维数的特征子集依次作为特征属性使用支持向量机学习出分类模型后验证其分类能力，进而评估选择的特征基因组合的识别能力。可以与SVM的结果进行对比。 5.2.4 问题二的结果及分析 本题我们首先使用FSSM方法从250个信息基因中搜索，生成34个不同大小的特征子集，子集大小从2到35不等。然后在62个样本集合上使用SVM（支持向量机）和人工神经网络算法训练出分类器，然后评估分类器的分类能力，进而评估所用基因特征子集识别 下面给出经过SVM100次10-fold交叉验证之后，样本分类准确率在前五的特征子集以及它们的分类准确率。最好的一组特征子集的样本分类能力达到了95.79%，并且该组特征子集的大小只有5。符合我们之前对特征子集的分类能力以及子集大小的要求。 表5.2.1 特征子集分类准确率前五名 特征子集 大小 分类准确率 Hsa.37937, Hsa.710, Hsa.3016, Hsa.5392, Hsa.6080 5 95.79% Hsa.37937, Hsa.710, Hsa.3016, Hsa.5392 4 94.32% Hsa.37937, Hsa.549, Hsa.710, Hsa.3016, Hsa.5392, Hsa.6080 6 91.10% Hsa.37937, Hsa.549, Hsa.710, Hsa.3016, Hsa.5392, Hsa.6080, Hsa.2058, Hsa.43331, Hsa.8214, Hsa.823, Hsa.957, Hsa.33965, Hsa.816, Hsa.490, Hsa.732, Hsa.36689, Hsa.2928, Hsa.8147, Hsa.6814, Hsa.2250 20 89.56% Hsa.37937, Hsa.549, Hsa.710, Hsa.3016, Hsa.5392, Hsa.6080, Hsa.2058, Hsa.43331, Hsa.8214, Hsa.823, Hsa.957, Hsa.33965, Hsa.816, Hsa.490, Hsa.732, Hsa.36689, Hsa.2928, Hsa.8147, Hsa.6814, Hsa.2250, Hsa.7048, Hsa.58 22 89.27% 图5.2.2 SVM下不同特征子集的样本识别情况 为了验证最优特征子集的准确性，我们使用Matlab工具绘制了Hsa.37937，Hsa.710，Hsa.3016，Hsa.5392，Hsa.6080这五个基因在62个样本中的表达水平。说明一下，下面的五张图中，*代表正常人样本的基因表达水平，+代表癌症病人样本的基因表达水平。 图5.2.3 基因标签Has.37937在不同样本的表达水平 通过这张图，我们可以明显地看出Has.37937基因在正常人和结肠癌病人两类样本的表达水平差异非常明显。 图5.2.4 基因标签Has.710在不同样本的表达水平 这张图显示基因Has.710的分类效果其实并不理想，无论从均值还是从方差来看，都不显著。它被选进特征子集的原因既有可能是噪声的影响，这点我们会在第三问中继续讨论。 图5.2.5 基因标签Has.3016在不同样本的表达水平 基因标签Has.3016的在两类样本中的表达水平差异也很大。 图5.2.6 基因标签Has.5392在不同样本的表达水平 图5.2.7 基因标签Has.5392在不同样本的表达水平 同样的，我们也利用人工神经网络20次循环验证特征子集的分类能力，下表给出了样本分类准确率在前三的特征子集以及它们的分类准确率。通过与前面的SVM的样本分类能力对比，我们发现{Hsa.37937, Hsa.710, Hsa.3016, Hsa.5392, Hsa.6080}这组特征子集的分类能力确实最强，所以确定结肠癌的基因标签就是这组特征子集。 表5.2.2 人工神经网络中特征子集分类准确率前三名 特征子集 大小 分类准确率 Hsa.37937, Hsa.710, Hsa.3016, Hsa.5392, Hsa.6080 5 93.54% Hsa.37937, Hsa.549, Hsa.710 , Hsa.3016, Hsa.5392, Hsa.6080, Hsa.43331, Hsa.8214, Hsa.823 9 93.54% Hsa.37937, Hsa.549, Hsa.710 , Hsa.3016, Hsa.5392, Hsa.6080 6 91.94% 5.3 问题三的分析、建模与求解 5.3.1 问题三的分析 使用DNA微阵列（DNA Microarray）技术测量基因表达水平得到的数据具有噪声强、波动大的特点，同时在大量数据的背后还有很多相关变量不能被直接观察到[4]，所以我们有必要仔细研究表达谱中的噪音，并给出相应的理论观点，处理思路和方法。 首先，从提取基因表达谱的角度看，噪声可以分为两类：一种噪声可以认为测量误差引入的噪声，在任何测量过程中无法避免的；另一种噪声是无关基因，该类基因在正常样本和结肠癌样本中的表达水平非常接近，没有为肿瘤的判断提供有用的信息，该类基因的存在增加了提取信息基因的难度。 无关基因是本文研究的问题之一，已经在第一问中试图解决了。下面主要关注第一种噪音。这类噪音主要是测量误差引起的，我们知道测量误差可以分为以下三类[5]： 1. 系统误差，主要是由于测量设备的缺陷、测量环境变化、测量时使用的方法不完善、所依据的理论不严密或采用了某些近似公式等造成的误差。 2. 随机误差，在同一测试条件下，多次重复测量同一量时，误差大小、符号均以不可预定的方式变化着的误差上。 3. 疏失误差，是指在一定的测量条件下，测得的值明显偏离其真值，既不具有确定分布规律，也不具有随机分布规律的误差。 从误差的定义可以知道，系统误差和随机误差是必然存在的。而且在DNA微阵列技术中系统误差是比较大的。随机误差是时刻存在的，且其服从一种自然分布。系统误差与随机误差的划分是相对的，二者在一定条件下可以相互转化，即同一误差，既可以是系统误差，又可以成为随机误差。疏失误差是由于测试人员对仪器不了解、或因思想不集中、粗心大意导致错误的读，使测量结果明显地偏离了真值。 对这三种误差的处理方法不同：对于含有疏失误差的测量值应予以剔出；对于随机误差的影响用统计的方法来消除或减弱；对于系统误差则主要靠测量过程中采取一定的技术措施来削弱或对测量值进行必要的修正来减弱其影响。结合到基因图谱信息提取这项具体技术中，我们无法求证系统误差的大小，只能结合到系统误差和随机误差的关系将其中的一部分作为随机误差来处理。而疏忽误差和设备老化、环境突变等恶劣因素引发的系统误差导致的是异常数据，会影响到建模方法的有效性，所以是我们需要在数据中剔除或者修正的。 下面我们针对疏忽误差和随机误差引发的噪声分别建立模型。 5.3.2 问题三的模型建立 我们对两种噪声分别讨论与之对应的模型： 图5.3.1 去噪模型 5.3.2.1异常数据滤波模型和主成分分析去噪模型 去除噪声数据是数据预处理的一项基本处理过程。针对从基因表达谱这类特殊的数据中发现基因标签这类特点问题，我们讨论一般的滤波去噪方法并提出一种主成分分析的去噪模型。 考虑到基因表达谱中的数据主要是来源于在癌症患者和正常人的基因上通过微观实验和生物化学原理等复杂技术获取到的。鉴于微观实验操作难度，系统误差可能导致突兀数据。同时考虑到癌症这种本身就是生理及其异常的疾病，即使是患有同一种癌症的病人也可能存在某些局限于个别人的基因突变的情况。考虑到研究中主要是发现问题的一般规律，我们做如下假设： 1. 某些癌症病人的基因突变引发的数据突变视为异常数据 2. 由于数据具有高噪声、波动大的特点，将系统误差和疏忽误差引起的突兀数据作为异常数据 3. 数据中一定存在噪声，主要是由测量的随机误差和固有的系统的误差共同决定，并且噪声数据符合高斯分布。 1）异常数据滤波模型 我们建立去除异常数据带来的噪声的模型。主要采用两种方法：均值滤波和中值滤波。通过观察，我们发现基因表达谱数据中存在突兀的数据项，即基因的样本数据严重偏离该样本的均值，下面统一称为异常数据。我们分别使用均值和中位数作为数据的参考基准，实验中设定阈值并通过调整其大小来设定过滤异常数据的滤波窗口的大小。 均值滤波算法过程如下： Step 1. 设定滤波窗口阈值a和调整数据比例阈值b=5%； Step 2. 对每个正常人和癌症病人分别计算Step3； Step 3. 计算每个基因维度数据的均值，并根据窗口阈值调整落在窗口外面的奇异的到窗口边缘； Step 4. 统计被调整的奇异数据项的数量，并计算其在整个数据中的比例。调整比例如果约为5%则终止程序，否则跳转到Step 1按照一定步长调整a的值。 中位数滤波算法类似，不再赘述。下面主要从基因谱数据挖掘基因标签这个具体应用出发构建主成分分析去噪模型。 2）主成分分析去噪模型 主成分分析是一种采用组合特征的方法将多维数据降维的方法。方法的目标是寻找在最小化重构误差的意义下最能够代表原始数据的投影方法。降维后的数据能够比较好的代表原始数据。主成分分析的主要思想是： 1. 用一维向量表示高维样本 2. 将一维投影量扩展到相对低维的空间 3. 低维空间是由高维空间数据的散布矩阵的最大几个特征值向量构成 通常情况下高维空间数据的散步矩阵的最大几个特征值占据了特征之和的绝大部分，所以可以认为少数几个最大特征值对应的特证向量即可表示原数据中 的绝大部分信息，而剩下的小部分（即对应较小的特征值的特征向量所表示的信息），通常可以认为是数据噪声而丢掉。考虑到基因谱数据的高噪音、多异常、大波动的特点，我们可以通过主成分分析的方法在降低维度的过程中去除噪音数据。 主成分分析一般是对样本的特征属性维度进行降维，在基因谱数据中即对应于基因维度。一方面由于主成分分析降维中是将当前的维度空间映射到低维空间，映射后将会当前维的多个维度映射到低维空间某个维度上，从而降维后的数据维度并不对应于某个当前维度；另一方面考虑到基因的维度在选取特征基因问题要求下需要是不能在降维去噪的过程中被“坏掉”的。所以我们考虑从样本的维度使用主成分分析方法进行降维去噪。 我们将数据集合按照样本的种类分成多个数据集合（这里是两个类别）分别进行主成分分析，这样降维后的每个数据项可以看成一个新的样本，并且其类别保持降维前所属的类别。 5.3.2.2 随机噪音的高斯模型 由随机误差和固有的系统误差引入的噪音有很好的随机性，一般假设其符合某种自然分布，其中以高斯分布最为普遍。下面介绍两种使用概率方法基于噪音模型的建模方法。 下面首先介绍一种对噪音建模的高斯过程分类器（Gaussian Process Classification）[6]，其基于贝叶斯（Bayes）理论的概率学习算法，使用高斯过程模型对噪音建立模型效果很好。 图5.3.1 GPC图模型示意图[6] 高斯过程分类器（简称为GPC）可以看作是一个图模型（如图5.3.1），用随机变量表示输入、潜在变量表示函数值和类别标签。潜在函数值完全决定类别标签。有很多噪音模型用来建模类别标签的似然函数。该问题可以如下作形式化定义： 仅考虑二分类问题，假定有数据集 其中 表示实例， 表示两类标签。在训练数据集上，我们希望训练得到一个分类模型能通过计算新实例 的可能性最大的所属类别 。 高斯过程分类器的核心思想是假设可以通过一些与 关联的且真正有价值的隐藏变量 来判别出实例的类别 。下面阐述从贝叶斯框架的角度建立GPC模型的主要步骤。 首先，我们在函数 上设定一个先验概率，即给定一个有限集合 ，随机向量 服从高斯分布。不是一般性我们假设过程的期望为0，且有 和 的协方差为： 其中超参数 反映了隐藏变量方差的垂直波动， 表示隐藏变量离0期望的偏置， 表示隐藏噪音的方差（ 表示噪音，期望为0，方差为 ）， 表示第m个特征属性在模型中的贡献权重。隐藏函数值 服从多元高斯分布： 其次，确定似然值为给定样本和隐藏函数值的结合，是似然函数的乘积： 我们假定隐藏函数值被高斯噪音影响，并且和输入独立。其中考虑了高斯噪音后的似然函数为： 第三，可以得到后验概率： 第四，根据文献[6]中的推导可以预测给定的实例 的类别 的分布： 其中 是假设发现的最优超参数。 至此建立了贝叶斯框架下的基于噪声模型的高斯过程分类器。 文献[6]中使用Expectation Propagation算法求解GPC模型，并在colon cancer数据集（和本文实验的数据集十分类似）上面进行试验，得到如下结果： 表5.3.1 2000个基因的测试错误率 从这个结果可以看出基于噪音模型的GPC模型可以取得比较好的结果。 5.3.3 问题三的讨论 从上面的论述可以看出，噪音在建立优秀的分类器，尤其是基于概率模型的分离器中起到了很重要的角色。下面我们将从分类器的泛化能力的角度，阐述一下噪音扮演的重要作用。 我们在确定基因标签的过程中主要是要判别基因标签对癌症的识别能力，通常的方法是使用基因组合作为特征属性训练出一个分类器，并通过分类器的分类能力来鉴别基因组合是否为基因标签（参考第一、第二问的解决方案）。 而数据中存在噪音是不可避免的，同时也是有其优点的。在使用机器学习的方法训练一个分类器的过程中，我们需要避免的一个问题就是分类器对训练数据的过渡拟合。 过渡拟合：是训练获得的模型过于符合训练数据的特点，而泛化能力较弱，从而在未知的新的测试数据集上面的分类效果一般，甚至错误较多。 这样噪音的存在可以再很大程度上扰动了分类器的容忍能力，使他不可能完全拟合训练数据，从而对新的测试样本有比较好的识别能力。而且在确定肿瘤的基因标签的背景之下，未知类别的测试样本的数据异常情况比较多见，这时候分类器的泛化能力显得尤为重要。所以基于噪音较大的数据发现出来的基因标签的对这些异常癌症样本的识别能力在一定程度上得到保证。 5.4 问题四的分析、建模与求解 5.4.1 问题四的分析 问题一、二、三的基因图谱分析模型都是以基因图谱的统计数据为基础进行分析的，而基因图谱信息不可避免的含有噪声，而这些噪声会影响特征肿瘤信息基因的确定。实际在肿瘤的研究领域，根据临床经验会已知若干个基因与某种癌症的关系密切，因此将包含临床经验的知识库融入到基因图谱模型中更加有利于癌症信息基因的确定。 题目中已知信息临床生理学信息：大约有90%结肠癌在早期有5号染色体长臂APC基因的失活，而只有40%~50%的ras相关基因突变。根据这个信息可以建立结肠癌的知识库，在基因表达谱的分析上可以重视知识库中的基因，提高其重要性参数，因此得到的肿瘤信息标签在具有很好的结肠癌判别能力的同时，更加尊重了实际临床数据的重要性。 该部分首先分析我们提出的基于知识库的基因图谱分析模型（Knowledge-based FSSM VSM ，KFS模型），然后结合题目中给定的结肠癌数据计算结肠癌的信息基因，并分析KFS模型与问题二的基因图谱分析模型间处理结果的分类能力比较，该部分最后给出结果分析以及出现该结果的原因，并提出了本模型还需要解决的问题。 5.4.2 问题四的模型建立 基于上述对问题四的分析，本文提出一种基于知识库的基因图谱分析模型(KFS模型)，该算法引入信息基因知识库的概念，模型首先基于临床经验建立模型的知识库，对去噪后的数据进行样本评价函数增益，然后对剔除基因图谱中的无关基因，降维后的数据通过KFSSM（Knowledge-based FSSM）算法的处理得到分类特征子集空间，然后模型分别采用改进的支持向量机（ISVM）以及人工神经网络（IANN）对分类特征子集的分类能力进行考察，最终确定癌症的信息基因组合，KFS模型的结构图如图5.4.1所示。 图 5.4.1 基于知识库的基因图谱分析模型 5.4.2.1 癌症基因知识库 题目中提到肿瘤研究领域通常会已知若干个信息基因与某种癌症的关系密切，那么在基因图谱分析的中引入已知关系必然会提高肿瘤信息基因判别的实际准确率和有用性。本文提出一种肿瘤信息基因知识库概念，用于存储临床经验已经发现的肿瘤信息基因。 结合题意，本文提出的信息基因知识库具有可以表示如下： [基因名称， 基因变化，基因样本统计概率] （5.1） 其中，基因名称表示与该类肿瘤的信息基因名称（本文算法中采用编号标识）；基因变化分为突变和失活（本文算法中取值分别为1、-1），基因突变表示基因表达水平值明显上调，基因失活表示基因表达水平值明显下降，表现在基因图谱中是基因的表达水平出现较大波动；基因样本统计概念表示该基因在样本统计中发生变化（突变或失活）的概率。 假设知识库的一项 记录 混凝土 养护记录下载土方回填监理旁站记录免费下载集备记录下载集备记录下载集备记录下载 为[gi , 1 , k%]，参考表达式5.1可知该项记录意思是信息基因gi在临床经验中有k%的样本表现出基因突变特征。 引入信息基因知识库的概念就是为了帮助从基因信息图谱中判断出肿瘤的信息基因组，考虑信息基因与癌症的密切关系，本文提出的信息基因知识库在KFS模型中的作用主要有三个： 1）知识库在数据预处理中的作用 由于基因表达谱中不可避免的存在噪声等因素的干扰，并且由于基因图谱中样本数目相对于基因数目往往很少，那么即使是信息基因也可能因为噪声数据而没有被选为信息基因，因此需要知识库来修正这个结果。 由于知识库样本数据的多样性以及临床数据的可供参考性，参考知识库中的信息基因修正相应的样本值：通过增加该基因的评价参数（巴氏距离），则可以使得该基因表现的更像信息基因。 2）知识库在分类特征子集生成过程中的作用 FSSM算法的特征子集空间搜索过程中的评价函数是以Fi的Bhattacharyya距离为评价函数，这忽视了知识库中信息基因的重要性，在比较selectMax[i]与Fmax[i]的更新问题上就需要判断是否需要更新，具体算法见KFSSM算法的描述。 需要说明的是，虽然信息基因知识库中基因具有更重要的作用，但是不会出现在所有的特征子集中，这也是符合算法执行的情况的，具体分析见问题四的结果及分析部分。 3）知识库在分类器中的作用 本文同样将知识库的作用考虑进了支持向量机以及人工神经网络分类器中，在机器学习的过程中考虑知识库中基因的重要性，提升了知识库中基因在学习过程中的重要性。 综上所述，知识库的本质意义在于对基于基因图谱数据的处理过程中起到一个引导作用，使得结果在某种程度上偏向于具有较高参考价值的信息基因知识库。 5.4.2.2 基因图谱数据预处理 由于一般基因图谱中，样本相对于基因数目往往很少，如果直接用于分类会造成小样本的学习问题，因此需要对基因图谱的原始数据进行预处理，本文KFS模型中的图谱数据预处理主要包含三个功能：噪声处理、基于知识库的样本评价函数增益和无关基因剔除，该部分试图通过对原始基因图谱数据进行分析得到初步过滤的特征子集，其结构如下图5.4.2所示。 图 5.4.2 KFS模型中预处理结构图 1）噪声去除 由于在读取生物芯片时一些不可控因素导致了某些基因样本表达水平发生了很大的变化，因此信息基因的信号有可能被噪声淹没，因此在基因图谱数据分析前需要进行噪声处理。 本文提出的KFS模型中的噪声去除直接采用了问题三建立的噪声模型，去除了基因表达谱中的噪声数据，为后续的基因图谱分析提供较好的数据。 2）基于知识库的样本评价函数增益 本文采用Bhattacharyya距离作为衡量基因中蕴含的分类信息量的评价函数，基因的Bhattacharyya距离越大，该基因对于样本的分类能力就越强，其被选为信息基因的可能性就越大。 基于信息基因知识库的分析，知识库中的信息基因对判断正常样本和癌症样本的区分能力很强，体现在样本的评价函数上就是需要大的Bhattacharyya距离，因此本文对基因的Bhattacharyya距离做简单的增益，以增加基因图谱数据分析过程中该基因被选为信息基因的可能性。 假设[gi, x, k%]为知识库中的一个信息基因记录，增益后的评价函数为： Bhattacharyya(gi) = Bhattacharyya(gi) +α* （ k%* maxB） 其中，增益参数经α的选择取决于实际基因图谱数据,maxB表示样本中最大的Bhattacharyya距离值。经过该增益后，知识库中的基因在基因图谱中的相对表达水平获得一定的增益，增益的程度取决于知识库中基因样本统计概率。 因此在通过设定Bhattacharyya距离阈值θ来剔除冗余基因的过程中，由于知识库中的基因在选择过程中获得了相应的增益，所以其被选进信息基因集合SI的可能性就增加了。 3）无关基因剔除 本文的KFS模型中的无关基因的剔除采用了Bhattacharyya距离，在上一节增益的基础上，通过设定阈值θ来生成无关基因集合SN和信息基因集合SI。 5.4.2.3 基于知识库的FSSM的分类特征子集生成算法 通过KFS模型预处理过程，我们得到了一组约为样本规模10%的信息基因集合SI，本文提出一种基于知识库的FSSM算法（KFSSM）算法，对基因集合SI构成的特征子集空间进行搜索，最后得出34组具有不同维数的分类特征子集。 本文提出的KFSSM算法的动态搜索过程中，也采用了Bhattacharyya距离作为评价函数，同时基于知识库引导KFSSM算法的搜索过程。KFSSM算法的基本原理是在顺序搜索的过程中，当计算出由i+1个基因的构成的最大评价函数的基因集合Fmax[i+1]时，回溯的时候计算出其中i个基因构成的最大评价函数SelectMax[i]，此时在刷新Fmax[i]的时候需要考虑知识库中的信息基因，对于包含知识库中信息基因权重较大的特征子集最后被选中的可能性较大，特征子集在知识库中的权重函数的设定设定如下： 其中，KDB表示信息基因知识库，k表示基因的样本统计概率，gi表示特征子集Fi中包含的知识库中的基因。 算法评价的基本流程如图5.4.3所示，其中Fmax[i]表示含有i个基因的特征子集的且具有最大评价函数的基因组合，selectMax[i]表示从当期i+1个基因中选择i个取其评价函数最大的特征子集，B(Fmax[i])函数表示求特征子集的Bhattacharyya距离，W()函数表示求相应特征向量在知识库中的权重值。 图 5.4.3 KFSSM算法搜索过程原理图 本文采用matlab程序对KFSSM算法进行了实现，下面列出算法的伪代码如下： 设置分类特征子集空间大小target值 %target = 35 读取信息基因集合SI和知识库数据kdbArray 初始化Fmax[2]数据 %直接取B值最大的两个向量 while i < target SI =SI - Fmax[i]; %从候选基因中删除已选基因 计算Fmax[i+1] = max{Fmax[i],g} ; %选择Fmax[i]与SI中基因拥有最大评价函数的组合 j = i; %用于控制回溯 while true 计算selectMax[j]; %从Fmax[i+1]中选择i个基因组合的具有最大评价函数的组合 if B( selectMax[j] ) > B( Fmax[j] ) && W(selectMax[j]) > W(Fmax[j]) Fmax[j] = selectMax[j] ; %更新Fmax[j]值 if j==2 break; %回溯更新执行到Fmax[2]停止 else j = j - 1; end else break; % 不需要回溯更新 end end end 循环结束，获得target组不同维数的分类特征子集 5.4.2.4 改进型的SVM和ANN分类器 本文的KFS模型中的支持向量机和人工神经网络方法使用问题二的模型，本文同样将知识库的作用考虑进了支持向量机以及人工神经网络分类器中，并对其做了简要修改，在学习的过程中考虑知识库中基因的重要性，提升了知识库中基因在学习过程中的重要性。 知识库在支持向量机和人工神经网络的运用是本文需要改进的地方，后续工作是实现基于知识库的分类器模型。 5.4.3 问题四的模型求解 文章结合给出的基因图谱信息以及结肠癌的临床数据，利用本文提出的KFS算法处理基因表达谱数据，下面描述处理过程。 5.4.3.1 结肠癌信息基因知识库的建立 题目中指出临床经验：大约90%结肠癌在早期有5号染色体长臂APC基因的失活，40%~50的ras相关基因突变，因此根据这些已知关系可以建立结肠癌信息的知识库。 搜索文章中给出的基因图谱找出APC基因的位置是project_data1.xls文件第869行的Hsa.2238基因，经过去除重复基因得到该基因在本文算法中的编号823，因此可得关于 Hsa.2238基因的记录：[Hsa.2238 , -1 , 90% ]; 搜索ras相关基因可得其两个记录分别为Hsa.2974和Hsa.2868 ，它们在本文算法中对应的编号分别为1027和1044，因此记录分别为：[Hsa.2974 , 1, 40%~50% ]和[Hsa.2868, -1, 40%~ 50%]。 为方便本文算法计算，文章中的基因名称替换成在本文算法中的基因编号，因此根据题意，建立起关于结肠癌信息基因的