基于相对距离的复杂网络谱粗粒化方法*

杨青林 王立夫 李欢 余牧舟

(东北大学秦皇岛分校 控制工程学院，秦皇岛 066004)

1 引 言

从20世纪末开始,复杂网络理论逐步成为社会、经济、化学和信息系统等领域的重要研究方法,对复杂网络静态特性和动态特性的定量定性描述,已经成为网络时代一个极其重要的挑战性课题,被称为“网络的新科学(new science of network)”[1].复杂网络的同步作为一种重要的网络动态行为,在激光系统、超导材料和通信等领域有重要的意义.同步现象在自然界中也广泛存在,例如萤火虫有规律地闪光、心脏细胞同步震荡、挂在同一横梁上的钟摆同步摆动等.近20多年来,研究人员对复杂网络的同步研究取得了丰硕的成果[2-15].其中,对于一般连续时间耦合网络最常用的方法是将网络的状态方程线性化后转化成主稳定函数来研究[2,3].然而对于实际网络动辄成千上万的节点来说,判断他们是否同步需要求解大量耦合微分方程,并且小规模网络得到的同步判据或定理不能直接推广到大规模网络[16].因此,如何将复杂网络进行约简并保留其同步能力不变,对大规模网络的同步分析和仿真至关重要.Gfeller和Rios[17,18]提出的谱粗粒化(spectral coarse graining,SCG)算法通过分析网络的结构矩阵,把特征向量分量相同或相近的节点合并,从而保证了约简后网络的同步能力不变.周建等[19]提出的改进谱粗粒化(improved spectral coarse graining,ISCG)算法采用了不同的分类方式,改善了在特征向量分量分布极不均匀时的粗粒化效果,并且降低了算法的复杂度.Chen等[20]进一步研究了网络的聚类结构对谱粗粒化效果的影响.Wang和Xu[21]发现采用SCG算法对网络约简后,其可控性也能很好地保留下来.

SCG算法和ISCG算法都给出了一种简单可实现的谱粗粒化方法,但是在确定网络中哪些节点需要合并时都是以特征向量分量间的绝对距离为判断依据,没有考虑到特征向量分量间的相对距离,即绝对距离占分量绝对值的比例.本文经过理论推导证明了两个节点合并前后特征值的变化量与分量间的相对距离成正比例关系,SCG和ISCG算法在分类时都忽略了分量间相对距离的影响,从而导致节点的分类不准确.为了克服这一缺陷,本文提出了一种基于相对距离的谱粗粒化(improved spectral coarse graining based on relative distance,ISCGR)算法,在对节点分类时以特征向量分量间的相对距离作为判断标准,改善了网络的粗粒化效果,提高了网络约简后保持同步的能力.采用本文所提的算法对27种实际网络进行粗粒化仿真,发现互联网、生物、社交、合作等具有明显聚类结构的网络在约简后比电力、化学等网络更容易保持同步.

2 同步能力的刻画

主稳定方程法是研究一般连续时间耦合网络同步的常用方法[2,3].考虑由N个节点构成的复杂网络,其中第i个节点的动力学方程为

式中,xi∈Rm为节点i的m维状态变量；f(xi)为第i个节点自身的状态方程；常数µ＞0 是网络的耦合强度；H(·):Rm→Rm表示各个节点之间的内部耦合函数；Laplacian矩阵L=(lij)∈RN×N是外部耦合矩阵,表示网络的拓扑结构,满足耗散条件特别的,当Laplacian矩阵L表示一个无向无权的网络时,有如下定义：

对角元为

其中ki为节点i的度数.此时,外部耦合矩阵L是一个不可约的对称阵,具有非负特征根并且特征根满足 0=λ1≤λ2≤···≤λN.如果在t→∞ 时有x1(t)=x2(t)=···=xN(t)=s(t),那么就称动态网络(1)式能够达到完全同步,其中同步流形s(t)满足(t)=f(s(t)).对(1)式在同步状态s(t)处线性化,令ξi(t)=xi(t)-s(t),可以得到如下线性动力学方程：

其中Df(s)和DH(s)分别是f(s)和H(s)在同步s(t)处的Jacobi(雅可比)矩阵,令ξ=[ξ1,ξ2···ξN],则(4)式可以写为

记LT=JΛJ¯1为Laplace矩阵L的Jordan分解,Λ=diag[λ1,λ2,···,λN]为L的特征值矩阵.再令ψ=[ψ1,ψ2,···ψN]=ξJ,将(5)式对角化可得

(6)式是一个关于λk的N维微分方程组,展开后可以写成

将(7)式一般化后即可得到主稳定方程：

这里,α=µλk.(8)式的最大Lyapunov指数Lmax是变量α的函数,称为动力学网络(1)式的主稳定函数.给定一耦合强度µ,在坐标轴上可以相应地找到一点µλk,若该点对应的Lmax为负数时,则该特征模态为稳定态；反之,则该特征模态为不稳定态.如果与特征值λk(k=2,3,···,N)对应的最大Lyapunov指数Lmax均为负数,那么在该耦合强度µ下,整个网络的同步流形是渐近稳定的.我们把使得Lmax＜0 的α取值范围S称为动态网络(1)式的同步化区域,它由孤立节点的动力学函数f(xi)和内部耦合函数H(·)决定.如果耦合强度µ与Laplacian矩阵的每个非零特征值之积都在同步区域S内,即µλk∈S(k=2,3,···,N),那么动态网络(1)式的同步流形是渐近稳定的,网络就能达到同步[2,3].同步区域S可以分为以下4种类型：Ⅰ.S=(α1,+∞)；Ⅱ.S=(α1,α2)；Ⅲ.S=(α1,α2)∪(α2,α3)∪···；Ⅳ.S=∅ .其中,类型Ⅰ网络的同步化能力可以用Laplacian矩阵L的最小非零特征值λ2来刻画,λ2越大,µλ2越容易大于α1,网络的同步化能力越强；类型Ⅱ网络的同步能力可以用Laplacian矩阵L的最大特征值λN与λ2的比值λN/λ2来刻画,λN/λ2越小,网络的同步能力越强；类型Ⅲ网络需要所有特征值都满足一定的条件,很难达到同步；类型Ⅳ网络无论耦合强度µ和Laplacian矩阵L的特征值为何值,网络都不能达到同步.大多数网络都属于Ⅰ或Ⅱ两种类型,因此在分析网络的同步能力时,通常考虑λ2和λN/λ2两个指标[22].

3 网络同步的SCG算法和ISCG算法

3.1 SCG方法与ISCG方法的原理

由主稳定方程法可知,当动态网络节点间的耦合强度、孤立节点的动力学方程和内部耦合函数固定后,网络的同步能力由Laplacian矩阵的λ2或λN/λ2决定.所以要想保持网络的同步能力不变,需要尽量保持网络在粗粒化前后的λ2或λN/λ2不变.Dffller和Rios[17,18]提出的SCG算法较好地保留了约简前后网络的λ2和λN/λ2值,从而保持了约简前后网络的同步能力不变.这种算法主要解决了合并过程中的两个主要问题,一是如何合并节点及更新边,即如何更新约简后网络的Laplacian矩阵；二是确定哪些节点可以被合并.

首先,对于如何更新合并后的边,我们先看一个小型的网络.图1(a)是一个无向无权的网络,如果把节点c、d、e合并为图1(b)中的节点g,则称节点c、d、e组成了团Cg,记ωx→y是节点x到节点y的权值,我们有：

这里,|Cg| 表示团Cg中包含的节点个数.同理,我们也可以求出其他的边权.在大型网络中,记初始节点的标号i=1,2,···,N.约简后的节点标号为C=1,2,···,,为约简之后网络的节点数,C也可以表示原始网络中要约简的节点所在的团,同属于一个团的节点被合并成一个节点.约简之后网络的Laplacian矩阵可以用以下公式求解

图1 合并节点的过程Fig.1.The processing of merging nodes.

这里,K∈R×N和R∈RN×分别为

其中,Ci表示原始网络节点i所在团的标号；|C| 表示团C中节点的个数,δ是常用的Kronecker函数.

接下来考虑合并哪些边,分别对λ2和λN/λ2两种情况进行讨论.若要保持网络(1)式的Laplacian矩阵L的最小非零特征值λ2不变,SCG方法的原理是将λ2所对应特征向量P2中分量相同或相近的节点合并在一起,这样得到合并后网络的Laplacian矩阵的最小非零特征值2与原来网络的λ2相同或相近.通常,我们定义特征向量P2所对应的两个分量p2i和p2j相近是指‖p2i-p2j‖/‖p2max-p2min‖≪1,其中p2max和p2min分别是特征向量分量中的最大值和最小值.具体的操作步骤为,先把区间 [p2min,p2max] 平均分成i个区间,再把分在同一个区间内的特征向量分量对应的节点合并为一个节点,最后利用(11)式求解约简后网络的Laplacian矩阵.类似地,当要保持网络(1)的Laplacian矩阵L的最大特征值与第二小特征值之比λN/λ2时,需要同时保持λ2和λN不变.具体的步骤为：先把λ2所对应的特征向量分量分为i个等分区间,再把λN所对应的特征向量分量分为i个等分区间,找出P2和PN中同时落入等分区间的分量,将这些分量对应的节点合并,最后利用(11)式更新合并后的Laplacian矩阵,便可以得到约简后的网络.

图1表明,采用谱粗粒化思想把外部连接相同的节点合并,网络的第二小特征值λ2不变,即约简后网络的同步能力不发生改变.

ISCG方法是谱粗粒化方法的改进算法[19],也是通过合并特征向量P2中相同或相近的分量达到约简网络的目的.但与SCG算法不同的是,ISCG算法不是把特征向量分量等间距平分,而是采用分裂聚类的方法对节点进行分类.具体操作为,当要保持λ2不变时,先把λ2所对应的特征向量分量从小到大排列,算出两两相邻分量间的距离,找出最大的前-1 个间距,在这-1 个间距设置分割点,将特征向量分量分裂为个聚类,把属于同一聚类的节点合并为一个节点,最后更新Laplacian矩阵,就得到了节点数为的粗粒化网络.同理,当要保持λN/λ2不变时,先根据特征向量分量的间距分别将P2的分量和PN的分量分裂为个聚类,找出同时被分裂到同一聚类里的分量,把分量对应的节点合并,最后按照(11)式来更新合并后网络的Laplacian矩阵.ISCG算法的优点在于可以将原始网络约简为任意指定规模的粗粒化网络,复杂度为O(N2),要低于SCG算法的复杂度O(N3),并且在特征向量分量分布不均匀时,使用ISCG算法的粗粒化效果要明显优于SCG算法[19].

3.2 SCG算法和ISCG算法的局限

SCG算法和ISCG算法在保持网络的同步能力时都有着得天独厚的优势.现实世界中大部分特征向量分量分布非常集中,少数特征向量分量分布特别分散,对于这种特征向量分量分布不均匀的网络采用SCG算法一方面会出现距离很近的特征向量分量被分到两个聚类里面,另一方面由于少数特征向量分布较为分散,导致得到的聚类数随i变化不大,不容易得到任意规模的粗粒化网络[19].而ISCG方法采用分裂聚类的方式对节点分类,因此,在特征向量P2的分量极不均匀时,采用ISCG算法约简网络更能保持其同步能力.

图2 15个节点并为14个节点的两种方案Fig.2.Two schemes of merging 15 nodes into 14 nodes.

这两种算法分类的方式不同,但都是以特征向量分量的绝对距离作为依据对节点分类.如果λ2对应的特征向量P2的第i个分量和第j个分量相等,即p2i=p2j,按照SCG算法或ISCG算法把这两个节点合并,λ2的值在网络约简前后不发生改变.而当P2的两个分量有较小的距离时,两种算法的分类方式都遵循两个分量间距离越小,被分在一起的几率越大,而没有考虑分量间相对距离的因素.如果设分量p2i与p2j之间的间距为ε,则p2i与p2j间的相对距离就是指ε占 min{|p2i|,|p2j|} 的比例ρ=ε/min{|p2i|,|p2j|}.下面我们通过一个例子说明在分类时只考虑绝对距离是不合理的.

假设原始网络如图2所示,网络的Laplacian矩阵为

其最小非零特征值λ2=0.4130,对应的特征向量P2=[-4.19,2.87,-1.92,3.61,-2.80,2.87,-1.14,0.11,1.51,1.00,2.80,-3.85,1.71,0.70,-3.27]T,P2的第9个与第10个分量分别是1.51和1.00,相差0.51；第12个与第15个分量分别是 -3.85 和 -3.27 相差0.58.如果以绝对距离作为合并节点的依据,则合并节点9和10后λ2的变化量会更小.实际上,由(11)式得到合并节点9和10后的Laplacian矩阵1为

根据以上的讨论,我们提出约简后网络λ2的变化量与所合并节点对应的特征向量分量间相对距离的关系.

定理1：在进行谱粗粒化过程合并两个节点时,如果同一个特征值对应的两个特征向量分量p2i和p2j的间距ε远远小于 min{|p2i|,|p2j|},则合并节点i和节点j后网络特征值的变化量δ与p2i和p2j间的相对距离ε/min{|p2i|,|p2j|} 成正比.

证明：假设λ2对应的两个正的特征向量分量分别为为p2i与p2j(p2i＜p2j),第j个分量与第i个分量之间的距离为ε,即p2j=p2i+ε,且ε≪p2i.合并节点i和节点j后λ2的改变量设为δ,写出原始网络Laplacian矩阵的特征方程

把(13)式的第i个等式和第j个等式两边相加得到

因为分量之间的距离ε≪p2i,所以可以忽略ε的大小,把合并后得到新节点所对应的特征向量分量近似为p2i.按照(11)式求出合并后网络的Laplacian矩阵并写出特征方程

(15)式的第i个等式为

由(14)式和(16)式联立可以得到

(17)式表明,约简网络后λ2的改变量δ不仅与特征向量分量间的绝对距离ε有关,还与特征向量分量p2i的大小有关.在ε相对于p2i足够小时,δ与绝对距离占 min{p2i,p2j} 的比率ρ=ε/min{p2i,p2j}成正比例关系.同理,可以得到当p2i为负值时δ与ε/min{|p2i|,|p2j|}成正比.因此,合并节点i和j后网络特征值的变化量δ与p2i和p2j间的相对距离ε/min{|p2i|,|p2j|}成正比.证毕.

从定理1可知,被合并节点对应的特性向量分量间的相对距离会直接影响合并后特征值的变化量,从而影响到粗粒化的效果.因此,在对特征向量分量进行分类时,应该以特征向量间的相对距离为依据.

4 ISCGR算法

与SCG和ISCG算法不同,ISCGR算法以相对距离作为分类依据,从而更合理地对节点分类,改善了网络的粗粒化效果,使所得粗粒化网络保持同步的能力增强.

在实际操作中,不能直接以ε/min{|p2i|,|p2j|}的大小作为设置分割点的判断标准,因为接近0的分量不满足分量间的距离ε远远小于min{|p2i|,|p2j|}这一条件,从而ε/min{|p2i|,|pε2j|} 会得出很大的距离比率.这里用来代替ε/min{|p2i|,|p2j|}作为设置分割点的标准,其中是λ2对应的所有特征向量分量绝对值的平均值.以代替ε/min{|p2i|,|p2j|}一方面避免了靠近0的分量算出过大的比率值,另一方面此指标包含了相对距离的影响因素,比采用绝对距离ε作为分类标准有更好的粗粒化效果.与ISCG算法相比,ISCGR算法在对特征向量分量进行分类时,越靠近0的分量分得越细化,聚类区间的长度越小,越远离0的分量分得越粗糙,聚类区间的长度越大.

ISCGR算法的具体的步骤为：1)把λ2对应的特征向量分量从小到大进行排序为计算两两相邻分量之间的距离εi(i=1,2,···,N-1)；2)用这个距离除以对应分量绝对值与所有分量绝对值平均值的和,即计算出3)根据ηi的大小设置分割点,将分量分裂成若干个聚类；4)把同一聚类中的节点合并为一个节点,利用(11)式得到粗粒化网络的Laplacian矩阵.

我们举例来说明ISCGR的操作步骤.假设λ2对应的特征向量P2的分量为 {p21,p22···p2N},要将其分为3类.首先把特征向量分量从小到大排列为相邻分量之间的间距为 {ε1,ε2···εN¯1},再通过公式计算出{η1,η2···ηN¯1},选出最大的两个相对距离记为ηj,ηk(j＜k),ηj对应着2j和2,j+1之间的相对距离,ηk是2k和2,k+1之间的相对距离.接着在这两个间距处设置分割点,把分量分裂成3类,分别为[21,2j],[2,j+1,2k]和 [2,k+1,2N] .最后将落入同一区间内的点合并为一个节点,利用(11)式更新合并后的Laplacian矩阵,就可以获得节点数为3的粗粒化网络.ISCGR算法对节点分类时考虑了特征向量分量间的相对距离,在两组分量间距相同的情况下,由于靠近0的分量相对距离大,所以会优先在靠近0的两个向量间设置分割点,把相对距离小的分量所对应的节点合并.由此可见,ISCGR算法在区分不同节点的相似程度时更加合理,保证了分在同一聚类的节点相似程度高于不同聚类的节点,因此粗粒化效果较ISCG算法更好.

与λ2类似,对于类型II网络,当考虑λN/λ2作为同步能力衡量指标时,先用ISCGR算法分别对λN和λ2对应的特征向量分量进行分类,将同时分裂在同一聚类的节点合并为一个节点,最后更新Laplacian矩阵,得到约简后的粗粒化网络.

5 仿真分析

本节将对一般连续时间的动态网络模型和27种现实世界网络进行数值仿真,分别应用ISCG方法和ISCGR方法对网络进行约简,比较两种算法对网络约简后同步能力的保持效果.

5.1 模型网络仿真

首先对连续时间动态网络模型进行数值仿真,这里选取三种典型的网络模型,即BA无标度网络[23]、ER随机网络[24]、和NW小世界网络[25,26],分别应用ISCG算法和ISCGR算法对其进行约简,分析比较没有考虑相对距离和基于相对距离的算法的粗粒化效果.

对于类型Ⅰ网络,保持同步能力不变只需要保持λ2不变即可.图3展示了分别采用ISCG算法和ISCGR算法对1000个节点的BA、ER、NW网络约简后λ2的变化情况.图3(a)是对1000个节点的BA网络的粗粒化过程仿真,原始网络的λ2为2.95.采用ISCG方法将网络约简至200,102,72个节点时,网络的2分别变为2.96,3.13,3.26,比原网络增加了0.34%,6.10%,10.51%；采用ISCGR方法将网络约简至200,102,72个节点时,网络的2为别是2.95,2.95,2.98,变化了0,0,1.02%.由此可见,采用ISCGR算法保持BA网络λ2的能力要优于ISCG算法.对于图3(b)1000个节点的ER随机网络和图3(c)1000个节点的NW小世界网络也能得出同样的结论,即在网络约简至相同的节点时,采用ISCGR算法λ2值的变化量要小于ISCG算法,这表明采用ISCGR算法保持类型Ⅰ网络同步能力的效果要优于ISCG算法.

对于类型II网络,同步能力用λN/λ2衡量.图4展示了分别采用ISCG算法和ISCGR算法对1000个节点的BA无标度网络、ER随机网络、NW小世界网络约简后保持λN/λ2情况的对比.从图中可以看出,采用这两种算法对网络进行约简后,λN/λ2都是随着约简后网络节点数的减少而减小,但在相同的情况下,用ISCGR算法得到粗粒化网络的λN/λ2变化更小,这表明了ISCGR算法保持类型Ⅱ网络的同步能力比ISCG算法更强.

图3 采用ISCG与ISCGR算法获得谱粗粒化网络对λ2的保持情况 (a)BA无标度网络；(b)ER随机网络；(c)NW小世界网络Fig.3.The maintaining of λ2obtained by using ISCG and ISCGR algorithms in coarse graining metwork：(a)BA network；(b)ER network；(c)NW network.

由此可见,对于不同类型的网络和不同的网络模型,采用ISCGR算法约简网络后保持同步的能力都要优于ISCG算法.

图4 采用ISCG算法与ISCGR算法获得谱粗粒化网络对 λN/λ2的保持情况 (a)BA无标度网络；(b)ER随机网络；(c)NW小世界网络Fig.4.The maintaining of λN/λ2obtained by using ISCG and ISCGR algorithms in coarse graining metwork：(a)BA network；(b)ER network；(c)NW network.

5.2 实际网络仿真

下面应用一些真实网络进行粗粒化过程仿真,进一步比较ISCG算法和ISCGR算法保持实际网络同步能力的效果.我们选取了协作、沟通、计算机、基础设施、社会、生物、词汇和电力等不同领域的网络[27-29].这些网络的规模从几百个到几千个节点不等,既包括稀疏网络又包括致密网络.表1给出了这些网络的基本数据及分别采用ISCG和ISCGR算法约简到原始节点数N的30%,20%,10%,2%时,网络Laplacian矩阵L的最小非零特征值λ2的变化量.图5展示了分别采用ISCG和ISCGR算法对这些实际网络约简后λ2的变化情况,其中图5(a)是从表1中挑选了一些随着节点数减少,λ2变化比较大的网络,图5(b)挑选了一些随着节点数减少,λ2变化不太明显的网络.图5可以看出约简后网络节点数越小,λ2的变化量越大,但是在相同的情况下,采用ISCGR方法约简网络后得到的λ2误差比采用ISCG方法要小.这表明对于实际网络,采用ISCGR方法比ISCG方法更能保持网络的同步能力.

对比表1中平均度和约简后λ2值变化量的关系发现,在约简程度相同时,网络的平均度越大,约简后网络λ2的变化量越小,保持同步的能力越强.这是因为网络的度越大,越接近于全局耦合网络,网络中外部连接相同或相似的节点越多,所以合并这些节点后对网络的特征值影响不大；而平均度越小,网络越稀疏,网络中外部连接相同或相似的节点越少,约简后网络的特征值变化量就越大.

图5 分别采用ISCG与ISCGR算法对实际网络进行粗粒化后保持 λ2情况的对比图Fig.5.The maintaining of λ2obtained by using ISCG and ISCGR algorithms for real-world networks in coarse graining network.

表1 分别采用ISCG和ISCGR算法对实际网络约简后 λ2的保持情况统计表Table 1.The Statistics table of maintaining λ2 obtained by using ISCG and ISCGR algorithms for some real-world networks in coarse graining network.

另外,表1的数据表明,ISCGR算法对互联网、生物、社交、合作等网络的同步保持能力更强一些,对电力、化学等网络的同步保持能力则弱一些.这是因为互联网[30]具有模块性,生物网络(如新陈代谢网络)中的拓扑结构是按等级组织起来的[31],社交网络与合作网络都具有“物以类聚,人以群分”的特点[32],即这类网络都具有高聚类的特性,网络中包含大量外部连接相同或相似的节点,所以ISCGR算法对这类网络的约简效果更好.相反,电力网络、化学网络中的节点大都是链状式排列,这类网络聚类特性不明显,外部连接相同的节点很少,所以采用SCG算法对网络约简到一定程度时很难保持网络的同步能力.

6 结 论

小规模网络中的一些性质不能直接应用在大规模网络中,所以复杂网络的约简成为了复杂网络理论中的重要课题.通过对网络约简,不仅可以把小规模网络中已经研究过的性质应用到大规模网络中,且简化了分析和计算.SCG算法是一种保持类型I和类型II复杂网络同步能力不变的约简方法.本文在原始谱粗粒化方法基础上提出了一种基于相对距离的谱粗粒化方法,使应用特征向量分量分类时比原算法更加合理,约简后的网络保持同步的能力得到增强,为大型网络的谱粗粒化提供了一种更加精准的约简算法；并通过对3种经典网络模型和27种实际网络的仿真分析表明,ISCGR算法比原算法粗粒化效果更好,进一步发现SCG算法对具有明显聚类结构的实际网络(互联网、生物、社交、合作等)比模糊聚类结构的实际网络(电力、化学等)更容易保持粗粒化网络的同步能力.

需要指出的是,SCG算法仅考虑了保持网络的同步能力不变对网络进行约简,所得粗粒化网络的其他动态特性或拓扑结构却不一定和原始网络相同或相近.SCG算法除了能保持同步能力之外,是否还能保持约简后网络的其他动力学特性,例如可控性、一致性、稳定性等? 是否还能保持约简后网络的小世界特性、度的幂律特性等网络拓扑性质? 谱粗粒化后的网络还有哪些性质保留了下来,其他性质发生了什么改变? 有没有一种约简方法可以兼并保持多种性质不变? 这些问题都是值得我们进一步研究和思考的.网络的约简问题还有太多的未知,这些未知也吸引着我们不断去探索.