CN110890130B - 基于多类型关系的生物网络模块标志物识别方法 - Google Patents

Abstract

基于多类型关系的生物网络模块标志物识别方法，属于组学数据分析技术领域。发明结合当前组学数据变量多，噪音多，生物体内变量以各种通路相互关联，变量之间关系复杂多样等实际特点，构建成对变量的多种组合变量，使用统一指标，对成对分子的各个组合进行类间区分能力评价，使用评价指标构建网络，并使用贪婪搜索技术搜索网络模块。有助于组学数据研究者迅速从基因组、代谢组等组学数据中识别具有强区分能力的变量模块，来进行后续针对分子功能，调控过程的定性定量研究，是一种基于多类型关系的生物网络模块标志物识别方法。

Description

基于多类型关系的生物网络模块标志物识别方法

技术领域

本发明属于组学数据分析技术领域，具体涉及一种基于多类型关系的生物网络模块标志物识别方法。

背景技术

代谢产物是基因表达的最终产物，在代谢酶的作用下生成。通过对生物体内的代谢物进行定性定量分析，能够挖掘代谢物在生理病理过程中的变化规律，代谢组学已被广泛地应用于临床标志物的发现、疾病的早期诊断和治疗及药物研发等领域。目前代谢组学研究涉及以下几个步骤：(1)样品制备；(2)代谢物检测(用高通量、高分辨率、高重复性的谱学手段，如核磁共振NMR技术、液质联用LC/MS技术、气质联用GC/MS技术等进行检测)；(3)峰匹配、基线矫正、积分、归一化、标度化；(4)模式识别和统计学分析技术，如PCA、PLS-DA、OPLS-DA等。

代谢通路是指细胞中的代谢物质在酶的作用下转化为新的代谢物质过程中所发生的一系列生物化学反应，代谢流的存在促成细胞的某种生理状态的出现。代谢网络则是指由代谢反应以及调节这些反应的调控机制所组成的描述细胞内代谢和生理过程的网络，代谢网络一直处于对环境的变动的响应之中，是生物体生理病理状态的综合反映，对代谢网络的研究有助于理解生物体代谢反应机制，以及疾病生物标志物的发现。在代谢网络的建模研究中，分子浓度及其转化速率尤为重要，代谢反应的情况可以通过研究代谢物含量间的相互关系得出，这种相互关系在代谢物含量间的组合运算上有所体现。

无论是原核细胞还是真核细胞，都有一套精确的基因表达和蛋白质合成调控机制。芯片技术实现了在同一时间测定大量基因的表达水平，从而进行功能富集，挖掘具有统计学意义的差异表达，或共表达基因的功能类别等分析。基因表达调控模式的分析，是基因组学数据挖掘的重要方向。基因表达调控是一个复杂的过程，基因和基因之间存在复杂的相互作用，形成网络结构，即基因调控网络。研究基因调控网络对于分析和研究生物体中的现象，寻找疾病生物标志物以及在生物分子学层面开发治疗靶基因等具有重要意义。

本发明从网络的角度分析基因组学、代谢组学数据，通过多种类型关系衡量组学数据变量间的关联性，构建多类型关系网络，并识别能够用于不同生理状态的样本判别的潜在网络模块标志物。该技术对每一对原始变量进行缩放后，构建四种类型组合变量，并分别对其进行类间区分能力评分，选取最高评分所对应的组合形式作为此对变量的最终组合形式，并将区分能力评分确定为网络中对应结点间连边的权重。以构建的网络中权重最小的边为搜索起点，逐步加入与当前模块相连，且分类性能提高最多的边，直至没有边使得分类性能提高为止，寻找具有判别能力的潜在网络模块标志物。

发明内容

本发明的目的是基于基因组学、代谢组学数据变量多、数据噪音大、变量间关系复杂多样并形成关联网络这一特点，采用一定的指标构建起能反映生物体内基因表达间或代谢物含量间相互关系的生物网络，并迅速定位具有判别能力的差异性组合变量模块，以此寻找能够反映生物生理状态发生发展的重要特征，筛选潜在网络标志物。该方法适用于组学数据的变量筛选和变量提取，可以用于代谢组学数据分析、基因组学数据分析以及转化医学等领域。该方法的核心技术为：基于多类型关系的融合生物网络构建与分析。

为了实现上述目标，本发明采用的技术方案如下：

一种基于多类型关系的生物网络模块标志物识别方法，步骤如下：

步骤一、获取数据

对于代谢组学分析：检测所研究问题的样本中的小分子代谢成分(相对分子量1000以下)，并确定这些成分在不同样本中的含量。将各种成分看作变量，将成分的含量看作变量的值。对于基因组学分析：测定不同样本中的microRNA或基因表达情况，将芯片上的各个信号看作变量，将信号强度数值看作变量的数值。

对数据集进行抽象：数据集含有若干样本，每个样本根据其所处生理病理状态的不同将其划分为不同的类别。每个样本拥有若干变量，样本在变量上具有一定的值。对数据集进行预处理：删除每类样本上缺失值数量均超过本类样本总数的百分之二十的变量，剩余变量的缺失值填补为同类样本在这个变量上的均值。F＝{f₁,f₂,…,f_m}定义为变量集合，m表示变量的个数；C＝{0,1}定义为类标集合，其中，0代表样本归属于0类，1代表样本归属于1类；S＝{s₁,s₂,…,s_n}定义为样本集合，n表示样本的个数。

步骤二、计算并筛选所有组合变量；

为了保证组合运算的有效性，对数据集的变量进行归一化，构成归一化变量集

变量最大最小归一化如公式(1)所示：

其中，

为变量i在样本j上经过归一化后的取值，f_ij为变量i在样本j上的原取值，

为变量i在所有样本上的原始最小值，

为变量i在所有样本上的原始最大值，min为转换后的下界，max为转换后的上界，此方法中分别设为1和2。经过调整后，缩放后的变量在所有样本上可能的最小值为1，最大值为2。

对于任意两个经过缩放的变量

1≤i<k≤m。加法、减法、乘法、除法四种组合方式定义如公式(2)-(5)：

按照公式(2)-(5)构建

四种新组合变量，通过单因素方差分析(ANOVA)技术，计算两类样本在四个新组合变量上的显著性差异p值，记为{p_ik ⁺,p_ik ^-,p_ik ^×,p_ik ^÷}，以此评价各个新组合变量的类间区分能力。

单因素方差分析技术常用于生物信息学数据统计分析，用于两个及两个以上样本均数差别的显著性检验。方差分析的基本原理是认为不同处理组的均数间的差别基本来源有两个，分别为实验条件(即不同的处理造成的差异，称为组间差异，记作SSA)，以及随机误差(如测量误差造成的差异或个体间的差异，称为组内差异，记作SSE)。组间差异SSA、组内差异SSE除以各自的自由度，得到其均方MSA和MSE。一种情况是处理没有作用，即各组样本均来自同一总体，MSA/MSE≈1。另一种情况是处理确实有作用，组间均方是由于误差与不同处理共同导致的结果，即各样本来自不同总体，那么，MSA>>MSE(远远大于)。通过将MSA/MSE构成F分布，用F值与其临界值比较，推断各样本是否来自相同的总体，对某变量x进行单因素方差分析的公式如(6)-(10)所示：

其中，m为不同类别的组数，s为总样本数量，s_d为第d组的样本数量，

为第d组的均值，

为所有样本的均值。x_dw为第d组第w个样本的值。利用上述公式求出的F值，通过查询统计学表可以得到变量x在两类样本所代表的总体上均值相同的概率p，若p<0.05，则认为均值相同是小概率事件，可以认为其均值不同，换言之此变量对于区分两类样本具有重要意义。

方法中，p值越小者将被认定为类间区分能力越强。对于一对变量

区分能力最强的组合变量将被认定为最佳组合变量，记作

同时将相应的显著性p值记作p_bestkindik，将区分能力最强的组合方式记为bestkind。全体最佳组合变量构成的集合记为F^best，对应的显著性p值构成的集合记为P^best，即有集合：

步骤三、网络构建与网络分析；

(3.1)构建生物网络。具体方式为：将变量集F^scaled中的变量看作生物网络中的结点，F^best中的组合变量看作用于计算组合变量的两个变量所对应结点的边，组合变量的类间差异显著性p值看作边的权重。若某条边权重大于等于0.05，则从网络中删除此边。经过筛选后，网络边集合、边权重集合与点集合分别为netEdge、netWeight、netNode，如公式(11)-(13)所示：

netNode＝F^scaled (13)

(3.2)差异模块搜索。具体方式如下：

(3.2.1)模块扩张种子边选取。设构成模块的结点集为moduleNode，

构成模块的边集为moduleEdge，

初始moduleNode＝φ，moduleEdge＝φ。从构建的生物网络中选取权重最小的边，不妨记为

作为模块扩张的起点，加入到模块边集合moduleEdge中，并将构成边的结点加入到模块点集合moduleNode中，加入后，当前模块结点集合

模块边集合

由步骤二可知，此边对应的组合变量是边结点对应的两个变量的四种组合中区分能力最强的，同时也是所有保留下来的组合变量中区分能力最强的。使用当前边所代表的组合变量构建临时数据集，设该临时数据集变量集合为

样本集合不变，仍为S＝{s₁,s₂,…,s_n}。

(3.2.2)模块搜索。当前模块分数记为moduleScore，初始moduleScore＝0。

(3.2.2.1)所有与当前模块相连的边构成边候选集candiEdge,

(3.2.2.2)逐个尝试将边候选集candiEdge中的边

加入模块中，1≤p<q≤n，则构成的临时数据集的变量集合变为

样本集合为S＝{s₁,s₂,…,s_n}。将临时数据集进行十倍交叉验证划分，针对每一倍产生的训练数据训练线性核SVM分类器，并在此倍交叉验证对应的测试数据上测试，将预测结果与真实结果对比，可以计算出每一倍上测试结果的受试者工作曲线(ROC)的曲线下面积(AUC)，此指标为机器学习分类技术常见指标，常用于评判分类器的分类性能，AUC值越高说明当前条件下分类器性能越好。对十倍交叉验证产生的共计十个AUC值求平均，作为此次加入的候选边的分数，记作score_p,q。删除临时数据集中加入的变量，恢复到加入之前的状态，即

继续尝试加入边候选集candiEdge中下一个候选边，按照上述相同的方式计算此候选边的分数。重复上述过程，直至遍历完边候选集candiEdge，构成候选边分数集合SCORE＝{score_p,q|f_pf_bestkindq∈candiEdge}。

(3.2.2.3)选取边候选集candiEdge中分数最高的边，记作

其中，

且score_best1,best2＝max(SCORE)。

(A)若score_best1,best2>moduleScore，则将此边加入到当前模块中，更新当前模块边集合

同时更新模块结点集合

更新临时数据集的变量集合

更新模块分数moduleScore＝score_best1,best2。重复步骤(3.2.2.1)-(3.2.2.3)。

(B)若score_best1,best2≤moduleScore，说明与当前模块相连的任意边均不能使得分类器性能有所提升，则停止模块更新，当前模块边集合moduleEdge以及模块点集合moduleNode为最终的识别结果。

本发明的有益效果：结合生物体内基因调控，代谢反应的特点，从网络的角度系统地分析基因组学、代谢组学数据。方法中，为消除变量的量级差异所带来的对组合变量有效性的影响，本发明使用缩放(scale)技术将变量缩放到相同的量级，同时考虑到生物体内各成分间关系的多样性，采用四种组合方式选最优的方式来刻画复杂的生理过程，并迅速找出差异性过程，拓宽了识别生物标志物的方法。基于公共miRNA数据以及公共代谢数据的结果表明，通过本发明所找出的差异性模块相比常用的单特征筛选方法，分类性能优越，变量集规模小。经过理论与实验分析，本发明能够为基因组学、代谢组学的研究提供切实有效的数据处理手段，具有较强的应用价值。

附图说明

图1为使用该方法，在结直肠癌代谢公共数据集上找到的网络模块。

图2为使用该方法，基于在结直肠癌代谢公共数据集上找到的模块，绘制的PLS-DA图。

图3为组合变量:赖氨酸/葡萄糖在两类样本上的箱线图。

图4为组合变量:己糖醛酸+赖氨酸在两类样本上的箱线图。

图5为组合变量:赖氨酸/肌醇在两类样本上的箱线图。

图6为组合变量:羟脯氨酸/肌醇在两类样本上的箱线图。

图7为组合变量:己糖醛酸/1_5-脱水山梨糖醇在两类样本上的箱线图。

图8为组合变量:羟脯氨酸/吡嗪-2,5-二醇在两类样本上的箱线图。

图9为组合变量:己糖醛酸/海藻糖1在两类样本上的箱线图。

图10为组合变量:谷氨基酸/吡嗪-2,5-二醇在两类样本上的箱线图。

具体实施方式

实施例：基于人类代谢的结直肠癌判别潜在网络标志物筛选。

(1)人类结直肠癌代谢谱数据

本例中所使用的人类结直肠癌代谢谱数据集为公共数据集，利用UPLC-QTOFMS和GC-TOFMS技术分别对脂肪组织的脂质组、代谢组和血清进行分析，经过QC处理后，数据集中包含物质99种，样本分为内脏脂肪(subcutaenous SAT)和皮下脂肪(visceral VAT)两类。两类样本数分别为57和49。此数据集的研究类型为生物标志物研究，研究目的包括(1)比较结直肠癌病人的VAT脂肪组织代谢组和脂质组和SAT脂肪组织代谢组和脂质组。(2)描述患者脂肪组织(VAT/SAT)中脂质组和代谢组之间的关系。

(2)计算并筛选所有的组合变量。

本例中变量共计99个，对于每对变量均构建4种组合变量，利用单因素分析公式(6)-(10)计算各个组合变量在两类上的类间差异显著性p值，组合变量个数共计(99×(99-1)/2)×4＝19404。

(3)网络构建与网络分析

对于每对变量的四种组合变量，仅保留单因素方差分析显著性p值最小者作为生物网络中对应边的权重，未筛选前网络边数共计99×(99-1)/2＝4851，将权重大于0.05的边从网络中移除，经过移除后，网络剩余边数为3226。网络构建完毕后，从中找出权重最小的边，此边对应的两个变量为葡萄糖、赖氨酸，对应运算方式为赖氨酸/葡萄糖，对应权重为2.35^-11。将此边作为模块搜索的起点进行搜索，寻找与当前模块相连，且使得十倍交叉验证下线性核SVM分类性能指标AUC的均值最大的边，若此边并未使得模块分数有所提升，则停止搜索，此数据集上停止搜索时，模块分数为0.9767，搜索结果为{(赖氨酸,葡萄糖,÷),(己糖醛酸,赖氨酸,+),(赖氨酸,肌醇,÷),(羟脯氨酸,肌醇,÷),(己糖醛酸,1_5-脱水山梨糖醇,÷),(羟脯氨酸,吡嗪-2,5-二醇,÷)}其中(赖氨酸,葡萄糖,÷)的含义为：将赖氨酸以及葡萄糖在所有样本上的值缩放到区间[1,2]，然后将各个样本在这两个缩放后的变量上的取值相除，形成新组合变量。如附图1所示，最终搜索出来的模块包含9种物质，8种组合变量。附图2为基于搜索出来的模块，在两类样本上构建的PLS-DA图，从图中可以看出，两类样本具有清晰的分离趋势，表明搜索出来的模块具有较强的区分能力。图3-图10为模块内各边所代表的组合变量上，两类样本的箱线图，图中可以看出，模块中各条边均具有强区分能力，信息丰富。

上表格为此方法(FusioNet)与常见的基于单变量的特征降维技术SVM-RFE、Lasso、ElasticNet方法在三个miRNA公共数据集，三个代谢公共数据集上30次hold-out的平均分类AUC比较结果(AUC(平均模块结点数/平均模块边数))，表格中加粗位置为数据集上几种方法的最高AUC值，使用分类器为线性核函数SVM分类器。从结果可以得出，本技术所确定的差异性模块具有较强的区分能力，且使用变量数远少于其余技术。

Claims (3)

1.基于多类型关系的生物网络模块标志物识别方法，其特征在于，步骤如下：

步骤一、获取数据

对于代谢组学分析：检测所研究问题的样本中的小分子代谢成分，所述小分子代谢成分的相对分子量为1000以下，并确定这些成分在不同样本中的含量；将各种成分看作变量，将成分的含量看作变量的值；

对于基因组学分析：测定不同样本中的microRNA或基因表达情况，将芯片上的各个信号看作变量，将信号强度数值看作变量的值；

对数据集进行抽象：数据集含有样本，每个样本根据其所处生理病理状态的不同将其划分为不同的类别；每个样本拥有变量，样本在变量上具有数值；

对数据集进行预处理：删除每类样本上缺失值数量均超过本类样本总数的百分之二十的变量，剩余变量的缺失值填补为同类样本在这个变量上的均值；将FS＝{f₁,f₂,…,f_m}定义为变量集合，m表示变量的个数；将C＝{0,1}定义为类标集合，其中，0代表样本归属于0类，1代表样本归属于1类；将S＝{s₁,s₂,…,s_n}定义为样本集合，n表示样本的个数；

步骤二、计算并筛选所有组合变量；

为了保证组合运算的有效性，对数据集的变量进行最大最小归一化，构成变量集

对于任意两个经过缩放的变量

加法、减法、乘法、除法四种组合方式定义如公式(2)-(5)：

按照公式(2)-(5)构建四种新组合变量

通过单因素方差分析ANOVA技术，计算两类样本在四个新组合变量上的显著性差异p值，记为{p_ik ⁺,p_ik ^-,p_ik ^×,p_ik ^÷}，以此评价各个新组合变量的类间区分能力；

方差分析的基本原理是认为不同处理组的均数间的差别基本来源有两个，分别为组间差异，记作SSA；以及组内差异，记作SSE；组间差异SSA、组内差异SSE除以各自的自由度，得到其均方MSA和MSE；一种情况是处理没有作用，即各组样本均来自同一总体，MSA/MSE≈1；另一种情况是处理确实有作用，组间均方是由于误差与不同处理共同导致的结果，即各样本来自不同总体，那么，MSA>>MSE；通过将MSA/MSE构成F分布，用F值与其临界值比较，推断各样本是否来自相同的总体，对某变量x进行单因素方差分析的公式如(6)-(10)所示：

其中，m为不同类别的组数，s为总样本数量，s_d为第d组的样本数量，

为第d组的均值，

为所有样本的均值；x_dw为第d组第w个样本的值；利用上述公式求出的F值，通过查询统计学表得到变量x在两类样本所代表的总体上均值相同的概率p，当p<0.05时，则认为均值相同是小概率事件，认为其均值不同，换言之此变量对于区分两类样本具有重要意义；

方法中，p值越小者将被认定为类间区分能力越强；对于一对变量

区分能力最强的组合变量将被认定为最佳组合变量，记作

同时将相应的显著性p值记作p_bestkindik，将区分能力最强的组合方式记为bestkind；将全体最佳组合变量构成的集合记为F^best，对应的显著性p值构成的集合记为P^best，即有集合

步骤三、网络构建与网络分析；

(3.1)构建生物网络；具体方式为：将变量集F^scaled中的变量看作生物网络中的结点，将F^best中的组合变量看作用于计算组合变量的两个变量所对应结点的边，组合变量的类间显著性p值看作边的权重；当某条边权重大于等于0.05，则从网络中删除此边；经过筛选后，网络边集合、边权重集合与点集合分别为netEdge、netWeight、netNode，如公式(11)-(13)所示：

netNode＝F^scaled (13)

(3.2)差异模块搜索；具体方式如下：

(3.2.1)模块扩张种子边选取；设构成模块的结点集为moduleNode，

构成模块的边集为moduleEdge，

初始moduleNode＝φ，moduleEdge＝φ；从构建的生物网络中选取权重最小的边，记为

作为模块扩张的起点，加入到模块边集合moduleEdge中，并将构成边的结点加入到模块点集合moduleNode中，加入后，当前模块结点集合

模块边集合

由步骤二可知，此边对应的组合变量是边结点对应的两个变量的四种组合中区分能力最强的，同时也是所有保留下来的组合变量中区分能力最强的；使用当前边所代表的组合变量构建临时数据集，设该临时数据集变量集合为

样本集合不变，仍为S＝{s₁,s₂,…,s_n}；

(3.2.2)模块搜索；当前模块分数记为moduleScore，初始moduleScore＝0；

(3.2.2.1)所有与当前模块相连的边构成边候选集candiEdge,

(3.2.2.2)逐个尝试将边候选集candiEdge中的边

加入到当前模块中，1≤p<q≤n,则构成的临时数据集的变量集合变为

样本集合为S＝{s₁,s₂,…,s_n}；将临时数据集进行十倍交叉验证划分，针对每一倍产生的训练数据训练线性核SVM分类器，并在此倍交叉验证对应的测试数据上测试，将预测结果与真实结果对比，计算出每一倍上测试结果的受试者工作曲线ROC的曲线下面积AUC，AUC值越高说明当前条件下分类器性能越好；对十倍交叉验证产生的共计十个AUC值求平均，作为此次加入的候选边的分数，记作score_p,q；删除临时数据集中加入的变量，恢复到加入之前的状态，即

继续尝试加入边候选集candiEdge中下一个候选边，按照上述相同的方式计算此候选边的分数；重复上述过程，直至遍历完边候选集candiEdge，构成候选边分数集合SCORE＝{score_p,q|f_pf_bestkindq∈candiEdge}；

(3.2.2.3)选取边候选集candiEdge中分数最高的边，记作

其中，

且score_best1,best2＝max(SCORE)；

(A)若score_best1,best2>moduleScore，则将此边加入到当前模块中，更新当前模块边集合

同时更新模块结点集合

更新临时数据集的变量集合

更新模块分数moduleScore＝score_best1,best2；重复步骤(3.2.2.1)-(3.2.2.3)；

(B)若score_best1,best2≤moduleScore，说明与当前模块相连的任意边均不能使得分类器性能有所提升，则停止模块更新，当前模块边集合moduleEdge以及模块点集合moduleNode为最终的识别结果。

2.根据权利要求1所述的基于多类型关系的生物网络模块标志物识别方法，其特征在于，步骤二中，最大最小归一化如公式(1)所示：

其中，

为变量i在样本j上经过归一化后的取值，f_ij为变量i在样本j上的原取值，

为变量i在所有样本上的原始最小值，

为变量i在所有样本上的原始最大值，min为转换后的下界，max为转换后的上界，此方法中分别设为1和2；经过调整后，缩放后的变量在所有样本上最小值为1，最大值为2。

3.根据权利要求1或2所述的基于多类型关系的生物网络模块标志物识别方法，其特征在于，步骤一中，检测所研究问题的样本中的小分子代谢成分所用方法为LC-MS或GC-MS；测定不同样本中的microRNA或基因表达情况所用方法为microRNA芯片或基因表达芯片。

Images