CN111739582B - 一种基于协同作用网络的生物组学数据分析方法 - Google Patents

Abstract

本发明提供一种基于协同作用网络的生物组学数据分析方法，属于生物组学数据分析技术领域。本发明针对基因组学、转录组学与代谢组学数据维数高、样本量小的特点，通过构建协同作用网络，来充分挖掘组学分子间的协同互补关系，基于模块搜索的方式来选定重要模块，从而确定潜在的模块标志物。本发明为基因组学、转录组学与代谢组学数据的前期分析处理提供了切实有效的方法，对后续的医疗诊断具有较强的应用价值。

Description

一种基于协同作用网络的生物组学数据分析方法

技术领域

本发明属于生物组学数据分析技术领域，通过分析评价基因组学、转录组学或代谢组学数据特征间的协同作用，利用信息论的技术，从网络层面来挖掘组学数据重要特征，从而确定与疾病相关的模块标志物，对后续医学研究具有重要意义。

背景技术

随着人类基因组计划的完成，高通量测序技术的日臻完善，这极大的促进了基因组学、转录组学数据的发展，近些年来，代谢组学成为继基因组学、转录组学后新近发展起来的热门学科，成为生物组学数据重要组成部分之一。生物组学数据通常具有维数较高、样本量较小的特点。受实验条件、实验设备等限制，往往会导致数据所含噪音杂质较多，因此如何从生物组学大数据中挖掘出重要的信息，对疾病诊断、药物研发、个性化治疗等具有重要的意义。

目前大多数方法只停留在分子层面来筛选生物标志物，而生物系统是一个复杂的网络，代谢反应机制十分复杂，分子以各类通路反应的形式相互关联，相互作用，共同反映生命活动变化，因此基于网络方式来筛选模块标志物是尤为必要的。生物分子间往往会存在协同交互作用，两个区分能力较弱的分子联合所提供的交互信息可能与类别高度相关。协同作用网络通过计算信息增益得分，充分考虑分子与分子的关联关系，对具有协同交互作用的分子之间构建连边。本发明在综合分析基因组学、转录组学与代谢组学的基础上，提出了一种基于网络方式的模块标志物筛选技术即MSIG算法。该技术通过构建协同作用网络，综合考虑分子间的交互作用，在所构建的网络上搜索特征模块，从而来确定潜在的模块标志物。本发明基于协同作用网络来筛选具有重要生物意义的模块标志物，在生物学领域中有着广泛的实际应用价值。

发明内容

本发明的目的是在基因组学、转录组学、代谢组学大数据中提取具有生物学意义的模块标志物，对后续医疗诊断提供实际应用价值。本发明提出了一种基于协同作用网络确定模块标志物的数据分析方法。该方法的核心技术为：对于二分类问题，令F＝{f₁,f₂,…,f_m}表示特征集合，对于集合内的每一对特征f_i，f_j∈F(1≤i≠j≤m)，计算信息增益值IG，该指标可以衡量组学数据分子之间复杂的线性关系和非线性关系，只保留具有显著性IG的特征对来构建协同作用网络。在所构建的网络上采用贪心策略搜索模块标志物，取网络中与类标互信息值最大的节点作为初始种子节点，依次检索在网络中与当前模块所含节点相邻且不属于该模块的节点集作为候选节点集合，综合衡量每一候选点与当前模块协同作用的紧密程度以及自身所包含的与类标的互信息值的大小，将网络的拓扑结构与边的权值相结合评价每一个候选节点的重要程度，每次尝试将得分最高的节点加入到当前模块，并计算模块的整体分类性能是否较之前有所提高，若提高，则将其加入并继续检索当前模块的邻接节点集合，直到模块的性能不再增加或无邻接点为止，则该模块搜索结束并从网络中删除其所包含的节点及其连边。重复此搜索过程，直至网络中所剩节点数量为空。最后根据模块AUC值对所得模块降序排序，利用内部交叉验证的方式，选取准确率最高的若干模块，作为所选模块标志物。

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

一种基于协同作用网络的生物组学数据分析方法，步骤如下：

将高通量测序技术所得的基因组学或转录组学数据中的基因序列视为特征，基因序列在不同样本上的表达值视为特征值；利用LC-MS或GC-MS技术得到代谢组学数据中的小分子代谢成分(相对分子量1000以下)，将实验所得的小分子代谢成分视为特征，小分子代谢成分在不同样本上的含量视为特征值。

令F＝{f₁,f₂,…,f_m}为原始特征集合，包含m个特征；X＝{x₁,x₂,…,x_n}代表样本集合，n为样本总数；Y＝(y₁,y₂,…,y_n)为n个样本的类标向量，其中y_d∈C＝{c₁,c₂},d∈{1,2,…,n}，C为类标集合；

(1)计算特征f_i与f_j的信息增益得分IG(f_i；f_j；Y)，1≤i≠j≤m，公式如下：

IG(f_i；f_j；Y)＝I(f_i,f_j；Y)-I(f_i；Y)-I(f_j；Y)(1)

其中I(f_i,f_j；Y)表示特征f_i、f_j与类标向量Y的联合互信息大小，I(f_i；Y)、I(f_j；Y)为特征f_i、f_j与类标向量Y的互信息；IG(f_i；f_j；Y)反映f_i(f_j)与类别Y的依赖关系由于f_j(f_i)的加入所发生的变化；若IG(f_i；f_j；Y)>0，则表示特征f_i与f_j联合所提供关于Y的信息量比它们各自与Y的信息量之和大，进而认为f_i与f_j具有协同作用；

(2)获得具有显著性的信息增益得分；

当信息增益得分接近于零时，需确定具有显著性交互增益的特征对；对此保留步骤(1)大于零的信息增益值，并对其升序排列，得到

利用公式(2)计算组距：

再计算序列IG⁺中相邻两个变量的差值{p₂-p₁,p₃-p₂,…,p_t-p_t-1}，只保留结果大于θ₁的差值变量，并对其求平均值记为θ₂,将满足条件

的最小下标为i₀所对应的

值设定为边的阈值，如果

表明f_i与f_j具有较强的协同作用，IG值越大表明协同作用越强，因此在具有较强协同作用的特征对之间构建连边，此时网络中的节点代表特征，边的权值即为IG值；

(3)在所构建的协同作用网络上进行模块搜索；

由上述步骤(2)所得的协同作用网络，执行模块搜索的过程，过程如下：

(3.1)初始化：特征模块集合

特征模块子集

互信息集合

候选节点集合

(3.2)计算每个特征与类标向量的互信息，加入集合MIList；

(3.3)取MIList中互信息最大的点作为初始种子节点，加入到特征模块子集S中，并计算此时模块受试工作特征曲线ROC下的面积AUC值大小；

(3.4)定义候选节点集合CS＝Adj(S)-S，其中Adj(S)为当前特征模块子集S内所含节点在网络中的邻接点集，若当前特征模块子集S的候选集合为空集，则该模块搜索结束，继续向下执行步骤(3.5)；若候选集合不为空，则计算候选集合内的每个候选节点得分，公式为：

其中

分别表示点f_i、f_j在网络中的一近邻节点集，

则代表点f_i与f_j属于S相同的一近邻节点数量，|S|代表当前模块中节点个数，

与|S|相除表示f_i与f_j属于S相同的一近邻节点数量占当前模块节点总数的比例，若该比例越大且f_i与f_j构成的连边权值越大，则表明f_i与f_j连接越紧密，从而有效的衡量了候选节点与当前模块中节点连接的紧密程度。MI_i为特征f_i的互信息值，代表特征自身的分类性能，因此公式(3)将网络拓扑结构与特征自身分类性能相结合，综合衡量了每一个候选节点的重要性，得分越高表明该节点越重要，进而每次尝试将综合得分最高的节点加入到当前特征模块子集S中，检查模块的AUC值是否增加，若增加，则将该点加入到S中，继续重复步骤(3.4)的搜索过程，直至模块AUC值不再增加为止，此时该模块搜索结束，继续向下执行步骤(3.5)；

(3.5)S即作为本次搜索所得到的模块加入到特征模块集合MS中，从MIList中删除S所包含的特征，并从网络中删除S所包含的特征及其中的连边，以保证每个特征属于且仅属于一个特征模块，将S置空；

(3.6)重复上述步骤(3.3)-(3.5)的搜索过程，直至网络中节点数量为空，此时所有模块搜索结束；

(4)根据特征模块集合MS，确定最终的模块标志物；

对步骤(3)集合MS中所得模块，根据AUC值降序排列且只保留AUC值大于等于0.7的模块,利用内部交叉验证的方式，计算每加入一个新模块后已选集合的准确率，选择准确率最高的模块集合作为最终所确定的模块标志物。

本发明的有益效果：

本发明针对基因组学、转录组学与代谢组学数据维数高、样本量小的特点，通过构建协同作用网络，来充分挖掘组学分子间的协同互补关系，基于模块搜索的方式来选定若干重要模块，从而确定潜在的模块标志物。本发明为基因组学、转录组学与代谢组学数据的前期分析处理提供了切实有效的方法，对后续的医疗诊断具有较强的应用价值。

附图说明

图1为本发明的整体流程图。

图2为协同作用网络图。

具体实施方式

下面结合技术方案在模拟数据上进一步说明本方法MSIG的具体实施方式，样本集X＝{x₁,x₂,…,x₁₀}，共有10个样本，特征集为F＝{f₁,f₂,…,f₅}，共5个特征。

(1)计算对特征的信息增益值：首先将数据在每一维特征上进行等宽离散化，划分为5等份，然后计算每一对特征的信息增益得分，以上述模拟数据为例，数据集中共有5个特征，则会产生10组对特征，分别是：(f₁,f₂)，(f₁,f₃)，(f₁,f₄)，(f₁,f₅)，(f₂,f₃)，(f₂,f₄)，(f₂,f₅)，(f₃,f₄)，(f₃,f₅)，(f₄,f₅)，利用公式(1)分别计算每组对特征与类别的信息增益值，可得：IG(f₁；f₂；Y)＝0.04，IG(f₁；f₃；Y)＝0.06,IG(f₁；f₄；Y)＝-0.01，IG(f₁；f₅；Y)＝0.03，IG(f₂；f₃；Y)＝0.05，IG(f₂；f₄；Y)＝0.02，IG(f₂；f₅；Y)＝-0.03，IG(f₃；f₄；Y)＝0.02，IG(f₃；f₅；Y)＝-0.01，IG(f₄；f₅；Y)＝0.01，只保留IG>0的特征对进行下一步分析；

(2)利用信息增益得分构建协同作用网络：求得具有显著性信息增益得分的阈值θ₂＝0.02，则特征f₄与f₅的连边被删除，最后特征f₁与f₂，f₁与f₃，f₁与f₅，f₂与f₃，f₂与f₄，f₃与f₄构成连边；

(3)在所构建的网络上进行模块搜索；

(3.1)初始化：特征模块集合

特征模块子集

互信息集合

(3.2)计算每个特征与类别的互信息值，得到集合MIList＝{I(f₁；Y)＝

0.05,I(f₂；Y)＝0.02,I(f₃；Y)＝0.02,I(f₄；Y)＝0.01,I(f₅；Y)＝0.01}；

(3.3)取MIList中互信息值最大的点f₁作为初始种子节点，加入到当前特征模块子集S中，则此时S＝{f₁},计算AUC＝0.3；

(3.4)判断此时S是否有邻接点，由图2可知，此时S的邻接点不为空，候选节点集合为{f₂,f₃,f₅}。利用公式(3)计算出Score(f₂)＝0.06，Score(f₃)＝0.08，Score(f₅)＝0.04，可得特征f₃分数最高，尝试将其加入S中，并计算此时S＝{f₁,f₃}的AUC值，假设此时AUC＝0.5>0.3，AUC值增加则将其加入S中，继续重复步骤(3.4)的搜索过程。由图2可知当前S的邻接点不为空，候选集合为{f₂,f₄,f₅},利用公式(3)计算出Score(f₂)＝0.11，Score(f₄)＝0.02，Score(f₅)＝0.025，因此需将f₂加入S中，并计算此时S的AUC＝0.7>0.5，继续重复上述搜索过程。考虑此时S的邻接点仍不为空，候选集合C＝{f₄,f₅}，利用公式(3)计算此时,Score(f₄)＝0.0367，Score(f₅)＝0.02，因此需将f₄加入到S中，计算S＝{f₁,f₂,f₃,f₄}的AUC＝0.8>0.7，需继续执行模块搜索的过程，此时特征模块子集S的邻接点为f₅,计算Score(f₅)＝0.04，尝试将其加入S中，得到S的AUC＝0.75<0.8，AUC值减小，因此不再将其加入，当前特征模块子集S搜索结束，所包含的特征集合为{f₁,f₂,f₃,f₄}。继续向下执行步骤(3.5)；

(3.5)从集合MIList中删除S所包含的全部特征，在网络中删除S中节点及其连边，此时网络中仅剩特征f₅，MIList＝{f₅}；

(3.6)此时网络中节点数量不为空，因此仍需继续执行步骤(3.3)-(3.5)的模块搜索过程。由于只剩一个节点，因此只能选择f₅作为初始种子节点，计算S＝{f₅}的AUC＝0.72，由于f₅无邻居节点，因此该模块搜索结束。从MIList中删除f₅，并在网络中删除f₅及其连边，此时网络为空，因此所有模块搜索结束，那么最后所得特征模块集合MS＝{(f₁,f₂,f₃,f₄),(f₅)}；

(4)根据所得模块，确定最终的模块标志物。

对上述步骤(3)所得两个模块，其AUC值均大于等于0.7，因此两个模块均需保留，再对所得模块按AUC值降序排列，将排名靠前的模块顺次加入到已选集合中，通过内部交叉验证的方式，计算每加入一个新模块后已选集合的准确率。假定本次模拟数据中，第一个模块所得准确率最高，因此模块集合{f₁,f₂,f₃,f₄},即为最终确定的模块标志物。

表1 MSIG和FAST、FCBF之间的比较

Claims (1)

1.一种基于协同作用网络的生物组学数据分析方法，其特征在于，步骤如下：

将高通量测序技术所得的基因组学或转录组学数据中的基因序列视为特征，基因序列在不同样本上的表达值视为特征值；利用LC-MS或GC-MS技术得到代谢组学数据中的小分子代谢成分，将实验所得的小分子代谢成分视为特征，小分子代谢成分在不同样本上的含量视为特征值；

令F＝{f₁,f₂,…,f_m}为原始特征集合，包含m个特征；X＝{x₁,x₂,…,x_n}代表样本集合，n为样本总数；Y＝(y₁,y₂,…,y_n)为n个样本的类标向量，其中y_d∈C＝{c₁,c₂},d∈{1,2,…,n}，C为类标集合；

(1)计算特征f_i与f_j的信息增益得分IG(f_i；f_j；Y)，1≤i≠j≤m，公式如下：

IG(f_i；f_j；Y)＝I(f_i,f_j；Y)-I(f_i；Y)-I(f_j；Y) (1)

其中I(f_i,f_j；Y)表示特征f_i、f_j与类标向量Y的联合互信息大小，I(f_i；Y)、I(f_j；Y)为特征f_i、f_j与类标向量Y的互信息；IG(f_i；f_j；Y)反映f_i(f_j)与类别Y的依赖关系由于f_j(f_i)的加入所发生的变化；若IG(f_i；f_j；Y)>0，则表示特征f_i与f_j联合所提供关于Y的信息量比它们各自与Y的信息量之和大，进而认为f_i与f_j具有协同作用；

(2)获得具有显著性的信息增益得分；

当信息增益得分接近于零时，需确定具有显著性交互增益的特征对；对此保留步骤(1)大于零的信息增益值，并对其升序排列，得到IG⁺＝{p₁,p₂,…,p_t}，p₁≤p₂,…,≤p_t，

利用公式(2)计算组距：

再计算序列IG⁺中相邻两个变量的差值{p₂-p₁,p₃-p₂,…,p_t-p_t-1}，只保留结果大于θ₁的差值变量，并对其求平均值记为θ₂,将满足条件

的最小下标为i₀所对应的

值设定为边的阈值，如果

表明f_i与f_j具有较强的协同作用，IG值越大表明协同作用越强，因此在具有较强协同作用的特征对之间构建连边，此时网络中的节点代表特征，边的权值即为IG值；

(3)在所构建的协同作用网络上进行模块搜索；

由上述步骤(2)所得的协同作用网络，执行模块搜索的过程，过程如下：

(3.1)初始化：特征模块集合

特征模块子集

互信息集合

候选节点集合

(3.2)计算每个特征与类标向量的互信息，加入集合MIList；

(3.3)取MIList中互信息最大的点作为初始种子节点，加入到特征模块子集S中，并计算此时模块受试工作特征曲线ROC下的面积AUC值大小；

(3.4)定义候选节点集合CS＝Adj(S)-S，其中Adj(S)为当前特征模块子集S内所含节点在网络中的邻接点集，若当前特征模块子集S的候选集合为空集，则该模块搜索结束，继续向下执行步骤(3.5)；若候选集合不为空，则计算候选集合内的每个候选节点得分，公式为：

其中

分别表示点f_i、f_j在网络中的一近邻节点集，

则代表点f_i与f_j属于S相同的一近邻节点数量，|S|代表当前模块中节点个数，

与|S|相除表示f_i与f_j属于S相同的一近邻节点数量占当前模块节点总数的比例，若该比例越大且f_i与f_j构成的连边权值越大，则表明f_i与f_j连接越紧密，从而有效的衡量了候选节点与当前模块中节点连接的紧密程度；MI_i为特征f_i的互信息值，代表特征自身的分类性能，因此公式(3)将网络拓扑结构与特征自身分类性能相结合，综合衡量了每一个候选节点的重要性，得分越高表明该节点越重要，进而每次尝试将综合得分最高的节点加入到当前特征模块子集S中，检查模块的AUC值是否增加，若增加，则将该点加入到S中，继续重复步骤(3.4)的搜索过程，直至模块AUC值不再增加为止，此时该模块搜索结束，继续向下执行步骤(3.5)；

(3.5)S即作为本次搜索所得到的模块加入到特征模块集合MS中，从MIList中删除S所包含的特征，并从网络中删除S所包含的特征及其中的连边，以保证每个特征属于且仅属于一个特征模块，将S置空；

(3.6)重复上述步骤(3.3)-(3.5)的搜索过程，直至网络中节点数量为空，此时所有模块搜索结束；

(4)根据特征模块集合MS，确定最终的模块标志物；

对步骤(3)集合MS中所得模块，根据AUC值降序排列且只保留AUC值大于等于0.7的模块,利用内部交叉验证的方式，计算每加入一个新模块后已选集合的准确率，选择准确率最高的模块集合作为最终所确定的模块标志物。

Images