CN111145830A

CN111145830A - 基于网络传播的蛋白质功能预测方法

Info

Publication number: CN111145830A
Application number: CN201911369728.0A
Authority: CN
Inventors: 胡赛; 赵碧海; 熊慧军; 王皓鋆
Original assignee: Changsha University
Current assignee: Changsha University
Priority date: 2019-12-26
Filing date: 2019-12-26
Publication date: 2020-05-12

Abstract

本发明公开了一种基于网络传播的蛋白质功能预测方法，包括以下步骤：根据蛋白质相互作用网络拓扑结构建立无向加权网络，并将无向加权网络表征为无向图邻接矩阵M；根据蛋白质结构域信息和蛋白质复合物信息计算得到初始传播距离矩阵MI；根据无向图邻接矩阵M和初始传播距离矩阵MI，迭代地计算蛋白质间的网络传播距离以得到网络传播矩阵MP；根据网络传播矩阵MP逐步添加与测试蛋白质的邻接节点构成功能模块，多个功能模块构成候选功能列表；对候选功能列表中的候选功能评分、排序，选择排名靠前的若干项功能注释测试蛋白质。本发明改进了蛋白质功能预测方法研究中多源生物数据的融合方式，降低了相互作用网络中假阳性和假阴性对预测造成的负面影响，大大提高了蛋白质功能预测的准确率。

Description

基于网络传播的蛋白质功能预测方法

技术领域

本发明涉及蛋白质功能预测技术领域，具体涉及一种基于网络传播的蛋白质功能预测方法。

背景技术

蛋白质功能的准确注释是从分子水平理解生命的关键，对生物医学和制药产生巨大的影响。随着很多物种基因组测序工作的完成，序列数据与他们功能注释之间的差异日益增大。UniProt数据库中有功能注释的序列数据不足14％，序列数据还在成倍增长。早期预测蛋白质功能都是基于实验方法，如基因敲除，有针对性的突变与基因表达的抑制等。这些实验方法需要付出大量的实验成本和人力成本，也不适合注释蛋白质的重要部分，从而使得序列与功能之间的差距日益扩大。因此，提出高效的蛋白质功能注释的计算方法成为人们关注的热点问题之一。

目前蛋白质功能预测的计算方法大致可以分为基于序列结构的方法、基于网络的方法和多元生物信息融合的方法三种。

(1)基于序列结构的方法

早期的计算方法大都基于蛋白质氨基酸序列相似性或编码基因融合原理以及蛋白质系统发生谱研究。Henning等首先将GO注释应用于他们的GOblet的BLAST搜索结果，Groth只是简单地将序列匹配到有向无环图表示的GO。虽然没有结合相似分数，用户可以直观地了解有多少个序列匹配到共同的GO条目。Zehetner通过加权检索的GO条目列表，将这种匹配扩展到OntoBlast方法。Khan等提出的GOFigure方法和Martin等提出的GOtcha方法采用了相似的加权方式或排名方案，但是所有的分值从GO条目传播到GO DAG中的父条目层次。任何预测GO条目的最后分数根据各自的贡献被归一化。Stavros等使用SSP(序列相似性概要)将蛋白质与一组带注释的训练蛋白质进行比较，提出了两种新的LSDR方法，一种是基于GO结构的LSDR方法，另一种是基于术语语义相似性的LSDR方法。

(2)基于网络的方法

通过分析蛋白质间的相互作用，利用已有的蛋白质功能信息，根据未知蛋白质与已知蛋白质间的关联性来注释未知的蛋白质。Schwikowski等提出NC(neighbor counting)方法，将在邻居蛋白质中出现次数最多的功能预测为蛋白质的功能。Hishigaki等利用卡方统计量对待测蛋白质的近邻蛋白质中功能出现的次数进行评估，用期望出现的次数和实际出现的次数的差异来预测蛋白质是否具有某一功能。Chua等调查分析功能相似性与网络距离的联系，他们关注蛋白质的一级和二级邻居，设计了一种功能相似性计算方法，根据与目标蛋白质的距离为不同的蛋白质赋予不同的权值。

Vazquez等为每一个未注释的蛋白质赋予一个功能，使得连接相同功能蛋白质(未注释或已注释)的边数量最大。Karaoz等采用一种局部搜索程序轮流地搜索每一个节点，直至全部覆盖。节点的状态将根据大部分邻居的状态而改变，这个过程保证至少一半得到最优。Sun等也使用了类似的方法，只是他们使用了不同的目标函数和求解方法。Brun等提出基于蛋白质距离的PRODISTIN算法，采用BioNJ算法产生层次化的功能树，并进一步创建了在线的功能注释工具。类似的，Samanta等采用超几何分布计算P-value值来度量蛋白质间的距离，而预测功能时则仍采用多数表决的方法。Nabieva等提出一种相关的方法，将功能注释问题描述为一个多路切割问题。

(3)多元生物信息融合的方法

考虑到蛋白质相互作用网络的不完整性，研究者们通过在基于网络的方法基础上融合多元生物数据，提高预测的准确率。一类做法是将蛋白质相关的属性通过网络的形式表示和关联起来。网络有两种形式，一种是诸如蛋白质相互作用网络、代谢网络、基因调控网络等通过实验技术测量得到的网络；另一种是基于蛋白质的某些特性，用计算方法构建的相似性网络。如Liang等利用蛋白质共结构域特性构建蛋白质重叠网络(PON,ProteinOverlap Network)，并应用于功能注释。PON网络中，如果两个蛋白质共享相同的结构域，则两者之间存在相互作用。Peng等构建了三个网络：蛋白质作用网络、结构域共发生网络和功能相互关系网络，并通过在三个网络上采用不平衡随机游走算法实现功能预测。集成多个网络常用的方法是将多个网络中的蛋白质映射到同一个网络中去，构建一个超级的高可信网络。多个网络可以是采用不同实验手段检测到的蛋白质相互作用网络，也可以是基于基因共表达、蛋白质结构域相似性、序列相似性创建的网络。加权网络的方法也有多种。人们通常的做法是统计相互作用在不同实验结果中出现的次数来给网络加权，或者设计概率模型综合多个实验的结果或者不同生物信息创建的网络来构建加权的蛋白质相互作用网络。然后在加权网络上利用扩散的方法或者聚类的方法注释蛋白质功能。

还有些方法通过将蛋白质相互作用数据和其它生物数据相融合来计算蛋白质之间的相似性。Zhang等利用结构域信息计算两个蛋白质之间的功能相似性值，并将相似性最高的已知功能的蛋白质的全部功能对功能未知的蛋白质进行功能注释。Wang等提出DCS(Domain Combination Similarity)方法，改进Zhang的方法，通过加入二级邻居的结构域信息改进相似性计算公式，从而使得预测性能得到提高。Damiano等结合序列相似性、结构域体系结构搜索和蛋白质相互作用网络提出名为INGA的功能预测方法。Makrodimitris等结合蛋白质序列数据和GO条目的语义相似性，改善了蛋白质功能预测的性能。

蛋白质相互作用网络具有无幂律特性和小世界特性；高通量技术获得的相互作用数据包含大量的假阳性和假阴性。这些问题制约了基于计算方法的蛋白质功能预测方法准确率的提升。

发明内容

本发明提供了一种基于网络传播的蛋白质功能预测方法以解决现有技术中存在的蛋白质功能预测性能差的问题。

为解决上述技术问题，本发明所采取的技术方案如下。

一种基于网络传播的蛋白质功能预测方法，包括以下步骤：

根据蛋白质相互作用网络拓扑结构建立无向加权网络，并将无向加权网络表征为无向图邻接矩阵M；

根据蛋白质结构域信息和蛋白质复合物信息计算得到初始传播距离矩阵MI；

根据无向图邻接矩阵M和初始传播距离矩阵MI，迭代地计算蛋白质间的网络传播距离以得到网络传播矩阵MP；

根据网络传播矩阵MP逐步添加与测试蛋白质的邻接节点构成功能模块，多个功能模块构成候选功能列表；

对候选功能列表中的候选功能评分、排序，选择排名靠前的若干项功能注释测试蛋白质。

优选地，建立无向加权网络时，通过边聚集系数计算蛋白质间的权值W(P_i，P_j)，

其中，

为蛋白质P_i邻居节点的集合，

为蛋白质P_j的邻居节点的集合，

为蛋白质P_i和蛋白质P_j的共同邻居节点的集合，蛋白质P_i和蛋白质P_j为蛋白质相互作用网络中的任意两个蛋白质。

优选地，

其中，Sim_D(i，j)为蛋白质P_i和蛋白质P_j之间基于结构域的功能相似性，Sim_C(i，j)为蛋白质P_i和蛋白质P_j之间基于复合物的功能相似性，n为蛋白质的数量，i∈[1，n]，j∈[1，n]。

优选地，

其中，DS_i为蛋白质P_i的不同类型结构域构成的集合，DS_j为蛋白质P_j的不同类型结构域构成的集合，DS_i∩DS_j为蛋白质P_i和蛋白质P_j的相同结构域类型构成的集合。

优选地，

其中，CS_i为包含蛋白质P_i的复合物的集合，CS_j为包含蛋白质P_j的复合物的集合，CS_i∩CS_j为同时包含蛋白质P_i和蛋白质P_j的复合物构成的集合。

优选地，网络传播矩阵MP的构建步骤如下：

S301、从蛋白质列表中选择一个蛋白质V_i；

S302、初始化迭代步伐t＝0,Dis₀＝[1/n，1/n，…1/n]，其中，n表示蛋白质数量；

S303、迭代地计算其他蛋白质与蛋白质V_i的网络传播距离，计算公式为：

Dis_t+1＝(1-α)·M·Dis_t+α·MI_i，

其中，参数α为重启概率，且0≤a≤1，MI_i为初始传播距离矩阵MI第i行元素构成的向量；

S304、如果||Dis_t+1-Dis_t||≥ε，使t＝t+1，并返回S302继续进行迭代，否则，迭代终止；

S305、返回S301，从蛋白质列表中选择下一个蛋白质继续计算网络传播距离，直至蛋白质列表中所有蛋白质节点遍历完毕得到网络传播矩阵MP，

优选地，通过聚类形成功能模块的步骤如下：

S401、将测试蛋白质v的邻居节点以与测试蛋白质v关联从强到弱的顺序依次加入队列；

S402、蛋白质v作为种子加入功能模块G；

S403、从队首取节点u，并计算节点u对功能模块G的适应度f_Gu，如果f_Gu>0，将节点u加入功能模块G；

S404、从队列中取下一个节点，返回至S503，直至队列为空。

优选地，节点u对功能模块G的适应度计算公式为：

其中，G+{u}和G-{u}分别表示功能模块包含节点u和不包含节点u，适应度函数f_G定义如下：

其中，

是功能模块G内所有边的权值之和，

是功能模块G中节点连向功能模块G外的边权值之和，β是控制功能模块规模大小的调节参数。

优选地，候选功能f_j的得分计算方法如下所示：

其中，F＝{f₁,f₂,…,f_m}是聚类得到的功能模块中所有蛋白质的功能集合；

MP(v,G_i)表示蛋白质v与蛋白质G_i之间的网络传播距离；若蛋白质G_i包含功能f_j，则t_ij＝1，否则t_ij＝0。

由于采用了以上技术方案，本发明所取得技术进步如下。

由于蛋白质相互作用数据是不完整的，且有些蛋白质之间的相互作用没有直接体现，但是这些蛋白质可能共享相同的结构域信息或复合物信息，因此，通过结合蛋白质相互作用网络拓扑结构、蛋白质结构域信息及蛋白质复合物信息计算蛋白质之间的网络传播距离得到的网络传播矩阵，能够挖掘隐含的、未知的蛋白质相互作用和关联，同时，蛋白质相互作用数据具有小世界特性，通过网络传播，能够让网络更加稠密，从而也能挖掘隐含的、未知的蛋白质相互作用和关联。基于这三种信息挖掘的蛋白质相互作用不仅更加全面和可靠，而且消除了一些不真实的蛋白质相互作用，从而使得蛋白质的功能预测准确率更高。

附图说明

图1是本发明的流程图；

图2为本发明的基于网络传播的蛋白质功能预测方法、FPM、Zhang、DCS、NC、PON在DIP数据集分别对酵母蛋白质的功能进行预测的结果；

图3为本发明的基于网络传播的蛋白质功能预测方法、FPM、Zhang、DCS、NC、PON在Krogan数据集分别对酵母蛋白质的功能进行预测的结果。

具体实施方式

下面将结合附图和具体实施例对本发明进行进一步详细说明。

本发明的基于网络传播的蛋白质功能预测方法可以用于酵母蛋白质的功能预测。

其中，蛋白质相互作用网络、蛋白质结构域信息及蛋白质复合物信息均源于公共数据库。

如图1所示，本发明的基于扩散距离的关键蛋白质识别方法，包括以下步骤：

S1：根据蛋白质相互作用网络拓扑结构建立无向加权网络，并将无向加权网络表征为无向图构建邻接矩阵M。

边聚集系数用于刻画网络中某个节点与其邻居之间的亲疏程度，也是复杂网络中最重要的拓扑特征之一。建立无向加权网络时，通过边聚集系数计算蛋白质间的权值W(P_i，P_j)，

其中，

为蛋白质P_i的邻居节点的集合，

为蛋白质P_j的邻居节点的集合，N_Pi∩N_Pj为蛋白质P_i和蛋白质P_j共同的邻居节点的集合，蛋白质P_i和蛋白质P_j为蛋白质相互作用网络中的任意两个蛋白质。

S2：根据蛋白质结构域信息和蛋白质复合物信息计算得到初始传播距离矩阵MI。

初始传播距离矩阵MI的计算公式如下：

其中，Sim_D(i，j)为蛋白质P_i和蛋白质P_j之间基于结构域的功能相似性(或者说，Sim_D(i，j)为基于结构域信息得到的蛋白质P_i和蛋白质P_j具有相同或相似功能的可能性)，Sim_C(i，j)为蛋白质P_i和蛋白质P_j之间基于复合物的功能相似性(或者说，Sim_C(i，j)为基于复合物信息得到的蛋白质P_i和蛋白质P_j具有相同或相似功能的可能性)，n为蛋白质的数量，i∈[1，n]，j∈[1，n]。

蛋白质P_i和蛋白质P_j之间基于结构域的相似性Sim_D(i，j)的计算公式如下：

其中，DS_i为蛋白质P_i的不同类型结构域构成的集合，DS_j为蛋白质P_j的不同类型结构域构成的集合，DSi∩DSj为蛋白质P_i和蛋白质P_j的相同结构域类型构成的集合。

蛋白质P_i和蛋白质P之间基于复合物的相似性Sim_C(i，j)计算公式如下：

其中，CS_i为包含蛋白质P_i的复合物的集合，CS_j为蛋白质P_j的所在复合物的集合，CS_i∩CS_j为同时包含蛋白质P_i和蛋白质P_j的复合物构成的集合。

S3：根据无向图邻接矩阵M和初始传播距离矩阵MI，迭代地计算蛋白质间的网络传播距离以得到网络传播矩阵MP。

网络传播矩阵MP具体的构建步骤如下：

S301、从蛋白质列表中选择一个蛋白质V_i。

S302、初始化迭代步伐t＝0,Dis₀＝[1/n，1/n，…1/n]，其中，n表示蛋白质数量。

Dis_t+1＝(1-α)·M·Dis_t+α·MI_i，

其中，参数α为重启概率，且0≤a≤1，1-α表示移动到相邻节点的概率，MI_i为初始传播距离矩阵MI第i行元素构成的向量。

S304、如果||Dis_t+1-Dis_t||≥ε，使t＝t+1，并返回S302继续进行迭代，否则，迭代终止。

其中，ε是一个控制迭代终止的参数，该参数值越大，迭代收敛越快。本发明中，默认地设置为10^-5。

S4：根据网络传播矩阵MP逐步添加与测试蛋白质的邻接节点构成功能模块，多个功能模块构成候选功能列表。

换而言之，每个蛋白质的邻接节点构成邻接节点集合(即功能模块)，多个蛋白质的邻接节点集合(即功能模块)构成候选功能列表。

通过聚类形成功能模块的具体步骤如下：

S402、蛋白质v作为种子加入功能模块G；

S403、从队首取节点u，并计算节点u对功能模块G的适应度

如果

将节点u加入功能模块G。

节点u对功能模块G的适应度计算公式如下：

G+{u}和G+{u}分别表示功能模块包含节点u和不包含节点u。适应度函数f_G定义如下：

其中，

是功能模块G内所有边的权值之和，

是功能模块G中节点连向功能模块G外的边权值之和，β是控制功能模块规模大小的调节参数。通过实验发现，β设置为时，预测准确性最佳。

S404、从队列中取下一个节点，返回至S503，直至队列为空。

S5：对候选功能列表中的候选功能评分、排序，选择排名靠前的若干项功能注释测试蛋白质。

候选功能f_j的得分计算方法如下所示：

其中，F＝{f₁,f₂,…,f_m}是聚类得到的功能模块中所有蛋白质的功能集合；MP(v,G_i)表示蛋白质v与蛋白质G_i间的网络传播距离；若蛋白质G_i包含功能f_j，则t_ij＝1，否则t_ij＝0。

需要说明的是，在本发明中，步骤S1和步骤S2可以依次进行，也可以步骤S2和步骤S1依次进行，还可以步骤S1和步骤S2同时进行。

在本实施例中，优选地，对蛋白质相互作用网络拓扑结构、蛋白质结构域信息及蛋白质复合物信息进行数据预处理。预处理主要时进行数据格式转换，确保蛋白质相互作用网络、蛋白质结构域信息及蛋白质复合物信息等采用相同的格式。

本发明的基于网络传播的蛋白质功能预测方法、FPM、Zhang、DCS、NC、PON在DIP和Krogan两个数据集分别对酵母蛋白质的功能进行预测，预测结果如图2和图3所示，其中，图2为本发明的基于网络传播的蛋白质功能预测方法、FPM、Zhang、DCS、NC、PON在DIP数据集分别对酵母蛋白质的功能进行预测的结果；图3为本发明的基于网络传播的蛋白质功能预测方法、FPM、Zhang、DCS、NC、PON在Krogan数据集分别对酵母蛋白质的功能进行预测的结果。图2和图3中的NGF为本发明的基于网络传播的蛋白质功能预测方法的预测结果。其中，蛋白质功能预测从以下三个方面体现：

特异性(Specificity)衡量的是在预测的功能中有多少比例的功能与已知的功能之间能够匹配；

敏感性(Sensitivity)衡量的是已知的蛋白质功能有多少能够与被预测的功能所匹配；

F-measure则能较好地反映算法的综合性能，它是Specificity和Sensitivity的调和平均值。

从图2和图3中可以看出，本发明基于网络传播的蛋白质功能预测方法相比于其他方法具有更高的准确度。

其中，蛋白质相互作用网络来源于酿酒酵母(面包酵母)DIP数据库，其由4023个蛋白质和22570条边组成；已经通过基因敲除实验被很好地特征化。

用于验证算法性能的蛋白质功能注释数据是从GO官方网站下载的最新版本。为了避免太特殊或者太一般化，仅仅使用那些至少注释了10个或者最多注释了200个蛋白质的GO Term来进行实验验证，处理后的GO Term数量为267个。此外，GO注释数据中，利用Uniprot网站将蛋白质格式从UniProtKB转换为Ensemble Genomes Protein，以便与蛋白质相互作用网络中蛋白质的格式匹配。

蛋白质结构域数据从Pfam数据库下载得到，包含1107个不同的结构域，涉及蛋白质相互作用网络中的3056个蛋白质。

蛋白质复合物数据采用CYC2008数据集，CYC2008包含408个通过生物方法预测得到的复合物，并被作为标准的已知复合物集合，广泛应用于蛋白质复合物预测方法评价。