徐 静,钱江海
(上海电力大学数理学院,上海 200090)
空间网络是指节点和边嵌入几何空间,并与某些实际的物理或几何效应(比如连边长度、构建代价、延时和能量耗散等)相联系着的一类复杂网络。这类网络包括航空网络[1-4]、道路网络[5-7]、Internet[8-9]、电力网络[10-12]和社会网络等等[13-16]。与那些只具有拓扑性质的系统不同[17-21],空间网络的结构、演化和动力学行为极大地受几何效应的制约[1-16]。这种制约带来的对拓扑结构的一个重要影响体现在网络度分布的标度律上。大量实证研究表明,随着网络空间效应的增强,其度分布会偏离幂律形态,而向较窄的指数型或高斯型分布过渡[1-7,22]。这一过渡在定性的层面是可以直观理解的。由于空间网络的连边和构建成本、时延以及能量耗散密切关联,为了降低这些代价,一个节点的连接目标自然局限于其附近的节点。这使得节点可连目标的数量大大减少,从而抑制了大度值节点的产生。
为了从理论上定量地解释这种空间效应对度分布的影响,科学家提出了各种空间网络模型。这些模型引入空间约束的方法大致可分为两种。一种方法是赋予每个节点一个有限边界的约束范围,只有当一方落在另一方的约束范围内时才可形成连边。这种约束模式及其相关的各类变种组成空间网络模型的一大类,即几何图模型[23-25]。另一种方法是引入距离衰减因子来修正连边的概率,这类模型包含Waxman模型[26]、空间小世界模型[27-28]、空间无标度模型和引力模型等等[4,29-30]。当节点在网络演化中互相竞争着获取连边时,距离衰减函数将导致每个节点在其自身周围产生一个有效影响范围。当目标落入该范围内时,节点与之连边的概率高于其他节点。因此,距离衰减函数本质上具有与有限边界相类似的约束范围的思想和内涵。事实上,有些基于距离衰减函数的空间网络模型正是通过引入有限边界的概念来帮助实现模型解析的[31]。但不管用哪种约束方式,通过调节约束强度,基本都能实现网络度分布从幂律到指数再到高斯型分布的转变。然而对于所有这些模型,有一种分布始终无法通过调节约束强度的方式产生,它就是在各类航空网络中出现的双段幂律度分布[1-3,32-34]。
双段幂律分布在其分布的上下端分别具有两种不同的标度律特征,其互补累积分布表达为
(1)
式(1)中的γ1与γ2表示分布指数,kc表示分段点。除了航空网络,这类分布还出现在语言网络[35]、森林火灾的规模[36]、个人收入[37]、网站点击次数[36]、龙卷风的规模[36]、恐怖袭击破坏程度以及人类行为等诸多复杂系统中[36,38]。实证得到的双段幂律分布一般具有γ1<γ2,这里称之为正双段幂律。而γ1>γ2的情形,即反双段幂律,似乎很少有见报道。关于双段幂律分布的生成模型并不多,目前主要有Reed提出的基于几何布朗运动的模型和Dorogovtsev、Mendes等人提出的语言网络模型[35,39]。Reed发现若一随机变量依据几何布朗运动变化,且其演化时间本身还服从指数分布,则演化终态得到的变量服从双段幂律分布[39]。该模型后来被Mitzenmacher详述并进一步衍生出其他变种[37,40-42]。但这些模型都不是网络模型。此外Reed自己指出模型结果对指数分布的演化时间假设是敏感的。Dorogovtsev和Mendes在BA模型的基础上引入了线性加速内部连边机制[35]。该模型可以产生双段幂律度分布,此模型的机制与空间约束无关,且只针对语言网络,其机制也未经实证检验。
既然空间约束可以导致度分布的过渡性转变,那么是否也存在某种约束模式可以产生双段幂律度分布?实际网络的双段幂律分布是否具有几何学上的起源?解决这些问题不仅有助于理解这类分布的形成原因,还能拓展对空间约束运作机理的认识。为此,本文提出了一种基于约束面积双值跳变的空间网络模型。它的核心机制可以简述为当一个节点的度演化到预设的kc时,该节点的约束面积便会从初始S1跳变为S2。换言之,模型要求约束强度本身具有动力学演化的性质。在现有的研究中,人们一般认为约束强度只和网络内在空间属性有关[43],通常是不随网络演化而变化的。然而一些实证研究明确给出了空间约束具备演化行为的证据[25,44-45]。考虑到基于静态空间约束的网络模型都未能产生双段幂律度分布,改用动态约束似乎是一种必然的尝试。事实上,我们认为这可能就是解决双段幂律问题的关键。
本文的内容安排如下:第1节给出约束面积双值跳变的空间网络模型的描述;第2节中对度分布进行解析计算,并证明该模型产生的度分布会在演化初期呈现单段幂律,而当有节点约束面积发生跳变后,网络逐渐演化为双段幂律;第3节对中国航空网络进行实证分析,验证模型机制对该真实网络的适用性,并探讨空间网络中约束面积减小的可能原因;最后一节对本文的工作予以总结。
本文的模型构建在[0,10]×[0,10]的二维方形区域中,网络通过不断添加节点和边而演化。我们遵循主流空间网络模型的共同思想内涵,赋予每个节点一个面积为S的约束范围,其形状规定为圆盘形。只有当一个节点的约束范围覆盖了另一个节点u的位置时(即其影响能够被u感知时),它才成为u连边时的候选目标。节点u仅在其候选目标集合内择优连一条边。每个节点的初始约束范围为S1,当其度增长到预设的kc时,其约束面积立刻跳变为S2。模型具体演化规则描述如下:
1)在每时间步u,一个具有初始约束范围S1的新节点(用u来标记)加入网络,它被随机分配一个位置(x,y)。
3)检查在此时间步中被u连接到的点,若连接后该节点的度达到预设的kc,令其约束面积等于S2。然后回到1)重复整个过程。
令模型所在空间区域的总面积为A,节点i被新加入节点u选中连接的概率为
(2)
(3)
在模型演化的第一阶段,所有节点度均未超过kc,此时网络中只存在唯一的约束面积S1。将其代入式(2),则模型直接退化为经典的BA模型[46],其互补累积度分布为P(k)=k-2。因此在这一阶段,度分布仅表现为单段幂律。
当有任意一个节点度值达到kc时,该节点的约束面积跳变到S2,网络从此时起就存在两种约束强度。式(3)因而对于不同度值范围内的节点具有两种不同的表达
(4)
(5)
将式(5)代入f(t)的定义式中并用积分运算替代求和,即
(6)
式(5)若为方程的真实解,式(6)关于f(t)的计算必须满足解的假设条件f(t)=Ct。从而可得到自洽方程
(7)
由式(7)可以确定C,从而完成对ki(t)的求解。然后采用文献[47]的方法可导出网络在第二阶段的互补累积度分布
(8)
为了验证以上理论分析,我们对模型进行了蒙特卡罗模拟。图1a~图1b是取kc=15、S1=30和S2=10的网络分别在演化初期(N=2 000)以及演化稳定后(N=15 000)的累积度分布图。可以看出度分布在两个阶段具有不同标度律表现,早期网络为单段幂律,而最终其演化为正双段幂律。当逐步增大S2时,模拟得到的稳态累积度分布如图1b~图1d所示。在S2 图1 kc=15时,不同模型参数下的累积度分布Fig.1 The cumulative degree distributions for different model parameters under kc=15 表1 双段幂律度分布指数Tab.1 Exponents of double power law degree distribution 为评估模型的有效性,验证空间网络中的双段幂律度分布确有可能起源于空间约束面积的变化,本节将针对中国航空网络(CAN)进行实证分析。需要指出,由于本模型只提供了一个概念性的框架,并不涉及任何实际系统的细节,因此寻求模型与真实观测量之间定量的匹配显然超出模型的能力。然而,如果双幂律度分布确实是由动态空间约束引起的,那么我们理应观测到这种动态约束造成的效应,即连边长度在不同时间段的变化趋势,而相应的度分布标度律演化特征原则上也应该能观测到。 文中所用到的航线数据全部来自《中国交通年鉴》(1988-2014),主要参照年鉴中统计表格《民航国内主要航线运输完成情况》中的航线。本文收集了跨度接近30年的CAN数据,可供我们建立每一年的网络拓扑并进行长时间的动态分析。 网络数据处理方法: 1)以通航的城市作为网络的节点。如果城市对之间实现了通航,那么两节点间就存在连边。最终某一城市的连边数目也即通航城市数,就是它在CAN中的度。 2)经停航线纳入统计。以经停2个城市的A-B-C-D航线为例,其中字母代表城市,该航线可分解为A-B、B-C、C-D 3个航段,每段按直飞航线处理,不区分航路方向。需要说明的是,我们主要采集了大陆区域的航线数据,不包括港澳台地区。 航线距离采集方法:任意一对空港城市间的航线距离都可以从Great Circle Mapper网站获得(http://www.gcmap.com/)。在该网站中,只要输入航线中机场对的三字码(IATA码)或四字码(ICAO码),就能采集到两地之间的飞行距离。 到目前为止,已有不少针对CAN的拓扑结构、流量动态、结构鲁棒性和服务质量的研究[2-4,48]。这些研究都证实了CAN具有典型的双段幂律度分布。但几乎所有这些研究的数据都采集于2002年以后,那时双段幂律度分布已经形成,从而无法辨识其可能存在的标度律演化特征。本研究追溯至CAN演化的早期,将观测时间提前至1988年,这比以前的研究都要早得多。 1988-2014期间CAN度分布的部分结果如图2所示。从图2a~图2b可看出在20世纪90年代早期至中期,度分布并不是双段幂律而是呈现明显的单段幂律特征。而在大约1999年左右,网络从单段幂律转换为双段幂律,并一直维持这一特征至今。图2c~图2d显示了2002年和2010年的CAN度分布,其为明显的双段幂律。经测量分段点kc≈22,这与先前的实证研究结果吻合[3,41]。因此CAN的双段幂律度分布并不是一开始就形成的,而是首先出现单一标度律并逐渐演化成双段幂律。这一结果完全吻合我们的模型对双段幂律分布形成过程的预言,但却无法由已有的其他相关模型解释。Reed和Mitzenmacher等人的模型[39-40]没有给出任何关于分布演化的动力学信息,且就其推导来看,他们的模型不具有实证中看到的幂律段接连出现的过程。Dorogovtsev和Mendes提出的语言网络模型[35]的度分布是含时的,但含时项并不处于指数项,即其标度律特征依然是稳定,这同样无法吻合目前的观测结果。此外Dorogovtsev和Mendes指出,在他们的模型中,度值大于分段点的节点集是恒定的。这与CAN情况不符合,因为我们明显可以看到不断有新节点城市加入到第二段幂律中去,比如深圳城市就是后续加入到第二段幂律段中的。 图2 中国航空网在不同年份的度分布Fig.2 Degree distributions of Chinese aviation network for different years 依照模型,度分布为正双段幂律的网络,其节点在跨越分段点kc后,约束面积将变小。为验证这一机制是否存在,本研究分析了每个空港城市在通航数经历kc≈22前后的约束范围的变化情况。我们将平均连边距离作为对约束空间范围的估计,如果分段点前后该距离在统计上有变小的趋势,这种趋势结合系统表现出的正双段幂律分布就完全吻合我们模型的预言。 对于每个空港城市,选取其演化过程中度值未超过kc的最后一年(记作Y,该年对不同节点是不同的),计算此时它所拥有的所有航线的距离平均值。然后选取该节点度值超过kc后的最后一年,即2014年,将该年里其航线中与Y共有的航线去除,计算剩下航线的平均距离。这样就能分别获得该节点在kc前后的约束空间范围的估计。然后通过比较各节点在到达kc前后连边平均距离的大小,便能看出约束空间范围变化的趋势。需要说明的是,CAN度分布的特征决定了大部分节点从没能演化到超过kc,因此这些节点无法纳入我们的统计和比较中。剔除这些节点后,还剩下29个空港城市。相比于网络的总节点数虽然少了很多,但若这些节点表现出统计上明显的趋势,则结果依然具有可信度。 将各个空港城市在达到kc前与后的平均连边距离作比,所得结果如图3所示。图中横坐标为赋予节点的标号,纵坐标为计算得到的比值。若比值大于1,则表明该节点的空间约束范围在到达kc后变小。图3里29个数据点中有22个数据点比值大于1,仅7个点小于1。这个结果表明,节点的约束空间在分段点前后具有统计上明显减小的趋势。而在该趋势下网络形成了正双段幂律度分布,这完全吻合了本模型的预言。 图3 29个空港城市在kc前后的平均连边距离的比值Fig.3 Ratios of average distance of edges of 29 airport cities before and after kc 如果空间约束是跳变减小的,那么是什么因素造成了这一结果?本文从预算和成本的角度出发,提出一种可能的解释。 我们假设空间网络中的连边可以分为两类,一类是与成本成正比的长程连边,另一类是地理上相邻节点之间无需代价即可形成的局部连边。这意味着连边长度有两种不同的尺度,即与整个空间的线性尺度成比例的长边L以及无需代价但只能局部连接的短边尺度l。事实上,一些推广的Kerlinberg模型就是采用类似的假设[49-51],而一些实证研究也证实了空间网络中连边长度具有的两种尺度的现象[52]。假定每个节点用来构造连边的预算B~kcL,在初始阶段,充足的预算允许它原则上将任何其他节点作为潜在的连接目标,因而这一时期获取的连边平均长度的规模约为L。随着边数的增多,该节点剩余的预算不断减少且每次平均减少大约L,该过程一直持续到预算降为B~L,然后只再需一个长程连边便可耗尽预算。而在随后的演化中,只有尺度为l的局部连接可用了。这便导致了本模型所描述的当节点达到某一临界度时,空间约束面积陡然缩减的效果,而分段点kc则刻画了节点构建连边的总预算。 本文提出了一种基于约束面积双值跳变的空间网络模型,该模型在演化的初期产生了单段幂律度分布,而后逐渐转变为双段幂律。通过改变跳变前后约束面积的相对大小,模型可以产生正反双段两种幂律度分布。这些结果通过解析方式得到,并经模拟仿真定量验证。为检验模型的适用性,本文对中国航空网络进行了实证分析。分析结果不仅表明度分布的演变过程恰如模型所预言,还直接观测到了模型所要求的约束面积减小的趋势。最后我们从预算和成本的角度出发,探讨了空间约束跳变减小的可能原因并指出分段点kc的物理意义。这些结果为双段幂律分布的起源提供了一种新的解释,同时拓展了对空间约束运作机制的认识。3 中国航空网实证
3.1 航空网络数据说明
3.2 中国航空网中非平稳的度分布演化特征
3.3 航线距离的演化特征
3.4 关于约束面积变化的原因
4 结论