基于复杂网络模型的基因调控网络的计算模拟

第37卷　第1期计算机科学Vol.37No.1基于复杂网络模型的基因调控网络的计算模拟

张律文1,2　谢　江1,2　陈建娇1　张　武1

(上海大学计算机工程与科学学院　上海200072)1　(上海大学系统生物技术研究所　上海200444)2

摘　要　随着大通量基因芯片数据的产生,基因调控机制的网络化研究需求日趋迫切。提出了基于复杂网络理论的基因调控网络的模拟方法,构建了基因调控网络模拟器GN2Simulator。通过分析真实基因调控网络的拓扑特性,提出了对应的矩阵模型,并充分考虑了网络的生物学鲁棒性和动力学稳定性,给出了人工基因网络的生成过程和计算模拟方法。计算实验表明,GN2Simulator能高效地模拟出与真实基因调控网络高度相似的大规模人工网络,并可为不同算法提供无偏验证的多样化人工模拟数据。关键词　复杂网络,无标度特性,动力学稳定性,模拟,基因调控网络

ArtificialGeneRegulatoryNetworksConstructionBasedonComplexNetwork

ZHANGLu2wen1,2　XIEJiang1,2　CHENJian2jiao1　ZHANGWu1

(SchoolofComputerEngineeringandScience,ShanghaiUniversity,Shanghai200072,China)1

(InstituteofSystemsBiology,ShanghaiUniversity,Shanghai200444,China)2

Abstract　Alongwiththeexplosionofhigh2throughputgenomicinformationproducedfDNAmicroarrays,modelingthemechanismofgeneregulationsandconstructinggenetask.Wepro2posedanovelmethod,GN2Simulator,tosimulategeneregulatoryontheory.TherealgenenetworksweremodeledexplicitlyinmatrixesAccordingtotherobustbiologicalmechanism,weconstructedartificialgenetothedynamicstability.WeconsideredthefeaturesofrealgenenetworkThecomputingexperimentsillustratethatGN2Simulatorcansimulatewithhighconfidencelevel.Moreover,GN2Simulatorgenera2tesvariousdifferentalgorithmsperformanceandprovidesreasonableestimationonthem.

Keywords　Complexnetwork,Scale2freetopology,Dynamicstability,Simulation,Generegulatorynetwork　

究成果之后,复杂网络被广泛地应用于生物网络,包括基因调控网络的构建和模拟等各方面,使得基因网络的研究取得了极大进展。

用来模拟真实基因调控的复杂网络模型经历了几个阶段的发展。最初,简单的随机网络模型[2]被用来描述基因网络的结构,但在1998年,Watts和Strogats在结合规则网络和随机网络特点的基础上,建立了小世界(Small2world)网络模型[3],接着1999年,Barabase和Albert发现了真实网络的无标度(Scale2free)性质[4]。在对大规模基因网络进行数据采集和统计分析后,无标度被认为是最接近真实基因网络的一种拓扑结构[5],基因相互作用网络的结构特性可能与其他复杂系统网络(比如Internet网)在很大程度上是一致的。

在构建人工基因网络的过程中,仅仅关注网络的大体拓扑结构是远远不够的。为了尽可能地模拟出接近真实基因网络的人工网络,构建了基因调控网络模拟器Genenetwork

simulator(GN2Simulator),从以下几个方面系统地考虑了网

1　引言

复杂生命现象是大量基因互相调控、协同作用的结果。随着基因组测序的完成、基因组学的快速发展,生物学积累了大量的基因调控实验数据,如何挖掘出海量数据所蕴藏的生命现象和生物规律已成为生命科学的研究热点。如今,复杂生物系统的研究开始从对细胞内个别基因功能的局部性描述发展到对复杂基因调控网络(GeneRegulatoryNetwork)的定量刻画。由于基因与其他生物小分子倾向于成组地通过网状的相互作用而影响生物系统的功能,因此对基因功能的研究必须分析其调控网络。

关于基因调控机制的系统性研究,最行之有效的方法就是结合计算机科学和数学等理论把基因之间的复杂相互作用进行整合,并简化成网络,然后建立理论模型,从而研究基因调控网络的结构、网络模块及动力学性质。自2000年Nature上发表了利用复杂网络就理论研究生物网络拓扑特性[1]的研

到稿日期:2009202220　返修日期:2009205208　　本文受上海市重点学科建设项目(项目编号:J50103),上海大学系统生物研究基金(200825)资助。

张律文　博士生,主要研究方向为高性能计算、生物信息学等,E2mail:[email protected];谢　江　博士,副教授,主要研究方向为高性能计算、生物信息学等;陈建娇　博士生,主要研究方向为生物信息学等;张　武　教授,博士生导师,主要研究方向为高性能计算、生物信息学等。

络特征。

1)网络度分布:在随机一致的原则下挑选出的节点的度数为k的概率为p(k)。

2)网络拓扑结构:无标度网络模型虽然最符合现实复杂网络的论断仍有争议,但依然被认为是目前最接近现实网络特性的模型,因此我们的模拟器采用的是无标度模型。

3)人工网络的鲁棒性:模拟出的网络必须符合动力学稳定性。

GN2Simulator的意义在于,理论上将复杂网络模型与经

或是连线图。这包括两方面的工作:一是决定网络中节点(基因)的数量,二是确定对每个节点产生激活或抑制效应的是哪些节点,即每个节点的出入度。在对43个物种的基因网络的度分布进行统计分析[7]后发现,这些网络的入度和出度分布均服从幂率分布:p(k)～k-r,幂指数r约为2.2,为基因网络的无标度拓扑特性提供了又一有力证据。同时,研究表明,当网络的平均连接度在2～3时,网络中的大部分具有生物学意义的动力学行为在此区间产生[8]。如果仅关注某一执行具体功能的子模块,其平均连接度相对高于整个网络的平均连接度,约在2～4之间。因此,在构建基因网络模拟器(GN2Simu2lator)时,平均连接度的默认区间为k=2～4。另外,人工基因网络的构建也考虑了无标度网络演化的BA模型(Barabase&Albertmodel)的内在机制:网络规模的增长和择优连结。根据这种无标度特性,网络不断有新节点加入,且倾向于与度数较高的节点连接,那么,网络中连接度高的节点往往是网络形成过程中较早加入的节点。基于真实基因网络以上特性的分析,模拟真实基因网络的过程概括如下:

步骤1　基因网络初始状态:m0个节点,e0条边,m0个节点被随机分为l个子集,那么每个子集中各有m01,m02,…,

m0l个节点,e01,e02,…,e0l条边过大量实验证明的真实基因调控网络特征很好地结合起来,

使得模拟网络具有较高的可信度;技术上它为基因调控网络甚至生物网络的模拟提供了有效的一体化解决途径。在更为重要的实际应用中,许多构建基因网络的算法需要得到详细全面的测试,虽然在真实基因网络上可以进行直接而简单的测试实验,但由于已知真实网络的规模有限、类型很少、动力学机制不明等缺点,很容易对特定算法产生偏袒,那么由GN2

Simulator模拟产生的大规模、多类型、参数可控的高可信基

因调控网络无疑成为公平系统地验证算法的金标准(goldstandard)。本文将开展构建人工基因调控网络的研究,并对基因调控网络的模拟结果进行动力学稳定性分析,给出其应用于算法评估的步骤。

2　基因调控网络模拟器构建

2.1　人工网络的矩阵表示

步骤2　mm条出边,其mm01,,0l)。

。当新节点以概率p产生m条入,,然后在该子集中偏好选取m-1个节点与之相连。原则为:择优选择连接节点的概率取决于该节点的入度。与此同时新节点的另一条边也在择优连接的原则下与另一个子集的节点相连。

步骤4　产生新的出边。当新节点以概率q产生m条出边时,也随机为新节点在初始网络中选择一个子集,并在该子集中偏好选取m-1个节点与之相连。原则为:择优选择连接节点的概率取决于该节点的出度。同时新节点的另一条边也在择优连接的原则下与另一个子集的节点相连。步骤5　重复步骤2-4,同时统计节点数,直到网络增长到需要的规模为止。

通过以上步骤,初步构建了GN2Simulator。

2.3　人工网络的鲁棒性和动力学稳定性

即把基因简化为节点(vertex),的连线(edge),及其节点间边的集合E构成的:

Network=(V,E)

由于邻接矩阵可以用来描述一个网络中节点与节点之间的关系,因此网络的拓扑结构用邻接矩阵A表示:

AN×N=

a21

a11

…a1N…a2Nω…

…0

…

aN1

…

aN2

其中,A中的每一个元素aij代表基因之间的调控关系。考虑到真实调控网络中每个基因的状态改变主要取决于其他基因对它的作用,为了简化分析,基因自身的自调控可以被忽略,即A中对角元素均为0。基因之间的调控关系有两种:激活和抑制,真实调控网络中基因之间的抑制作用通常强于激活作用,且aij的具体值并不影响网络的整体动力学性质,因此以aij=-∞表示抑制作用,aij=1表示激活作用。

此外,在网络的存储方面,由于真实的生物网络大多是稀疏网络,如芽殖酵母的蛋白质相互作用网络的平度连结度仅约为2[6],在储存对应的稀疏矩阵时,如果直接用邻接矩阵表示法,对于一个有n个节点的网络,将占用n2的存储空间,当模拟大规模的全局基因网络时,将产生过高的空间复杂度;因此,采用边列表示法建立了一个以边为元素的数组来记录每条边的出节点和入节点,从而使得整个网络需要占用的存储空间下降为2n。2.2　给定度分布及基于无标度拓扑模型的网络构建

为了定义一个基因网络模型,必须建立网络的拓扑结构,

生物系统和调控网络的基本性质是鲁棒性和稳定性[9],这两种重要的性质有利于生物应对复杂多变的外界环境和不断受到扰动的内部环境,甚至在某些基因缺失的情况下,整体

的生物状态和重要的基本功能还能够保持稳定。因此,如何构建出具有动力学稳定性的人工网络是一项重要的研究内容。我们必须首先清楚地分析生物调控网络是如何实现鲁棒性和稳定性的。一方面通过网络的结构性质[10],例如基因的冗余性、功能模块化、网络中的负反馈机制,可以实现稳定性。另一方面,生物系统可以通过网络的整体结构和动力学性质来实现鲁棒性[11]。这种由网络的整体性质产生的鲁棒性是本文关注的重点。

在建立好基因调控网络的理论模型后,将进一步研究其网络的动力学稳定性,包括状态稳定性和结构稳定性。状态稳定性指研究网络中不同节点所代表的基因数量发生变化

时,网络所执行的生物学过程能否继续。结构稳定性指网络中不同节点的具体数值所代表的基因相互作用的强度发生改变时生物学状态的稳定性。

在模拟基因调控网络的过程中,为了构建具有动力学稳定性的人工网络,我们进行了具有全局稳定不动点的骨架网络的设计。在一个网络中,具有最少的边数且对整个网络的整体动力学性质起决定性作用的一个连通子集被称为该网络的骨架网。为了保证网络的整体动力学性质,这个特殊的子网是不能随意改变的。在网络的除骨架网之外的拓扑结构上进行随机扰动,例如,随机删除或增加一条边,抑或改变一条边的属性,如果网络仍然定性地保持其整体动力学性质不变,就成功构建了一个拓扑结构上具有鲁棒性的人工网络[12]。

此外,为验证模拟网络的稳定性,可以计算其相关矩阵的特征根,如果所有特征根实部为负或最大特征根实部为负,系统必将趋于稳定。同时,最大负特征根的绝对值越大,系统在经历微小的扰动后趋于稳定的速度将越快。

3　数值试验及应用分析

如前所述,GN2simulator

建立的一个重要目的在于为全面评估基因调控关系的预测算法提供“金标准”人工网络。首先,在应用于算法评估之前我们验证了人工基因调控网络的动力学稳定性。基于GN2simulator,本文模拟产生了10个有着不同平均连接度的基因调控网络样本,其平均连接度k=3.9,抑制作用关系与激活作用关系之比r=0.7,节点数n=100,图1一个微小的扰动后,10收敛并达到平衡状态下。

1)网络参数设置:根据用户选择,生成n个具有不同网络规模、平均度(average2degree)、平均路径长度(averagepathlength)、群聚系数(clusteringcoefficient)和介数(vertexbe2tweennessandedgebetweenness)的人工基因调控网络。

2)网络模拟:根据每个网络生成相应的邻接矩阵,记为A1,A2,…,An,以稀疏矩阵的形式表示,并以边列表的形式储存。

3)稳定性检验:随机选取矩阵中一定数量的节点,将其数值增加0.1(保证扰动足够微小,以免使整个系统偏离平衡点过大),检测整体网络在一定时间步内的收敛性,如果不符合动力学稳定性,可进行修正:将对角元素减去矩阵的最大特征值,这样既不会改变网络的拓扑结构和度分布等特性,同时又保证了系统的稳定。

4)结果比对:用户以矩阵形式提交其算法预测出来的基因调控网络,记为A1′,A2′,…,An′,将A与A′进行相似性比对,给出正确率、召回率和覆盖率等统计性比较结果。表2给出了不同网络参数下几种典型算法的表现情况。

表2　不同人工网络测试样本下的算法表现

AlgorithmsSWNI

EO[15]

Coverage(Numberofnodes,averagedegree)(50,2)65%8476%(200,2)9577%

,2)%58%53%

(500,3)91%64%58%57%

(500,4)85%68%59%60%

图1　微扰后网络中所有节点随时间收敛图

同时,从著名的基因调控和基因组数据库KEGG[13]及

NCBI[14]数据库中收集了人类转录因子调控网络、酵母基因调控网络和蛋白质组网络,将它们的网络拓扑性质与GN2Simulator模拟的网络进行对比,如表1所列,本文模拟的基因调控网络与真实网络相比具有高相似性。

表1　拓扑性质比较

GNSimulator

AveragedegreeAverageclusteringAveragepathlengthAssortativemixing

30.134.7-0.16

HTFN3.70.174.5-0.18

Yeast3.70.154.15-0.05

Proteome2.40.76.81-0.15

2可以看出,不同的算法均对网络规模和网络平均

度敏感,但敏感趋势和敏感度不同。SWNI在大规模稀疏网络的预测方面相对其它算法表现较优,但其预测覆盖率却随着网络稠密化而降低。另外,可以明显看出ARACNE等3种算法的预测覆盖率随着网络规模的增大而下降,但似乎在对较稠密网络的预测中有较好的表现。

结束语　本文提出了专门针对基因调控网络的人工网络建立模型。与其他的复杂网络分析软件Pajek,Netdraw和Ucient等比较,GN2Simulator具有很强的针对性和生物信息学应用性,为系统生物领域的研究人员提供了基因网络模拟工具。数值试验结果显示,GN2Simulator产生的网络与真实基因调控网络相比,在整体结构和动力学性质上有很高的相似性,能弥补真实网络规模不足等缺点,为各种基因调控关系预测算法提供了大量多种类的测试网络数据集;同时,它集成了算法评估、比较和分析功能,对算法在不同网络数据集上的表现进行打分,为生物网络结构及其预测算法提供了一个高效、合理的计算实验平台。

参考文献

[1]

JeongH,TomborB,AlbertR,etal.Thelarge2scaleorganizationofmetabolicnetworks[J].Nature,2000,407:6512654[2]

Erd

第37卷　第1期计算机科学Vol.37No.1基于复杂网络模型的基因调控网络的计算模拟

张律文1,2　谢　江1,2　陈建娇1　张　武1

(上海大学计算机工程与科学学院　上海200072)1　(上海大学系统生物技术研究所　上海200444)2

摘　要　随着大通量基因芯片数据的产生,基因调控机制的网络化研究需求日趋迫切。提出了基于复杂网络理论的基因调控网络的模拟方法,构建了基因调控网络模拟器GN2Simulator。通过分析真实基因调控网络的拓扑特性,提出了对应的矩阵模型,并充分考虑了网络的生物学鲁棒性和动力学稳定性,给出了人工基因网络的生成过程和计算模拟方法。计算实验表明,GN2Simulator能高效地模拟出与真实基因调控网络高度相似的大规模人工网络,并可为不同算法提供无偏验证的多样化人工模拟数据。关键词　复杂网络,无标度特性,动力学稳定性,模拟,基因调控网络

ArtificialGeneRegulatoryNetworksConstructionBasedonComplexNetwork

ZHANGLu2wen1,2　XIEJiang1,2　CHENJian2jiao1　ZHANGWu1

(SchoolofComputerEngineeringandScience,ShanghaiUniversity,Shanghai200072,China)1

(InstituteofSystemsBiology,ShanghaiUniversity,Shanghai200444,China)2

Abstract　Alongwiththeexplosionofhigh2throughputgenomicinformationproducedfDNAmicroarrays,modelingthemechanismofgeneregulationsandconstructinggenetask.Wepro2posedanovelmethod,GN2Simulator,tosimulategeneregulatoryontheory.TherealgenenetworksweremodeledexplicitlyinmatrixesAccordingtotherobustbiologicalmechanism,weconstructedartificialgenetothedynamicstability.WeconsideredthefeaturesofrealgenenetworkThecomputingexperimentsillustratethatGN2Simulatorcansimulatewithhighconfidencelevel.Moreover,GN2Simulatorgenera2tesvariousdifferentalgorithmsperformanceandprovidesreasonableestimationonthem.

Keywords　Complexnetwork,Scale2freetopology,Dynamicstability,Simulation,Generegulatorynetwork　

究成果之后,复杂网络被广泛地应用于生物网络,包括基因调控网络的构建和模拟等各方面,使得基因网络的研究取得了极大进展。

用来模拟真实基因调控的复杂网络模型经历了几个阶段的发展。最初,简单的随机网络模型[2]被用来描述基因网络的结构,但在1998年,Watts和Strogats在结合规则网络和随机网络特点的基础上,建立了小世界(Small2world)网络模型[3],接着1999年,Barabase和Albert发现了真实网络的无标度(Scale2free)性质[4]。在对大规模基因网络进行数据采集和统计分析后,无标度被认为是最接近真实基因网络的一种拓扑结构[5],基因相互作用网络的结构特性可能与其他复杂系统网络(比如Internet网)在很大程度上是一致的。

在构建人工基因网络的过程中,仅仅关注网络的大体拓扑结构是远远不够的。为了尽可能地模拟出接近真实基因网络的人工网络,构建了基因调控网络模拟器Genenetwork

simulator(GN2Simulator),从以下几个方面系统地考虑了网

1　引言

复杂生命现象是大量基因互相调控、协同作用的结果。随着基因组测序的完成、基因组学的快速发展,生物学积累了大量的基因调控实验数据,如何挖掘出海量数据所蕴藏的生命现象和生物规律已成为生命科学的研究热点。如今,复杂生物系统的研究开始从对细胞内个别基因功能的局部性描述发展到对复杂基因调控网络(GeneRegulatoryNetwork)的定量刻画。由于基因与其他生物小分子倾向于成组地通过网状的相互作用而影响生物系统的功能,因此对基因功能的研究必须分析其调控网络。

关于基因调控机制的系统性研究,最行之有效的方法就是结合计算机科学和数学等理论把基因之间的复杂相互作用进行整合,并简化成网络,然后建立理论模型,从而研究基因调控网络的结构、网络模块及动力学性质。自2000年Nature上发表了利用复杂网络就理论研究生物网络拓扑特性[1]的研

到稿日期:2009202220　返修日期:2009205208　　本文受上海市重点学科建设项目(项目编号:J50103),上海大学系统生物研究基金(200825)资助。

张律文　博士生,主要研究方向为高性能计算、生物信息学等,E2mail:[email protected];谢　江　博士,副教授,主要研究方向为高性能计算、生物信息学等;陈建娇　博士生,主要研究方向为生物信息学等;张　武　教授,博士生导师,主要研究方向为高性能计算、生物信息学等。

络特征。

1)网络度分布:在随机一致的原则下挑选出的节点的度数为k的概率为p(k)。

2)网络拓扑结构:无标度网络模型虽然最符合现实复杂网络的论断仍有争议,但依然被认为是目前最接近现实网络特性的模型,因此我们的模拟器采用的是无标度模型。

3)人工网络的鲁棒性:模拟出的网络必须符合动力学稳定性。

GN2Simulator的意义在于,理论上将复杂网络模型与经

或是连线图。这包括两方面的工作:一是决定网络中节点(基因)的数量,二是确定对每个节点产生激活或抑制效应的是哪些节点,即每个节点的出入度。在对43个物种的基因网络的度分布进行统计分析[7]后发现,这些网络的入度和出度分布均服从幂率分布:p(k)～k-r,幂指数r约为2.2,为基因网络的无标度拓扑特性提供了又一有力证据。同时,研究表明,当网络的平均连接度在2～3时,网络中的大部分具有生物学意义的动力学行为在此区间产生[8]。如果仅关注某一执行具体功能的子模块,其平均连接度相对高于整个网络的平均连接度,约在2～4之间。因此,在构建基因网络模拟器(GN2Simu2lator)时,平均连接度的默认区间为k=2～4。另外,人工基因网络的构建也考虑了无标度网络演化的BA模型(Barabase&Albertmodel)的内在机制:网络规模的增长和择优连结。根据这种无标度特性,网络不断有新节点加入,且倾向于与度数较高的节点连接,那么,网络中连接度高的节点往往是网络形成过程中较早加入的节点。基于真实基因网络以上特性的分析,模拟真实基因网络的过程概括如下:

步骤1　基因网络初始状态:m0个节点,e0条边,m0个节点被随机分为l个子集,那么每个子集中各有m01,m02,…,

m0l个节点,e01,e02,…,e0l条边过大量实验证明的真实基因调控网络特征很好地结合起来,

使得模拟网络具有较高的可信度;技术上它为基因调控网络甚至生物网络的模拟提供了有效的一体化解决途径。在更为重要的实际应用中,许多构建基因网络的算法需要得到详细全面的测试,虽然在真实基因网络上可以进行直接而简单的测试实验,但由于已知真实网络的规模有限、类型很少、动力学机制不明等缺点,很容易对特定算法产生偏袒,那么由GN2

Simulator模拟产生的大规模、多类型、参数可控的高可信基

因调控网络无疑成为公平系统地验证算法的金标准(goldstandard)。本文将开展构建人工基因调控网络的研究,并对基因调控网络的模拟结果进行动力学稳定性分析,给出其应用于算法评估的步骤。

2　基因调控网络模拟器构建

2.1　人工网络的矩阵表示

步骤2　mm条出边,其mm01,,0l)。

。当新节点以概率p产生m条入,,然后在该子集中偏好选取m-1个节点与之相连。原则为:择优选择连接节点的概率取决于该节点的入度。与此同时新节点的另一条边也在择优连接的原则下与另一个子集的节点相连。

步骤4　产生新的出边。当新节点以概率q产生m条出边时,也随机为新节点在初始网络中选择一个子集,并在该子集中偏好选取m-1个节点与之相连。原则为:择优选择连接节点的概率取决于该节点的出度。同时新节点的另一条边也在择优连接的原则下与另一个子集的节点相连。步骤5　重复步骤2-4,同时统计节点数,直到网络增长到需要的规模为止。

通过以上步骤,初步构建了GN2Simulator。

2.3　人工网络的鲁棒性和动力学稳定性

即把基因简化为节点(vertex),的连线(edge),及其节点间边的集合E构成的:

Network=(V,E)

由于邻接矩阵可以用来描述一个网络中节点与节点之间的关系,因此网络的拓扑结构用邻接矩阵A表示:

AN×N=

a21

a11

…a1N…a2Nω…

…0

…

aN1

…

aN2

其中,A中的每一个元素aij代表基因之间的调控关系。考虑到真实调控网络中每个基因的状态改变主要取决于其他基因对它的作用,为了简化分析,基因自身的自调控可以被忽略,即A中对角元素均为0。基因之间的调控关系有两种:激活和抑制,真实调控网络中基因之间的抑制作用通常强于激活作用,且aij的具体值并不影响网络的整体动力学性质,因此以aij=-∞表示抑制作用,aij=1表示激活作用。

此外,在网络的存储方面,由于真实的生物网络大多是稀疏网络,如芽殖酵母的蛋白质相互作用网络的平度连结度仅约为2[6],在储存对应的稀疏矩阵时,如果直接用邻接矩阵表示法,对于一个有n个节点的网络,将占用n2的存储空间,当模拟大规模的全局基因网络时,将产生过高的空间复杂度;因此,采用边列表示法建立了一个以边为元素的数组来记录每条边的出节点和入节点,从而使得整个网络需要占用的存储空间下降为2n。2.2　给定度分布及基于无标度拓扑模型的网络构建

为了定义一个基因网络模型,必须建立网络的拓扑结构,

生物系统和调控网络的基本性质是鲁棒性和稳定性[9],这两种重要的性质有利于生物应对复杂多变的外界环境和不断受到扰动的内部环境,甚至在某些基因缺失的情况下,整体

的生物状态和重要的基本功能还能够保持稳定。因此,如何构建出具有动力学稳定性的人工网络是一项重要的研究内容。我们必须首先清楚地分析生物调控网络是如何实现鲁棒性和稳定性的。一方面通过网络的结构性质[10],例如基因的冗余性、功能模块化、网络中的负反馈机制,可以实现稳定性。另一方面,生物系统可以通过网络的整体结构和动力学性质来实现鲁棒性[11]。这种由网络的整体性质产生的鲁棒性是本文关注的重点。

在建立好基因调控网络的理论模型后,将进一步研究其网络的动力学稳定性,包括状态稳定性和结构稳定性。状态稳定性指研究网络中不同节点所代表的基因数量发生变化

时,网络所执行的生物学过程能否继续。结构稳定性指网络中不同节点的具体数值所代表的基因相互作用的强度发生改变时生物学状态的稳定性。

在模拟基因调控网络的过程中,为了构建具有动力学稳定性的人工网络,我们进行了具有全局稳定不动点的骨架网络的设计。在一个网络中,具有最少的边数且对整个网络的整体动力学性质起决定性作用的一个连通子集被称为该网络的骨架网。为了保证网络的整体动力学性质,这个特殊的子网是不能随意改变的。在网络的除骨架网之外的拓扑结构上进行随机扰动,例如,随机删除或增加一条边,抑或改变一条边的属性,如果网络仍然定性地保持其整体动力学性质不变,就成功构建了一个拓扑结构上具有鲁棒性的人工网络[12]。

此外,为验证模拟网络的稳定性,可以计算其相关矩阵的特征根,如果所有特征根实部为负或最大特征根实部为负,系统必将趋于稳定。同时,最大负特征根的绝对值越大,系统在经历微小的扰动后趋于稳定的速度将越快。

3　数值试验及应用分析

如前所述,GN2simulator

建立的一个重要目的在于为全面评估基因调控关系的预测算法提供“金标准”人工网络。首先,在应用于算法评估之前我们验证了人工基因调控网络的动力学稳定性。基于GN2simulator,本文模拟产生了10个有着不同平均连接度的基因调控网络样本,其平均连接度k=3.9,抑制作用关系与激活作用关系之比r=0.7,节点数n=100,图1一个微小的扰动后,10收敛并达到平衡状态下。

1)网络参数设置:根据用户选择,生成n个具有不同网络规模、平均度(average2degree)、平均路径长度(averagepathlength)、群聚系数(clusteringcoefficient)和介数(vertexbe2tweennessandedgebetweenness)的人工基因调控网络。

2)网络模拟:根据每个网络生成相应的邻接矩阵,记为A1,A2,…,An,以稀疏矩阵的形式表示,并以边列表的形式储存。

3)稳定性检验:随机选取矩阵中一定数量的节点,将其数值增加0.1(保证扰动足够微小,以免使整个系统偏离平衡点过大),检测整体网络在一定时间步内的收敛性,如果不符合动力学稳定性,可进行修正:将对角元素减去矩阵的最大特征值,这样既不会改变网络的拓扑结构和度分布等特性,同时又保证了系统的稳定。

4)结果比对:用户以矩阵形式提交其算法预测出来的基因调控网络,记为A1′,A2′,…,An′,将A与A′进行相似性比对,给出正确率、召回率和覆盖率等统计性比较结果。表2给出了不同网络参数下几种典型算法的表现情况。

表2　不同人工网络测试样本下的算法表现

AlgorithmsSWNI

EO[15]

Coverage(Numberofnodes,averagedegree)(50,2)65%8476%(200,2)9577%

,2)%58%53%

(500,3)91%64%58%57%

(500,4)85%68%59%60%

图1　微扰后网络中所有节点随时间收敛图

同时,从著名的基因调控和基因组数据库KEGG[13]及

NCBI[14]数据库中收集了人类转录因子调控网络、酵母基因调控网络和蛋白质组网络,将它们的网络拓扑性质与GN2Simulator模拟的网络进行对比,如表1所列,本文模拟的基因调控网络与真实网络相比具有高相似性。

表1　拓扑性质比较

GNSimulator

AveragedegreeAverageclusteringAveragepathlengthAssortativemixing

30.134.7-0.16

HTFN3.70.174.5-0.18

Yeast3.70.154.15-0.05

Proteome2.40.76.81-0.15

2可以看出,不同的算法均对网络规模和网络平均

度敏感,但敏感趋势和敏感度不同。SWNI在大规模稀疏网络的预测方面相对其它算法表现较优,但其预测覆盖率却随着网络稠密化而降低。另外,可以明显看出ARACNE等3种算法的预测覆盖率随着网络规模的增大而下降,但似乎在对较稠密网络的预测中有较好的表现。

结束语　本文提出了专门针对基因调控网络的人工网络建立模型。与其他的复杂网络分析软件Pajek,Netdraw和Ucient等比较,GN2Simulator具有很强的针对性和生物信息学应用性,为系统生物领域的研究人员提供了基因网络模拟工具。数值试验结果显示,GN2Simulator产生的网络与真实基因调控网络相比,在整体结构和动力学性质上有很高的相似性,能弥补真实网络规模不足等缺点,为各种基因调控关系预测算法提供了大量多种类的测试网络数据集;同时,它集成了算法评估、比较和分析功能,对算法在不同网络数据集上的表现进行打分,为生物网络结构及其预测算法提供了一个高效、合理的计算实验平台。

参考文献

[1]

JeongH,TomborB,AlbertR,etal.Thelarge2scaleorganizationofmetabolicnetworks[J].Nature,2000,407:6512654[2]

Erd