基于复杂网络分析法的层级生态网络结构研究

王 戈 于 强 YANG Di 张启斌 岳德鹏 刘建华

(1.北京林业大学精准林业北京市重点实验室， 北京 100083； 2.弗罗里达大学地理系， 盖恩斯维尔FL 32611)

0 引言

高强度土地开发与土地利用方式的快速城市化进程已成为20世纪以来人类社会发展最显著特征[1]，资源环境过度利用与开发导致区域生态环境趋于恶化，严重影响区域景观格局与可持续发展[2]。尤其在西北半干旱区，伴随着人口增长，自然景观被人工景观代替，生境破碎，景观连通性变差，出现生物多样性加速下降等一系列问题[3]。基于先前层次生态网络构建方法研究,构建多层级空间生态网络是维持西北半干旱区生态安全的重要保障[4]。在层级生态网络中，低层级生态源地稳定依靠高层级生态源地，高层级生态源地对于维持层级生态网络稳定具有极其重要的意义[5]。高层级生态源地遭到破坏，易影响周围低层级生态源地，以至于影响低层级生态网络稳定，引发层级网络的级联失效，导致整个网络崩溃[6]。这种层级网络具有复杂的空间结构，层级之间具有复杂的影响关系，是典型的复杂系统，具有层次性、无序性与动态性等特征[7]。

复杂网络的研究方法已经被广泛应用到社交网络、专家网络、航空网络、交通网络等领域[8]。复杂网络的抽象方法可以将复杂系统简化为节点与边相连的集合，节点代表复杂系统的基本单元，边表示基本单元间的作用关系。结合景观生态学理论，可以理解为区域内生态廊道和生态节点两种景观格局要素所组成的复杂网络，结构、功能以及两者之间的联系一直是网络科学以及景观生态学的一个研究重点[9]。我国西北半干旱区具有生物多样性低、景观破碎、生态环境脆弱的特征，所构建的多层级生态网络具有多层生态网络，高层级生态网路是整个层级网络的“核心”，可以维持整个区域生态安全稳定性。低层级生态网络中存在的低层级生态源地，在小尺度上起到提升生态质量作用。层级生态网络在大尺度维持区域生态安全，并在小尺度上提升生态质量，形成完整的层级生态网络，从而全面保证区域生态安全。

本研究的复杂网络是具有层级性的复杂空间生态网络，并具有生态属性。层级复杂网络结构分析是网络优化的基础[10]，网络结构分析包括：网络节点度分析、节点与整个网络聚类系数分析、网络关联性分析、节点介数分析、网络连通性分析[11]。本文将复杂网络分析方法用于探究层级网络的空间拓扑结构[12]，选取典型西北半干旱城市包头市作为研究区，对包头市层级生态网络结构进行分析，并选取适合空间层级生态网络的鲁棒性评价指标对层级生态网络的结构鲁棒性进行分析。

1 材料与方法

1.1 研究区概况

包头市位于内蒙古自治区西部，面积27 768 km2，东经109°13′～111°26′，北纬40°13′～42°44′(图1)。包头市深处内陆，为典型的温带干旱、半干旱大陆性气候，冬季寒冷干燥、夏季炎热多雨[13]，年平均气温2.0～7.7℃，年均降水量175～400 mm，年均蒸发量为2 100～2 700 mm。可利用地表水总量为9×108m3，地下水补给量为8.6×109m3。黄河流经包头境内214 km，水面宽130～458 m，最大流量6 400 m3/s，年平均径流量为260×109m3，是包头市主要用水来源[14]。研究区海拔976～2 317 m，地势中间高南北低，北部丘陵、中部山地、南部平原分别占土地面积的14.49%、75.51%和10%[15]。干旱的气候条件与起伏的地貌特征使得包头市生态环境较为脆弱，面临较高的水土流失与土地沙漠化、荒漠化风险。近年来，包头市城市建设规模扩张迅速，房地产开发项目、工业园区等建设项目不断涌现，导致草地、林地、湿地等生态用地遭到占用与破坏，生态风险有所升高[11]。

1.2 数据来源与处理

选取包头市夏季且少云的Landsat OLI影像(2018年)为研究素材，空间分辨率30 m DEM等作为辅助数据。利用ENVI 软件对影像进行波段合成、图像增强和几何校正处理，选择最大似然监督分类法对遥感影像进行目视解译，提取包头市的景观类型信息，使用ArcMap软件进行细碎斑块处理，运用叠加分析工具进行空间数据分析，最终ArcInfo环境下完成拓扑和改错处理。自然保护区分布数据来自《包头市土地利用总体规划》；地面高程数据(DEM)来自地理空间数据云；根据研究需要，从遥感影像中提取的景观类型通过实地验证确保其精度。

1.3 层级生态网络构建

1.3.1层级生态源地节点提取

根据研究区景观格局分布数据，提取林地、草地、水体生态用地，将生态用地划分为绿色和蓝色生态用地。利用ArcGIS的分区统计工具，分别统计每个生态斑块所对应的植被指数(NDVI)和改进的归一化差异水体指数(MNDWI)值。计算所有生态斑块的能量因子Qj，并结合景观生态功能对生态用地进行层级划分[16]。能量因子Qj计算公式为

Qj=AjNjr

(1)

式中Aj——第j块生态用地斑块面积，km2

Njr为第j块生态用地斑块的第r个归一化指数，本文选择两种归一化指数描述生态用地，故r分别取1和2，Nj1为第j块生态用地的NDVI平均值，Nj2为第j块生态用地的MNDWI平均值。

1.3.2层级生态廊道提取

将生态源地、源地能量等级、源地间距离、基面阻力系数4方面因素考虑到模型中，得到修正层级生态廊道提取模型[17]，公式为

(2)

式中VMCRQj——最小累积阻力面值

fmin——一个土地单元的累积阻力最小值

Dij——从生态源地j到土地单元i空间距离

Ri——用地单元i运动过程的阻力系数[18]

Qj越大代表生态源地斑块的生态能量越大[19]。

根据包头市的实际情况，考虑从地形坡度、植被覆盖、水文分布、土地覆盖共4方面建立生态阻力的评价体系。进行叠加栅格计算得到生态阻力基面的综合评价结果。利用ArcGIS软件中的cost-distance模块完成层次生态源地的累积阻力面值的计算[20]，最终得到每层生态网络对应的最小生态累积阻力面值。利用ArcGIS软件中的cost-path模块完成层级生态廊道提取。

1.4 复杂网络模型

利用ArcGIS软件中空间关系建模等功能，把分层生态源地分别抽象成Ni个节点(i=1,2,3)，Ei条边构成包头市层级生态网络。利用空间统计中的工具把网络转换成表，再将表转换为Pajek能识别的.net格式的网络数据。

1.4.1网络基本静态统计特征

1.4.1.1度及度分布

生态节点度是与该生态节点相连廊道数量，节点度越大表明该生态节点重要性较高。网络平均度是网络中每个生态节点度平均值[21]，计算公式为

(3)

式中k——网络平均度

N——节点总数

k′i——生态节点度

在复杂生态网络中，节点的度在统计学上服从一定分布函数，生态节点度分布p(k)代表在一个生态网络中，度为k的生态节点占比，度分布函数p(k)是节点度为k的节点被抽到的概率。

1.4.1.2平均路径长度

对网络中任意两个节点最短路径边数求和后，对其取平均值，得到复杂生态网络平均路径长度[22]，计算公式为

(4)

1.4.1.3聚类系数

(5)

平均聚类系数C为所有生态节点聚类系数平均值，公式为

(6)

复杂生态网络平均聚类系数为0～1，C=0表明所有节点没有任何边连接；C=1表明网络中任意两节点均直接相连；C值越大，表明生态网络中节点联系越紧密，聚集系数越大的网络小世界特性越强，反之越弱。

1.4.2生态网络的关联性

1.4.2.1基于Pearson相关系数的度-度相关性

实际生态网络度与度之间存在相关性，网络中度大的节点连接概率高，表明该生态网络是度-度正相关生态网络，即同配生态网络[24]，反之为异配生态网络。一个生态网络度的Pearson相关系数r计算公式为

(7)

式中ki、kj——连接生态廊道eij两生态节点vi和vj的度

M——网络中生态廊道总数

E——所有生态廊道集合

1.4.2.2聚类系数分布及聚-度相关性

聚类系数分布函数P(C)代表在复杂网络中任意选择一个节点，其聚类系数为C的概率。一个生态网络的聚类系数和度之间也存在一定内在联系，即聚-度相关性。聚-度相关性计算中主要是局部聚类系数C(k)与k之间的关系。

1.4.3生态网络的节点介数

生态网络节点介数是生态网络中任意的两个生态节点之间最短生态廊道路径所通过的生态节点数量，生态节点介数可以用来反映生态网络中生态节点的重要性[25]。生态节点介数Bi计算公式为

(8)

式中njl——生态节点vj和vl之间的最短生态廊道数量

njl(i)——生态节点vj和vl之间最短生态廊道路径经过生态节点vi数量

1.4.4生态网络的连通性

1.4.4.1生态网络核数

生态网络核数是指不断去掉度小于k生态节点和其相互连接的生态廊道，最终剩余的联通子图中生态节点数目[26]。如果一个生态节点属于k核，但是不属于k+1核，这个生态节点的核数为k,该生态网络的核数同样为k。通过k核统计分析，复杂网络逐渐趋于核心的区域，越位于中心的核，连通性越强[27]。

1.4.4.2节点的连通度

网络节点的连通度主要反映了生态网络连通程度。连通复杂网络G的连通度k(G)定义为

k(G)=min{|S|:S⊆V(G)，ω(G-S)≥2}

(9)

式中V——生态网络G的生态节点组合

S——V的真子集

ω(G-S)——从生态网络G中删除生态节点集S后得到的子图的联通分支数

生态节点联通度就是指使G不联通或者平凡图(只有一个生态节点没有边的生态网络)所必须删除掉最少生态节点数量。对于不联通生态网络，定义k(G)=0;若G为N个生态节点完全生态网络[28]，则有k(G)=N-1。

1.5 分层生态网络鲁棒性

生态网络是维持区域生态环境稳定的重要保障，完整的空间结构能保证其发挥出正常功能。生态网络结构鲁棒性是用来衡量生态网络在遭受到外界干扰破坏时，维持其正常结构和功能的能力，即潜在生态网络抵抗能力，与之相对应的就是潜在生态网络恢复能力[29]。在其空间结构遭到破坏后，潜在生态网络能够恢复的能力是恢复鲁棒性。结构鲁棒性公式为

(10)

式中N′——初始潜在生态网络的节点数

Nr——从生态网络中去除的节点数

c——当节点被去除后生态网络中最大连通子网络中的生态节点数量

针对生态节点和生态廊道，恢复鲁棒性计算公式分别为

(11)

(12)

式中D——节点恢复鲁棒性指标

E——边恢复鲁棒性指标

Nd——通过某种策略恢复节点数

M′——初始网络中边的数量

Mr——从网络中去除的边数

Me——通过某种策略恢复的边数量

通过构建实际生态网络表示节点关系的邻接矩阵来评价生态网络的结构鲁棒性，制作无向无权生态网络拓扑图，并对网络进行随机攻击和恶意攻击[30]。随机攻击即从网络中随机去除若干个节点，恶意攻击即从潜在生态网络中，同时去除度最大的Nr个节点及其对应边。

2 结果与分析

2.1 层级生态网络构建

研究区域层级生态网络如图2所示，第1层生态网络由8个潜在生态源地、8条潜在生态廊道和7个生态节点构成。第2层生态网络由31个潜在生态源地、35条潜在生态廊道和28个生态节点组成。第3层生态网络由123个潜在生态源地、151条潜在生态廊道和47个生态节点组成。通过对包头市全域内生态网络进行分层研究，构建1、2、3层生态网络，有利于保障包头市生态环境安全。在市域尺度上构成了分层的点-线-面相互交织的层级空间生态网络。

图2 潜在生态网络Fig.2 Potential ecological networks

2.2 层级生态网络拓扑结构分析

2.2.1第1层生态网络

所提取的第1层生态网络共有15个节点，网络的平均度为2，第1层生态网络的度分布如图3所示。度最大为3，最小为1，度为1和3的节点分别有2个，度为2的生态节点数量最多，共有11个。

图3 第1层生态网络度分布Fig.3 Ecological network distribution of the first layer

在包头市第1层级生态网络中，生态节点与生态节点之间具有相关性。同样，度与度之间也存在一定相关性。计算得第1层潜在生态网络的度-度相关系数为-0.339 3，该生态网络负相关，即异配网络。在第1层生态网络中度大的节点倾向与度小节点相连。

通过计算该网络生态节点介数并绘制散点图(图4)，有两个生态节点的介数为0，编号为1和6。介数最高为101，编号为3，分布在红花敖包保护区附近。介数小于100的节点共有12个，编号10、11、12的生态节点并非生态源地且生态能量较低，但在网络中具有较高的介数。

图4 第1层生态网络介数分布Fig.4 Ecological network mediator distribution of the first layer

计算得第1层生态网络各节点核数均为2，得到该网络核数为2。对该网络节点连通性进行计算，得该网络所有节点连通性均为14，表明该网络节点连通性较低，网络结构简单，结构稳定性较差。

2.2.2第2层生态网络

第2层潜在生态网络共有60个节点，网络平均度为2.166 7，网络节点度分布如图5所示。该生态网络完全连通，不存在度为0的节点。生态节点的度最大值为4，最小值为1，数量最多的生态节点度为2，共有40个。

图5 第2层生态网络度分布Fig.5 Ecological network distribution of the second layer

包头市第2层生态网络度-度相关系数为-0.180，较第1层生态网络度-度相关系数降低，异配特性减弱，第2层生态网络节点分布较第1层潜在生态网络更加均匀。

计算第2层生态网络节点介数并绘制散点图(图6)，有7个生态节点介数为0，介数小于1 000的节点共有45个，介数大于1 000的共有8个，其中编号14、17、18、19、27、28为生态源地节点，节点的度较高，在网络起到能量循环与物质交换枢纽的作用。编号40、54的生态网络节点大于1 000并且两端连接均为介数大于1 000的生态源地节点，是第2层生态网络中最薄弱且最重要的生态节点，分布在下湿壕乡和都荣敖包苏木附近。

图6 第2层生态网络介数分布Fig.6 Ecological network mediator distribution of the second layer

计算得第2层生态网络各节点核数均为2，得到该网络核数为2。对该网络节点连通性进行计算，得该网络所有节点连通性均为59。由于在第2层网络中生态节点数量增加，导致生态网络节点连接性变好，整个网络连通性变好。

2.2.3第3层生态网络

图9 第3层生态网络的聚-度相关性散点和双对数分布Fig.9 Poly-degree correlation scatter and double logarithmic distribution of the third layer ecological network

第3层生态网络共有180个节点，网络平均度为2.549 5，第3层生态网络度分布如图7所示。度最大为12，生态节点编号为120，分布在九峰山自然保护区附近。度为2的生态节点数量最多，共有98个。度为1的生态节点有17个，第3层潜在生态网络整体散点分布既不是典型的幂律分布，也不是 Poisson分布，但是幂律分布特征更为明显，Poisson 分布特征不明显，表明该潜在生态网络无标度性特征要强于均匀性特征。

图8 第3层生态网络聚类系数分布Fig.8 Ecological network clustering coefficient of the third layer

第3层潜在生态网络的平均聚类系数为0.111 1，该网络的小世界特性不明显。该潜在生态网络节点的聚类系数分布如图8所示，聚类系数为0节点有142个，占该网络节点总数78%，这些节点与其他节点不具有集群特点。有10个节点的聚类系数为1，这些节点有明显集群特征。表明第3层潜在生态网络有明显非均匀性，生态网络结构不稳定性较大。

所提取的第3层潜在生态网络的度-度相关性为0.162 4，该生态网络是正相关，即同配网络。表明网络中度大的节点趋向和度大的节点相连。第3层生态网络的聚-度相关性散点和双对数分布如图9所示，随着生态节点度增加，节点聚类系数先增加后减少，该网络度小节点聚集特性较为明显，度大节点分布较分散。

图10 第3层生态网络介数分布Fig.10 Ecological network mediator distribution of the third layer

计算第3层生态网络节点介数并绘制散点图(图10)，有28个生态节点的介数为0，介数小于8 000的节点共有159个，介数大于8 000的共有5个，编号为62、70、79、110的节点距离较近，位于塔尔浑河附近，形成了一个较大的生态节点组团。生态节点136位于石宝镇附近，为非生态源地节点，由于两端连接介数较大的生态源地节点，是第3层生态网络中薄弱但十分重要的生态节点。

所提取的第3层生态网络中核数为2的生态节点最多，达176个。核数为4的生态节点有5个，核数最高值为12，共有2个节点。对该网络的节点连通性进行计算，得该网络所有节点的连通性均为180。相比第1层和第2层生态网络，第3层生态网络节点连接性变好，使整个网络的连通性变好。

2.3 层级生态网络鲁棒性分析

图11 第1层生态网络结构鲁棒性Fig.11 The first layer ecological network structure robustness

第1层生态网络初始连接鲁棒性为1(图11)，任意打击该层级网络中2个节点，网络的连接鲁棒性仍然为1，随着节点打击规模增加，连接能力下降迅速(图11a)。除去第9个生态节点时，网络的连接鲁棒性低于0.1。在恶意攻击下随着节点打击规模增加，连接鲁棒性下降迅速，除去第8个节点后连接鲁棒性低于0.1，生态网络连接能力极差，可见恶意攻击对第1层生态网络的连通能力破坏十分明显。网络节点恢复鲁棒性如图11b所示，在恶意攻击和随机攻击下，当去除两个节点时网络完全可以恢复。随着破坏节点数目增加，两种攻击下节点恢复鲁棒性均呈下降趋势，当攻击除去超过15个节点时，恶意攻击和随机攻击下丢失节点得不到恢复。对边恢复鲁棒性进行分析, 恶意攻击和随即攻击下的边恢复鲁棒性均呈现下降趋势。随机攻击下边恢复鲁棒性曲线为凸曲线，当除去边数大于6时，随机攻击下的边恢复鲁棒性优于恶意攻击。在恶意攻击下，当除去6条边时网络结构可以恢复。恶意攻击下随着边去除数目增加，边恢复鲁棒性下降呈现一定的“涌现”现象，对网络结构破坏非常显著。

第2层生态网络鲁棒性如图12所示，其初始连接鲁棒性为1，在随机攻击下，对2个节点进行随机打击，网络连接鲁棒性仍然为1，随节点打击规模增加，除去2～25个节点时，网络连接能力下降迅速。在恶意打击下，在除去第2个重要节点后网络连接鲁棒性下降显著，当除去第15个节点时，网络的连接鲁棒性低于0.1，生态网络的连接能力极差。在恶意打击和随机打击下，当除去的节点为17个时，网络可以完全恢复。去除17～60个节点时，恶意攻击和随机攻击时网络节点恢复鲁棒性均呈下降趋势，当除去55个节点后恶意攻击和随机攻击下丢失节点得不到恢复。网络的边恢复鲁棒性如图12c所示，除去5条边时，网络可以完全恢复。随着破坏边的数目增加，随机攻击下的边恢复鲁棒性均呈下降趋势，当攻击除去超过57条边时，随机攻击下去掉的边不能恢复。在恶意攻击下，当除去的边超过28条，边恢复鲁棒性下降更加迅速。

图12 第2层生态网络结构鲁棒性Fig.12 The second layer ecological network structure robustness

图13 第3层生态网络结构鲁棒性Fig.13 The third layer ecological network structure robustness

第3层生态网络鲁棒性如图13所示，其初始连接鲁棒性为1，在随机攻击和恶意攻击下，随着节点打击规模增加，连接能力呈下降趋势，在恶意攻击除去第10个节点和随机攻击除去第20个节点时网络的连接鲁棒性出现一定“涌现”现象。在恶意攻击除去第62个节点和随机攻击除去第80个节点时，网络连接鲁棒性低于0.1，网络连接性极差。网络节点恢复鲁棒性如图12b所示，在恶意攻击和随机攻击下，当去除20个节点时，网络完全可以恢复。破坏节点数目增加，两种攻击下节点恢复鲁棒性均呈下降趋势，当攻击除去超过170个节点时，恶意攻击和随机攻击下丢失的节点得不到恢复。对边恢复鲁棒性进行分析(图13c), 恶意攻击和随机攻击下的边恢复鲁棒性均呈现下降趋势。除去13条边时，网结可以完全恢复。随机攻击下的边恢复鲁棒性曲线为凸曲线，随机攻击下的边恢复鲁棒性优于恶意攻击。

3 结论

(1)根据包头市实际情况，在市域尺度上构建了点-线-面相互交织的层级生态网络。第1层由8个潜在生态源地、8条潜在生态廊道和7个节点构成，第2层由31个潜在生态源地、35条潜在生态廊道和28个生态节点组成，第3层由123个潜在生态源地、151条潜在生态廊道和47个生态节点组成。

(2)基于复杂网络中拓扑结构分析指标，分析提取第1、2、3层生态网络的拓扑性质。第1层网络平均度为2，度-度相关系数为-0.339 3，该层网络核数为2。第2层生态网络平均度为2.166 7，度-度相关性为-0.180，该网络核数为2，网络连通性为59。第3层生态网络平均度为2.549 5，网络聚类系数为0.111 1，度-度相关性为0.162 4，核数为4，网络连通性为180。第1层网络连通性最低，结构简单但重要性最高，影响整个区域内层级生态网络稳定。第3层生态网络连通性最高且结构复杂，可在小尺度上维护生境稳定。

(3)对层级生态网络进行结构鲁棒性分析。第1、2、3层生态网络的初始连接鲁棒性均为1，连接鲁棒性分别在除去9、15、20个生态节点时，连接鲁棒性明显下降。节点恢复鲁棒性分别在除去2、17、20个生态节点时，网络可以完全恢复。边恢复鲁棒性分别在除去6、5、13条边时，网络可以完全恢复。第2、3层生态网络中，低等级源地比例较高。由于高层级源地节点稳定性较高，可维持区域大尺度生态网络稳定。低层级源地可以提升小尺度生态安全，但对增强网络抗打击破坏能力与恢复能力效果不明显。