CN110610742B - 一种基于蛋白质互作网络的功能模块检测方法 - Google Patents

Abstract

本发明公开一种基于蛋白质互作网络的功能模块检测方法，在蛋白质互作网络中，使用LAR机制表示全部的解集，加快了算法的收敛速度。然后使用多目标粒子群优化算法在解空间搜索最优解。本发明基于粒子群算法使用帕雷托支配概念解决了多目标任务，并在在粒子群（PSO）算法框架基础上加入了遗传操作提高检测性能。本发明方法不仅显著优于最先进的方法，而且还能更有效、准确地检测功能模块。

Description

一种基于蛋白质互作网络的功能模块检测方法

技术领域

本发明涉及生物信息领域，尤其涉及一种基于蛋白质互作网络的功能模块检测方法。

背景技术

人类基因组计划的顺利进行标志着后基因组时代的到来。在后基因组时代，系统地分析和全面理解蛋内质之间通过相互作用完成生命活动的规律成为最热点的研究问题之一。蛋白质是构成一切细胞和组织结构必不可少的成分，是生命活动最重要的物质基础。但是，细胞中的每个蛋白质并不是孤立存在的，而是与其他蛋白质一起进行相互作用来行使其功能的蛋白质相互作用(Protein-Protein Interaction，PPI)在生命中起核心作用，如细胞的新陈代谢、信号传导、DNA复制、细胞周期控制等过程的基础，在许多疾病的发生和治疗过程中也起着至关重要的作用。一个生物体内所有蛋白质间的相互作用组成的生物分子关系网络称为蛋白质相互作用网络(PPI网络)。在该网络水平上对蛋白质相互作用的研究有助于对生命活动的全面和深入认识，是解决生命科学中重要难题的一条有效途径。其中，PPI网络中功能模块检测就是理解生命现象的一个突破口，是当前生命科学的研究前沿之一。

发明内容

本发明的目的在于提供一种基于蛋白质互作网络的功能模块检测方法。

本发明采用的技术方案是：

一种基于蛋白质互作网络的功能模块检测方法，其包括以下步骤：

S1.算法初始化：解集表示和适应度计算。使用LAR机制表示出全部的解集，后每个解对应一个粒子p，对于每个粒子分别计算他们的目标函数KKM、RC，使用非支配排序机制，得到非支配解集PS^*，从PS^*中选择NMI值最大的粒子作为全局最优解p’_gbest；

S2.移动局部最优策略：为了使每个粒子移动到局部最优位置，p_i和p_i,best之间做一次两点交叉操作，得到两个结果child₁、child₂，p_i指的是第i个粒子，p_i,best指的是第i个粒子的最优位置；分别计算对应的KKM、RC的值，若child₁可以支配child₂，则存储粒子的暂时位置temp_par＝child₁，若child₂可以支配child₁，则存储粒子的暂时位置temp_par＝child₂，若两者之间不存在支配关系，则选择NMI值最大的结果child₁或child₂存于temp_par；

S3.移动全局最优策略：为了使粒子移动到全局最优位置，粒子暂存位置temp_par和全局最优解p_gbest之间做一次两点交叉操作，得到两个结果child₃、child₄，分别计算对应的KKM、RC的值，若child₃可以支配child₄，则存储粒子的暂时位置temp_par＝child₃，若child₄可以支配child₃，则存储粒子的暂时位置temp_par＝child₄，若两者之间不存在支配关系，则选择NMI值最大的结果child₃或child₄存于temp_par；

S4.为了能够在整个解空间进行搜索最优解，粒子暂存位置temp_par做一次突变操作，得到结果p_i'，若p_i'可以支配p_i,best，则替换粒子最优位置p_i,best＝p_i'，否则不更新粒子最优位置；

S5.以粒子总数作为最大迭代次数，判断迭代次数是否超过最大迭代次数；当达到最大迭代次数时，则执行S6；否则，迭代执行S2至S5；

S6.等到所有粒子到达最优位置时分别计算每个粒子对应的KKM、RC、NMI值，利用非支配排序机制得到非支配解集，选择非支配解集中的粒子的NMI值最大的作为全局最优解p_gbest。

进一步地，步骤S1的具体步骤为：

S101.输入蛋白质互作网络数据邻接矩阵A，并使用LAR机制表示所有的解集，LAR机制将每个解看作长度为N的数组，其中每个位置存储的值为节点的任意一个邻居的对应的标号。因此每个位置的取值范围为1到N。第i解对应第i个粒子p_i，初始化粒子最优位置p_i,best＝p_i；

S102.每个解对应一个粒子p，对于每个粒子分别计算他们的目标函数

其中n为蛋白质个数，m为功能模块个数，V_j为第j个功能模块，|V_j|为V_j中蛋白质个数，为除V_j以外的所有功能模块，/>为/>中蛋白质个数，/> A_kl为A的第k行第l列元素；

S103.使用非支配排序机制，得到非支配解集PS^*；

S104.计算非支配解集中粒子对应的NMI值，选择NMI值最大的作为全局最优解p’_gbest，NMI的计算公式如下：

其中，PPI真实的划分为A＝{A₁,A₂,...,A_R},步骤S101得到的划分为B＝{B₁,B₂,...,B_D}，R和D分别为划分A和B中功能模块的个数，C为混淆矩阵C_st为同时出现在A_s和B_t中的节点个数，N为节点数。

本发明采用以上技术方案，使用的多目标粒子群优化算法，同时算法基于粒子群(PSO)算法使用帕雷托支配概念解决了多目标任务，同时在该算法中最小化两个目标函数Kernel k-means(KKM)和ratio cut(RC)。另外，为了提高检测性能，在粒子群(PSO)算法框架基础上加入了遗传操作，如交叉和突变操作。最后在真实PPI网络的实验结果表明所提出的方法不仅显著优于最先进的方法，而且还能更有效、准确地检测功能模块。

附图说明

以下结合附图和具体实施方式对本发明做进一步详细说明；

图1为本发明一种基于蛋白质互作网络的功能模块检测方法的流程示意图；

图2为本发明所用的LAR机制表示方法示意图；

图3为本发明搜索策略中的两点交叉操作示意图。

具体实施方式

以下结合附图和具体实施方式对本发明做进一步详细说明：

通常，蛋白质互作网络被建模为无向图G＝(V,E)，其中V和E分别代表蛋白质节点集合和蛋白质节点之间的相互作用集合。一个功能模块可以看作是其中的稠密子图。因此，功能模块检测方法可以看做图聚类问题。

为了能够更高效、更准确的检测蛋白质网络中的功能模块，本发明提出的一种基于蛋白质互作网络的功能模块检测方法。该方法使用的是一个新的多目标粒子群优化算法，算法主要有两大步骤：

首先为算法的初始化，初始化包括解集表示和适应度计算。本算法中解集表示采用的是基于位置的邻接表示(LAR)作为其表示方案。LAR机制将每个解看作长度为N的数组，其中每个位置存储的值为节点的任意一个邻居的对应的标号。因此每个位置的取值范围为1到N。每个解决方案可以被解码为一个图划分。其中得到的稠密子图便为该解决方案下对应的功能模块。使用这种表示形式有几个优点。首先，通过使用该方案，我们的算法不需要预先知道功能模块的数量。这是因为在解码过程中才能确定功能模块的数量。此外，解的解码过程可以在线性时间内完成。此外，使用LAR表示机制，可以方便地对解执行交叉操作。

在图2(b)中，给出了一个具有13个节点的PPI网络的LAR方案示例图。图2(a)绘制了PPI网络的图结构。一个任意由LAR格式表示的粒子(解)及其解码得到图，如图2(b)所示。可以看出，每个位置取值范围为1-13，对于每个位置存储的值为节点的任意一个邻居的对应的标号，如图2(b)，位置ID为6的值为7，表示在对应的图中有一个边从节点6到节点7。换句话说，这两个节点被放在同一个子图中，如图2(b)中的虚线圈所示。

在使用LAR机制表示全部解集后，每个解对应一个粒子p，对于每个粒子分别计算他们的目标函数对于每个粒子分别计算他们的目标函数KKM、RC。KKM为衡量功能模块内部连接密度的指标，RC衡量功能模块之间的连接密度的指标，最小化KKM、RC，能够确保检测得到的功能模块内部连接紧密，外部连接稀疏。使得检测到的功能模块更加准确。计算完之后使用非支配排序机制，得到非支配解集PS^*，选择其中NMI值最大的作为全局最优解p_gbest，NMI为测量真实蛋白质功能模块结构与算法得到功能模块结构的相似性指标。

其次为搜索策略。粒子群优化算法的搜索策略是将每个粒子分别移动到局部最优和全局最优。为了将每个粒子移动到最佳位置，采用了交叉和突变等遗传操作。首先为交叉操作，图3是一个两点交叉的例子。图3(a)和(b)分别给出了作为双亲的两个任意解P₁和P₂及其对应的图结构。图3(c)中确定了两个随机点i＝7,j＝10。为了产生第一个子解，从粒子开始到交叉点i的二进制串从父解P₁复制，从交叉点i到交叉点j的部分从父解P₂复制，其余部分从父解P₁复制。要生成第二个子解，此操作按相反的顺序执行。图3(d)展示了第一个子解及其相应的图划分。为了使每个粒子移动到局部最优位置，p_i和p_i,best之间做一次两点交叉操作，得到两个结果child₁、child₂，分别计算对应的KKM、RC的值，若child₁可以支配child₂，则存储粒子的暂时位置temp_par＝child₁，若child₂可以支配child₁，则存储粒子的暂时位置temp_par＝child₂，若两者之间不存在支配关系，则选择NMI值最大的结果child₁或child₂存于temp_par。

为了使粒子移动到全局最优位置，粒子暂存位置temp_par和全局最优解p_gbest之间做一次两点交叉操作，得到两个结果child₃、child₄，分别计算对应的KKM、RC的值，若child₃可以支配child₄，则存储粒子的暂时位置temp_par＝child₃，若child₄可以支配child₃，则存储粒子的暂时位置temp_par＝child₄，若两者之间不存在支配关系，则选择NMI值最大的结果child₃或child₄存于temp_par；

最后，为了对整个搜索空间进行搜索，对所有粒子进行基于一点邻居的突变，对于每个粒子，随机选取一个位置i，然后用其邻居节点的一个可能值替换现有值。这个限制是为了保证在解空间中产生可能的解。对所有粒子暂存位置temp_par做一次突变操作，得到结果p_i'，若p_i'可以支配p_i,best，则替换粒子最优位置p_i,best＝p_i'，否则不更新粒子最优位置。

以粒子总数作为最大迭代次数，判断迭代次数是否超过最大迭代次数。当达到最大迭代次数时，则等到所有粒子到达最优位置时分别计算其对应的KKM、RC、NMI值，利用非支配排序机制得到非支配解集，选择其中NMI值最大的作为全局最优解p_gbest，输出全局最优解p_gbest；否则，继续迭代。

在本发明中使用了LAR机制表示全部的解集，加快了算法的收敛速度。另外本算法基于粒子群(PSO)算法使用帕雷托支配概念解决了多目标任务，同时在该算法中最小化两个目标函数KKM和RC。另外，为了提高检测性能，我们在粒子群(PSO)算法框架基础上加入了遗传操作，如交叉和突变操作。最后在真实PPI网络的实验结果表明所提出的方法不仅显著优于最先进的方法，而且还能更有效、准确地检测功能模块。

如图1所示的流程图：本发明公开了一种基于蛋白质互作网络的功能模块检测方法，其包括以下步骤：

S1算法初始化：解集表示和适应度计算。使用LAR机制表示出全部的解集，后每个解对应一个粒子p，对于每个粒子分别计算他们的目标函数KKM、RC，使用非支配排序机制，得到非支配解集PS^*，从PS^*中选择NMI值最大的粒子作为全局最优解p_gbest；

S2移动局部最优策略：为了使每个粒子移动到局部最优位置，p_i和p_i,best之间做一次两点交叉操作，得到两个结果child₁、child₂，分别计算对应的KKM、RC的值，若child₁可以支配child₂，则存储粒子的暂时位置temp_par＝child₁，若child₂可以支配child₁，则存储粒子的暂时位置temp_par＝child₂，若两者之间不存在支配关系，则选择NMI值最大的结果child₁或child₂存于temp_par；

S3移动全局最优策略：为了使粒子移动到全局最优位置，粒子暂存位置temp_par和全局最优解p_gbest之间做一次两点交叉操作，得到两个结果child₃、child₄，分别计算对应的KKM、RC的值，若child₃可以支配child₄，则存储粒子的暂时位置temp_par＝child₃，若child₄可以支配child₃，则存储粒子的暂时位置temp_par＝child₄，若两者之间不存在支配关系，则选择NMI值最大的结果child₃或child₄存于temp_par；

S4为了能够在整个解空间进行搜索最优解，对所有粒子暂存位置temp_par做一次突变操作，得到结果p_i'，若p_i'可以支配p_i,best，则替换粒子最优位置p_i,best＝p_i'，否则不更新粒子最优位置；

S5.以粒子总数作为最大迭代次数，判断迭代次数是否超过最大迭代次数；当达到最大迭代次数时，则执行S6；否则，迭代执行S2至S5；

S6.等到所有粒子到达最优位置时分别计算每个粒子对应的KKM、RC、NMI值，利用非支配排序机制得到非支配解集，选择非支配解集中的粒子的NMI值最大的作为全局最优解p_gbest。

进一步地，步骤S1的具体步骤为：

S101.输入蛋白质互作网络数据邻接矩阵A，并使用LAR机制表示所有的解集，LAR机制将每个解看作长度为N的数组，其中每个位置存储的值为节点的任意一个邻居的对应的标号。因此每个位置的取值范围为1到N。第i解对应第i个粒子p_i，初始化粒子最优位置p_i,best＝p_i；

S102.每个解对应一个粒子p，对于每个粒子分别计算他们的目标函数

其中n为蛋白质个数，m为功能模块个数，V_j为第j个功能模块，|V_j|为V_j中蛋白质个数，为除V_j以外的所有功能模块，/>为/>中蛋白质个数，/> A_kl为A的第k行第l列元素；

S103.使用非支配排序机制，得到非支配解集PS^*；

S104.计算非支配解集中粒子对应NMI的值，其计算公式如下：

其中，PPI真实的划分为A＝{A₁,A₂,...,A_R},步骤S101得到的划分为B＝{B₁,B₂,...,B_D}，R和D分别为划分A和B中功能模块的个数，C为混淆矩阵C_st为同时出现在A_s和B_t中的节点个数，N为节点数。

本发明采用以上技术方案，在蛋白质互作网络中，使用LAR机制表示全部的解集，加快了算法的收敛速度。然后使用多目标粒子群优化算法在解空间搜索最优解，另外本算法基于粒子群(PSO)算法使用帕雷托支配概念解决了多目标任务，同时在该算法中最小化两个目标函数KKM和RC。另外，为了提高检测性能，我们在粒子群(PSO)算法框架基础上加入了遗传操作，如交叉和突变操作。最后在真实PPI网络的实验结果表明所提出的方法不仅显著优于最先进的方法，而且还能更有效、准确地检测功能模块。

本发明的有益效果在于：

1、初始化数据时使用LAR机制表示了对应的解集，提高了我们所提方法的收敛速度；

2、利用帕雷托概念解决了多目标任务，同时在该算法中最小化两个目标函数KKM和RC；

3、在粒子群算法框架基础上加入了遗传操作，如交叉和突变操作。大大地提高了检测性能，同时能够更准确地检测功能模块。

Claims (4)

1.一种基于蛋白质互作网络的功能模块检测方法，其特征在于：其包括以下步骤：

S1.算法初始化：使用LAR机制表示出全部的解集，每个解对应一个粒子p，分别计算每个粒子的目标函数KKM、RC，使用非支配排序机制得到非支配解集PS^*，从PS^*中选择NMI值最大的粒子作为初始全局最优解并存于暂存全局最优解p’_gbest；步骤S1的具体步骤为：

S101.输入蛋白质互作网络数据邻接矩阵A，并使用LAR机制表示所有的解集，LAR机制将每个解看作长度为N的数组，其中每个位置存储的值为节点的任意一个邻居的对应的标号，每个位置的取值范围为1到N；

S102.每个解对应一个粒子p，分别计算每个粒子的目标函数

其中n为蛋白质个数，m为功能模块个数，V_j为第j个功能模块，|V_j|为V_j中蛋白质个数，为除V_j以外的所有功能模块，/>为/>中蛋白质个数，/> A为蛋白质互作网络数据邻接矩阵，A_kl为A的第k行第l列元素；

S103.使用非支配排序机制得到非支配解集PS^*；

S104.计算非支配解集中粒子对应的NMI值，选择NMI值最大的粒子作为初始全局最优解并存于暂存全局最优解p’_gbest，NMI的计算公式如下：

其中，蛋白质互作网络真实的划分为A＝{A₁,A₂,...,A_R},步骤S101得到的划分为B＝{B₁,B₂,…,B_D}，R和D分别为划分A和B中功能模块的个数，C为混淆矩阵，C_st为同时出现在A_s和B_t中的节点个数，N为节点数；

S2.基于移动局部最优策略使粒子移动到局部最优位置：选择一个未寻优的粒子p_i，将该粒子p_i和该粒子的最优位置p_i,best之间做一次两点交叉操作，比较得到粒子的局部最优位置存储于暂存位置temp_par；

S3.基于移动全局最优策略使粒子移动到全局最优位置：将粒子暂存位置temp_par的数值与暂存全局最优解p’_gbest之间做一次两点交叉操作，比较得到粒子的全局最优位置存储于暂存位置temp_par；

S4.将粒子暂存位置temp_par做一次突变操作得到结果p_i'；当p_i'支配p_i,best时，则替换粒子最优位置p_i,best＝p_i'；否则，不更新粒子最优位置；

S5.以粒子总数作为最大迭代次数，判断迭代次数是否超过最大迭代次数；当达到最大迭代次数时，则执行S6；否则，迭代执行S2至S5；

S6.等到所有粒子到达最优位置时分别计算每个粒子对应的KKM、RC、NMI值，利用非支配排序机制得到非支配解集，选择非支配解集中NMI值最大的粒子作为最终全局最优解p_gbest。

2.根据权利要求1所述的一种基于蛋白质互作网络的功能模块检测方法，其特征在于：步骤S101中第i解对应第i个粒子p_i，初始化粒子最优位置p_i,best＝p_i。

3.根据权利要求1所述的一种基于蛋白质互作网络的功能模块检测方法，其特征在于：步骤S2的具体方法为：

将p_i和p_i,best之间做一次两点交叉操作，得到两个结果child₁、child₂，分别计算child₁、child₂对应的KKM、RC的值；

当child₁支配child₂时，则存储粒子的暂时位置temp_par＝child₁；

当child₂支配child₁时，则存储粒子的暂时位置temp_par＝child₂，

当两者之间不存在支配关系时，则选择child₁和child₂中NMI值最大的结果存于temp_par。

4.根据权利要求1所述的一种基于蛋白质互作网络的功能模块检测方法，其特征在于：步骤S3的具体方法为：

将粒子暂存位置temp_par和暂存全局最优解p’_gbest之间做一次两点交叉操作，得到两个结果child₃、child₄，分别计算child₃、child₄对应的KKM、RC的值，

当child₃支配child₄时，则存储粒子的暂时位置temp_par＝child₃；

当child₄支配child₃时，则存储粒子的暂时位置temp_par＝child₄；

当两者之间不存在支配关系时，则选择child₁和child₂中NMI值最大的结果存于temp_par。