荒漠绿洲交错区生态空间网络粒度效应研究

张启斌 刘冬悦 李 倩 于 强 宋宏利 李世冉

(1.河北工程大学地球科学与工程学院，邯郸 056038；2.河北工程大学园林与生态工程学院，邯郸 056038；3.北京林业大学精准林业北京市重点实验室，北京 100083；4.国家林业和草原局林产工业规划设计院，北京 100714)

0 引言

生态空间网络是景观生态学中重要的概念与方法[1]，是耦合格局、过程与功能的重要途经[2]。近年来，众多学者围绕生态空间网络展开研究，以解决城市化进程中引起的人类扰动增强、景观斑块破碎、生物多样性降低等问题[3-4]。在这一过程中，生态空间网络的相关研究方法不断丰富，地理信息系统、图论、复杂网络等多种分析方法均得到了广泛的应用。

袁少雄等[5]基于网络环通度指数、线点数和连接度对广东省自然保护区生态空间网络进行了评价并提出修复建议；PASCUAL-HORTAL等[6]提出了基于图论的整体连接度指数(IIC)、可能性连接指数(PC)以及等效连接性指数(EC)等，可以更好地对生态空间网络连接性与斑块、廊道重要性等指数进行分析。在生态空间网络的各种分析方法中，基于复杂网络的分析方法由于基本涵盖了图论分析方法中的各项指数，可以更全面地剖析复杂生态空间网络的各项拓扑特征[7]，正在成为领域内的热点方法。YU等[8]通过复杂网络的相关理论研究了磴口县生态空间网络的拓扑特征与鲁棒性，并对生态节点的空间布局进行了优化；牛腾等[9]通过多项复杂网络相关指数，研究了乌兰布和沙漠东北缘生态空间网络的拓扑特征。

在当前的各类研究中，相关学者往往在各自研究区域内采用各自不同的空间尺度展开，不同粒度下生态空间网络空间结构与拓扑结构的差异未见报道。尺度问题是景观生态学的核心问题之一[10]，生态空间网络在空间上表现为生态源地斑块、生态廊道与基质景观的镶嵌体[11-12]，具有典型的空间异质性和尺度依赖性[13]，进行生态空间网络的粒度效应研究对理解生态空间网络空间、拓扑结构与空间尺度的关系具有重要意义[14-15]。

本文以乌兰布和沙漠东北缘的磴口县为研究区，通过最小累积阻力模型构建县域生态空间网络，并分析其空间与拓扑特征的粒度效应，以期为景观生态学相关研究中适宜尺度的选取以及区域生态保护提供科学依据。

1 材料与方法

1.1 研究区概况

选取内蒙古自治区巴彦淖尔市磴口县为研究区，位于内蒙古自治区西南部(106°10′～107°10′E,40°10′～40°57′N)，如图1所示。磴口县县域西侧为狼山山脉，东部为黄河，北侧为河套平原，南侧为乌兰布和沙漠。县域内部从东南向西北逐渐由绿洲向沙漠过渡，属典型的荒漠绿洲交错区，生态区位极为关键。县域大陆性气候特征显著，降水稀少(144.5 mm)，蒸发强烈(2 397.6 mm)，年均风速3 m/s，风蚀强烈，土地沙化风险较高[16]。

图1 研究区边界Fig.1 Boundary of study area

1.2 数据来源与预处理

本研究原始数据包括磴口县土地利用、地下水埋深、植被覆盖、高程与坡度。其中，土地利用数据来自中国林业科学院沙漠林业研究中心，原始地类共有26种，经合并后分为耕地、林地、草地、水体、建设用地、未利用地，基于土地利用数据，提取了水网、建设用地以及路网密度。地下水埋深数据通过磴口县40个地下水监测点的埋深经空间插值得到。植被覆盖数据通过欧洲航天局下载的Sentinel-2影像提取得到。高程数据为ALOS卫星的DEM数据，空间粒度12.5 m，并基于此数据提取了坡度数据。

基于上述原始数据产生不同空间粒度的数据是本研究开展的前提。不同空间粒度数据的生产需要确定粒度的变化范围以及该范围内粒度的间隔。当前生态空间网络的研究粒度基本分布在10～1 000 m。因此，本研究也将粒度变化范围控制在10～1 000 m之间。本研究借鉴以往研究的粒度划分方法[17-19]，将各类数据的空间粒度统一为10 m×10 m，基于该初始粒度数据，在10～100 m每隔10 m，100～1 000 m每隔100 m进行重采样，得到研究区不同空间粒度数据。

由于本研究基础数据既包括离散数据(如土地利用)又包括连续数据(如DEM)，因此选取最邻近像元法进行重采样。重采样后，各粒度下数据采用相同的技术路径进行生态源地提取、生态网络构建以及各项指数的分析，这保证了不同粒度下仅有数据空间粒度的差别而其他条件保持一致。

1.3 生态源地提取

生态源地指提供重要生态功能，对维持生态过程连续性、防止生态退化发挥重要作用的景观斑块[20]。磴口县位于乌兰布和沙漠与河套平原交界地带，土地沙化面积较大，生态环境脆弱敏感，林地、草地、水体在区域中发挥着涵养水源、调节气候、防风固沙等重要的生态服务功能，因此本文选取林地、水体以及植被覆盖度大于0.5的草地作为生态源地。本研究所提取的生态源地初始粒度为10 m×10 m，其他各粒度生态源地数据通过初始数据重采样得到。

1.4 生态空间网络构建

采用最小累积阻力模型(Minimal cumulative resistance model)进行研究区生态空间网络的构建。该模型考虑源、距离和景观介面阻力三方面因素，基本公式可表示为[21-22]

(1)

式中fmin——未知正函数，代表景观中任何一个源地的最小累积阻力与其到其他所有源地的几何距离及景观基面特征的正相关关系

Dij——源地i和源地j之间的距离

Ri——景观i的阻力

M——最小累积阻力

景观阻力指景观对生态流流动速率的影响[23]。不同性质的景观元素产生阻力不同，这种阻力的差异将对物种的运动、物质能量的流动产生不同的影响[24]。磴口县位于干旱的乌兰布和沙漠边缘，地表水与地下水分布、土地利用格局、地形、植被覆盖度、人类活动等因素共同影响着该区域景观阻力。本文选取地下水埋深、水网密度、植被覆盖度、坡度、土地利用类型、建设用地密度、路网密度共7个因子对磴口县生态阻力进行评价，各因子采用等间距重分类为5类，分级标准及阻力评分如表1所示。各阻力因子分级结果如图2所示。

图2 阻力因子Fig.2 Resistance factors

表1 阻力因子分级Tab.1 Classification of resistance factors

各阻力因子数据初始空间粒度为10 m×10 m，各因子叠加后的阻力面具有同样的空间粒度，为了得到不同粒度的阻力面数据，本研究首先将各因子重采样，通过叠加分析得到不同粒度的阻力面数据。

在实际计算中，采用成本路径(Cost path)模型代替MCR模型中未知正函数的作用。该模型的输入包括源、目的地和阻力栅格数据。运行中，该模型利用累积阻力面为每个像元分配到最近像元的累积成本，通过对该像元临近像元的遍历，识别出成本最小的临近像元，不断迭代该过程便可得到最小成本路径，该路径就是生态源地间阻力最小的生态廊道，众多生态廊道互相交错，形成生态空间网络[25]。在不同空间粒度下分别运行该模型，得到不同粒度下的生态空间网络。

1.5 生态空间网络拓扑结构分析

生态空间网络由提供重要生态功能的生态源地及连接着这些生态源地的生态廊道构成[26]。生态源地间的这种连接关系可被看作复杂网络的一种[27]。利用复杂网络分析中的多种指数研究生态空间网络的拓扑结构，可揭示其连接特征、识别关键源地与关键廊道，为生态空间网络的构建与优化提供重要借鉴。然而当前往往基于某个固定的空间粒度进行生态空间网络拓扑结构的分析，分析过程中空间粒度影响程度仍然有待研究。因此，本文选取度、度分布、网络直径、平均路径长度、平均聚类系数、连通度、平均介数、模块化、度-度相关性、聚-度相关性等指数，计算各指数在不同空间粒度下的取值，分析各指数的粒度效应。各指数的计算方法详见文献[9,28]。

2 结果分析

2.1 生态源地粒度效应

磴口县生态源地提取结果如图3所示。由图3可知，生态源地主要分布于县域中东部，西部山区仅有少量河流。所有的生态源地中，面积较大的通常为水体，地类为林地与草地的生态源地空间形态比较破碎。

对磴口县生态源地进行重采样，得到不同空间粒度下的生态源地数据，其中粒度50、100、500、1 000 m的数据如图4所示。由图4可知，伴随像元粒度增大，生态源地呈现出越来越强的破碎性，众多带状生态源地消失，而团块状的生态源地得到保留。

图4 不同空间粒度下生态源地Fig.4 Ecological sources under different grain-sizes

统计不同空间粒度下生态源地的面积如图5所示。由图5可知，在粒度10～500 m，县域生态源地的面积基本保持稳定，然而在粒度600～1 000 m，生态源地面积出现了大幅波动，其总面积从544.75 km2下降到522.72 km2，之后又升高到561.33 km2，最后在粒度1 000 m回到540 km2附近。生态源地总面积的波动，将对生态空间网络整体的空间及拓扑结构产生明显的影响。

图5 不同空间粒度下生态源地面积Fig.5 Ecological source area under different grain-sizes

2.2 生态阻力粒度效应

将地下水埋深、水网密度、植被覆盖度、坡度、土地利用类型、建设用地密度、路网密度共7个因子的阻力进行栅格相加运算，得到县域阻力面数据。利用ArcGIS 10.2中的cost distance工具生成累积阻力面，结果如图6所示。

图6 累积阻力Fig.6 Accumulated resistance

由图6可知，磴口县生态阻力最小值为1，均分布在生态源地范围内，最大值为898 978，主要分布在狼山区域。磴口县粒度50、100、500、1 000 m下累积阻力如图7所示。由图7可知，磴口县累积阻力的最低与最高值并没有因空间粒度的变化产生变化，并且各粒度下累积阻力的空间变化趋势相同。空间粒度变化对县域累积阻力的影响仅仅使得较大粒度下产生了较明显的颗粒感。

图7 不同粒度下累积阻力Fig.7 Cumulative resistance under different grain-sizes

2.3 生态空间网络空间结构粒度效应

利用成本路径模型提取研究区潜在生态空间网络如图8所示。由图8可知，磴口县潜在生态空间网络在县域中部与东部较为密集，西部区域由于生态源地较为稀疏，生态空间网络密度较小。县域内部存在一些空间距离较近的生态源地之间没有廊道直接相连，例如县域西南方的部分源地。这主要是由于两个源地之间的累积阻力大于经过其他源地间接将二者相连的累积阻力的原因。

磴口县空间粒度50、100、500、1 000 m下潜在生态空间网络如图9所示。由图9可知，伴随空间粒度的升高，潜在生态空间网络呈现越发破碎混乱的状态，这主要是由于空间粒度的升高使得生态源地中的河流、沟渠等带状部分变得破碎，破碎的源地间互相连接，便产生了较为混乱的网络。伴随空间粒度的继续提高，原本错综复杂的网络又开始变得稀疏，这主要是由于空间粒度的继续升高使得众多面积较小的生态源地消失，而面积更大的源地得以加强，此时的网络仅需保持少数主要生态源地之间的互相连接，故变得较为稀疏。

图9 不同粒度下生态空间网络Fig.9 Ecological spatial network in different grain-sizes

不同空间粒度下生态空间网络中廊道总长度如图10所示。由图10可知，伴随空间粒度的升高，网络总长度逐渐增加，在空间粒度60 m时达到峰值，最终在空间粒度1 000 m下降为1 165.92 km。网络总长度的这一变化趋势与网络空间结构的变化相对应。

图10 生态空间网络总长度Fig.10 Length of ecological spatial network variation curve

2.4 生态空间网络拓扑结构粒度效应

度-度相关性、聚-度相关性度分布幂函数指数、模块化呈现上升趋势，如图11～14所示。由图11、12可知，网络的度-度相关性和聚-度相关性具有相似的变化趋势，在空间粒度10～300 m的范围内，其增加趋势较为缓和，尤其是度-度相关性，甚至出现了下降趋势，然而在空间粒度到达300 m之后，上升趋势突然陡峭，之后在空间粒度600～800 m范围内出现了轻微下降，在空间粒度900～1 000 m范围内又开始上升。度-度相关性的变化趋势表明，随着空间粒度的升高，磴口县生态空间网络的同配性不断升高，度大的节点更倾向于与度大的节点连接，网络不均匀性升高。由图12可知，随着空间粒度的增大，度大的节点更倾向于具有更高的聚类系数，网络中正在形成以度大节点为核心的模块，且出现了层次性。由图13可知，随着空间粒度的增大，网络模块化也增高，这与聚-度相关性的分析结果吻合。磴口县生态空间网络度分布幂函数的指数同样随着空间粒度的升高而升高，如图14所示，幂函数指数的升高表明网络度分布的变化曲线正变得越来越陡峭，度分布更加集中，网络不均匀性增强。

图11 度-度相关性变化曲线Fig.11 Degree-degree correlation variation curve

图12 聚-度相关性变化Fig.12 Clustering-degree correlation variation curve

图13 模块化变化Fig.13 Modularity variation curve

图14 度分布幂函数指数变化Fig.14 Exponent of power function variation curve

由图15、16可知，磴口县生态空间网络的直径与网络连通度呈现比较明显的下降趋势。其中，网络直径由最初的20降至13,表明网络中任意两个节点之间发生连接所需要经过的其他节点随着空间粒度的升高而减少。在不存在地理意义的复杂网络中，网络直径的减小可以减少网络中物质、能量、信息流动的成本，而在具有地理意义的生态空间网络中，网络直径的减小往往预示着生态廊道中踏脚石节点的消失，不利于网络结构与功能的完整性。磴口县生态空间网络的连通度在空间粒度10～30 m下出现了小幅上升，之后呈现明显的下降趋势，由最初的10 871减小到218，生态空间网络连通度的减少表明将网络变为平凡图需要移除的节点越来越少，网络的不稳定性随着空间粒度的增大而增强。

图15 网络直径变化Fig.15 Diameter variation curve

图16 连通度变化Fig.16 Connectivity variation curve

图17 平均介数变化Fig.17 Average betweenness variation curve

图18 平均路径长度变化Fig.18 Average path length variation curve

由图17、18可知，平均介数和平均路径长度呈现先升高后降低的变化趋势。其中，平均介数在空间粒度10～50 m由1.22×105升高到6.33×105，升高趋势明显，之后迅速下降，在粒度1 000 m下降低为2.20×103，网络平均介数的降低表明，随着空间粒度的增大，网络中各节点上所经过的最短路径在减少，单个节点相对于整个网络的重要性下降。磴口县生态空间网络的平均路径长度在空间粒度10～200 m下基本呈上升趋势，由5.247升高到7.262，之后不断下降，在粒度1 000 m下变为5.778，网络平均路径长度表示网络中任意两个节点发生连接所必须经过的其他节点的数量，在生态空间网络中，路径长度的升高表明廊道中踏脚石数量的增加，而减少则可能预示着踏脚石节点的消失。结合网络直径的变化趋势可以发现，平均路径长度的降低和网络直径的明显减小均发生在空间粒度200 m以后，这很可能是由于面积较小的生态斑块在较大空间粒度下被忽略引起的。

由图19、20可知，磴口县生态空间网络的平均聚类系数和平均度不存在明显的变化趋势。网络的平均聚类系数随着空间粒度的升高呈现上下波动的状态，其值在0.25左右徘徊。网络平均度的变化幅度明显高于聚类系数，然而经历几次上下波动之后，在空间粒度1 000 m下网络平均度达到3.04，这与空间粒度10 m下平均度3.02十分接近。

图19 平均聚类系数变化Fig.19 Average clustering coefficient variation

图20 平均度变化Fig.20 Average degree variation

3 讨论

本研究结果表明，空间粒度的变化对生态空间网络的相关研究具有明显影响，粒度效应是生态空间网络研究中必须考虑的重要问题。主要原因可能是采用数据的空间粒度与景观斑块不匹配。根据O’NEILL等[29]的研究，当使用栅格数据计算景观格局指数时，如果景观斑块面积小于数据的空间尺度，该斑块在计算过程中将被忽略，从而导致计算误差。

基于此，邱扬等[17]提出通过各指数变化的拐点粒度来确定最优粒度的方法。该方法认为，在所有粒度的计算结果中，最小的初始粒度有最高的精度，因为该粒度下被忽略的景观斑块最少。从初始粒度逐渐增大空间粒度的过程中，寻找各指数变化的拐点，该拐点出现之前的空间粒度为最优，即该粒度既可保证计算精度，又可免去小尺度下过大的计算量。

利用这一方法，本研究中度-度相关性、聚-度相关性两个指数计算的最优粒度分别可大致确定为200～300 m、80～90 m。本研究中其他指数由于初始波动较大，其最优粒度即为最小粒度。该方法的确可以在某一研究区识别出各指数计算的最优粒度，然而其局限在于，为了达到这一目的，不得不在每一粒度下计算各项指数，然后寻找拐点，这种做法显然工作量过大，不具备很高的应用价值。然而这种方法仍然为解决这一问题提供了很好的思路，可以在多个研究区开展生态空间网络粒度效应的研究，寻找拐点出现的粒度同研究区面积、研究区中最小斑块、最大斑块、平均斑块等指标之间的关系，构建一种基于上述各项指标的最优粒度推算方法，从而解决最优粒度选取的问题。

4 结束语

空间粒度的变化对生态空间网络的空间结构、拓扑结构以及生态源地的面积与空间分布均具有明显的影响，而生态阻力的粒度效应并不明显。随着空间粒度的不断增大，生态源地面积产生了明显的波动并且变得越发破碎；生态空间网络长度经历了从增加到减少的变化，其空间分布由最初的不断密集再到逐渐稀疏；生态空间网络的各项拓扑结构指数变化较为明显，随着空间粒度的不断升高，网络的度分布愈发集中，同配性、层次性、模块化增强，网络不均匀性不断提升，连接性不断降低，总体朝着更加不稳定的方向变化。