基于土地利用变化的黄河流域景观格局及生态风险评估

刘希朝，李效顺，蒋冬梅

（1. 中国矿业大学中国资源型城市转型发展与乡村振兴研究中心，徐州 221116；2. 自然资源部海岸带开发与保护重点实验室，南京 210095；3. 江苏自然资源智库中国矿业大学研究基地，徐州 221116）

0 引 言

流域是一种结构复杂的地理区域，隶属于生态、经济及社会多个系统，具有生产、生活、生态及文化等多项相互联系的功能[1-2]，过度利用其中一种或几种功能势必对流域整体功能产生负面影响[3-4]。在城镇化迅速发展的背景下，对流域的开发利用主要聚焦在其经济功能的发展，而对生态功能关注较少。研究流域生态风险能够为流域管理和风险防控提供决策依据，进而促进流域可持续高质量发展[5]。

生态风险评价是为生态风险管理服务的，美国环境保护署于1992年定义生态风险为评估由于一种或多种外界因素导致可能发生或正在发生的不利生态影响的过程[6-7]。早期的生态风险评价主要针对环境污染展开研究，评价规模多数是单一风险源和单一风险受体。近年来，有关生态风险评价的研究较多关注生态系统的整体影响以及生态风险的空间相关性，评价规模也扩展到区域尺度，如流域、城市群、经济带等[8-9]。随着对土地利用和生态风险研究的深入，从土地利用视角评估生态风险逐渐成为研究的主流，主要评价方法有景观指数法和相对风险模型（Relative Risk Model，RRM）[6]。景观指数法是选择与生态风险相关的景观格局指数构建评价模型[10]，如张玉娟等[3,7,11-13]选取不同类型的景观格局指数构建生态风险指数（Ecological Risk Index，ERI），分别评价省级区域、流域、海岸带等区域的生态风险。奚世军等[14-15]在计算ERI的基础上，分析ERI的空间相关性，并研究其驱动因素。相对风险模型是一种区域复合压力风险评价模型，应用于水域、陆地等多个生态系统的生态风险评价。如张天华等[4,16-18]运用RRM分别评价水系工程、流域和县级区域的生态风险，并提出相应的风险防控对策。王丽萍等[19]构建了基于RRM思想的多级模糊综合评价模型，得到生态风险模糊评价分值，为改善区域生态环境提供参考。综上，国内外学者主要通过建立评价模型来分析生态风险，评价方法和体系比较成熟，但从多尺度分析区域生态风险时空特征的研究较少。

黄河流域生态环境脆弱，地形复杂，水资源短缺，水土流失严重[20-22]。虽然近年来实施退耕还林、坡改梯等修复工程，但其生态保护工作仍面临挑战。在2019年召开的黄河流域生态保护和高质量发展座谈会上，习近平总书记提出既要合理开发利用资源，又要兼顾资源的可持续性，提升资源的循环利用[23]，自此，黄河流域生态保护正式上升为国家战略。分析长时间序列黄河流域景观格局并从多尺度诊断其生态风险，对流域高质量发展具有重要意义。基于此，本文采用黄河流域2000、2010和2018年的土地利用数据，运用Fragstats4.2计算景观格局指数，构建生态风险评价模型，从网格和县域不同尺度诊断生态风险，揭示生态风险的时空分异规律、转移特征和地类分布，并分析生态风险的空间关联性，通过划分生态风险管控区，提出黄河流域生态保护分区管控建议，以期为黄河流域风险防控和生态保护提供定量参考，并为其高质量发展提供决策支持。

1 材料与方法

1.1 研究区概况

黄河是中国第二长河，全长约5 464 km，流经青海、四川、甘肃、宁夏、内蒙古、陕西、山西、河南及山东9个省（自治区），约390个县（区）[20]，如图1所示。流域面积约79万km2，占国土面积的8.3%，横跨青藏高原、内蒙古高原、黄土高原和黄淮海平原4个地貌单元，地势西高东低。2018年黄河流域（九省区）GDP为2.39×105亿元，占全国GDP总量的26.5%。

1.2 数据来源与预处理

为深入诊断黄河流域景观及生态风险时序特征，本文将研究时段定为2000—2018年。数据均来源于中国科学院资源环境科学数据中心（http://www.resdc.cn），包括2000、2010、2018年全国土地利用遥感监测数据（栅格数据，分辨率为1 km）和中国行政区划数据（矢量数据）。按照研究区域的行政边界裁剪土地利用数据，得到黄河流域2000、2010、2018年的土地利用数据，运用ENVI5.3对其进行精度检验，三期影像的综合精度达91%以上。根据本文的研究需要，在ArcGIS中对其进行重分类等处理，将土地利用类型划分为耕地、林地、草地、水域、建设用地和未利用地6大类（如图2）。

1.3 研究方法

1.3.1 景观格局指数

景观格局包括景观单元的类型、数量和空间分布，用景观格局指数表达。为避免指标选取过多带来的重复性研究[24-28]，根据黄河流域土地利用情况，综合考虑各景观格局指数的生态学意义，本文从景观水平和类型水平两方面选取7个景观格局指数（表1）分析黄河流域景观要素的数量、形状和空间分布等特征。景观水平指数可以反映研究区整体的景观格局特征，类型水平可以反映研究区内不同斑块类型的结构特征。运用Fragstats4.2软件处理经预处理后的土地利用数据（栅格数据，1 km×1 km），可得到各景观指数的数值。

1.3.2 生态风险评价模型

从景观生态学角度出发，为建立景观结构和生态风险的关系[7,11]，选取能够反映区域生态损失的景观格局指数[25-27,29-30]：景观干扰度指数和景观脆弱度指数，构建生态风险评价模型，具体如下：

1）景观干扰度指数（Di）

Di反映不同景观生态系统受到外界干扰的程度，由景观破碎度指数（Fi）、景观分离度指数（Si）、景观分维数（FDi）计算得出[7,11]。Fi是指景观的空间分裂破碎程度，可以反映自然或人为干扰对景观的影响，其值越大，表明景观的稳定性越差，受到干扰的程度越大。Si是指景观中不同斑块的分离程度，其值越大，表明景观越分散，景观分布越复杂。FDi是用来测定斑块形状对内部斑块生态过程影响的指标，其值越大，表明斑块形状越复杂。指数计算公式如下：

表1 景观格局指数应用尺度、生态学意义及计算方法Table 1 Application scale, ecological significance and calculation method of landscape pattern index

式中a、b、c为权重，且a+b+c=1，根据其作用程度[25-27,29-30]，分别赋值为0.5、0.3、0.2；ni为第i种景观的斑块数；Ai为第i种景观的面积，km2；A为所有景观的总面积，km2；Pi为第i种景观的周长，km。

2）景观脆弱度指数（Vi）

Vi反映不同景观生态系统对外部干扰的敏感性，Vi值越大，生态系统越不稳定，越容易受到损害[25-27,29-30]。结合研究区的特点，对6类景观赋值（耕地-4，林地-2，草地-3，水域-5，建设用地-1，未利用地-6）[25-27]，并做归一化处理。

3）生态风险评价模型（ERI）

根据景观干扰度指数和景观脆弱度指数构建ERI，将景观空间结构转化为空间生态风险。ERI用于描述评价单元内的生态损失程度，其值越大，生态风险越高[7,11,29-30]。计算公式如下：

式中ERIk为评价单元k的生态风险指数；n为景观类型的数量；Aki为景观生态风险评价单元k中第i类景观的面积，km2；Ai为景观生态风险评价单元k的面积，km2。

为充分体现生态风险的空间分布情况，本文从网格和县域两个尺度展开研究，两尺度的研究结果可以进行对比，也可以互相印证，使生态风险评价结果更具有科学性。在网格尺度上，以黄河流域景观斑块平均面积的5倍（30 km×30 km）创建渔网[13-14]，将其划分成1 078个风险单元；在县域尺度上，黄河流域共包含391个县区，每个县区作为1个风险单元，如图3所示。

1.3.3 空间分析法

采用Moran’s I指数（I）和LISA指数（Ii）分析黄河流域生态风险的空间相关性。Moran’s I统计量是一种应用非常广泛的空间自相关统计量，可以反映空间相邻或相近单元属性值的相似程度，公式如下：

式中Yi、Yj为变量在相邻配对空间单元的取值；ωij为空间权重矩阵；¯为属性值的平均值。I的取值在[-1，1]之间，当I＞0时，表明研究单元的观测值趋于空间聚集，空间正相关；当I＜0时，表明空间呈离散分布状态，空间负相关；当I=0时，表明空间不相关。

LISA指数又称局部Moran’s I指数，可以反映某区域与相邻区域间的差异程度及显著性，公式如下：

式中n′为样本数量，即研究单元的数量，S2为统计量的方差。当Ii＞0时，表示一个观测值高（低）的区域被一个高（低）观测值的区域包围，即“高-高”（“低-低”）聚集；当Ii＜0时，表示一个观测值高（低）的区域被一个低（高）观测值的区域包围，即“高-低”（“低-高”）聚集；Ii=0时，表示观测区域与相邻区域无关联，即不显著。

2 结果与分析

2.1 黄河流域景观格局分析

2.1.1 土地利用变化

根据数据预处理的结果，统计2000、2010和2018年黄河流域各地类的面积及变化情况。如表2所示，黄河流域土地利用类型以草地和耕地为主，草地占比最大。从时序特征来看，2000—2018年耕地和未利用地面积分别减少14 243 km2和9 410.3 km2，林地、草地、水域面积在2000—2010年间维持相对稳定，2010—2018年间有小幅度增加，这与国家自2003年开始实行的退耕还林（草）政策相关。退耕还林（草）在一定程度上调整了土地利用结构，改善了农业生产条件。2000—2018年建设用地面积大幅度增加，增加面积为12 179.8 km2，尤其在2010—2018年，面积增长了68.36%，这是城镇化建设和经济发展的共同需求。现阶段，黄河流域仍处在快速城镇化发展的时期，对建设用地的刚性需求强劲。

表2 2000—2018年黄河流域地类面积及变化情况Table 2 Area and changes of land use types in the Yellow River Basin from 2000 to 2018

2.1.2 整体景观特征

由表3可以看出，2000—2018年，黄河流域整体斑块数（NP）呈上升趋势，增长幅度为5.83%，整体斑块密度（PD）呈上升趋势，增长幅度为6.12%。由此可知，黄河流域景观破碎化日益严重，破碎化进程加快，景观分布分散。最大斑块指数（LPI）呈上升趋势，增长幅度为7.14%，说明黄河流域的优势斑块类型控制景观的作用增强。景观形状指数（LSI）呈上升趋势，说明景观格局的形状变得更为复杂和不规则，形状种类多样化。香农多样性指数（SHDI）呈上升趋势，增长幅度为1.49%，说明景观类型丰富并趋于多样化。蔓延度指数（CONTAG）呈下降趋势，下降幅度为3.80%，说明流域景观聚集程度降低，景观破碎化加剧，景观的抗干扰能力下降。整体来看，黄河流域受人类活动影响严重，生态功能有所损害，景观格局趋于复杂化、破碎化和分散化。

表3 黄河流域景观水平的景观格局指数计算结果Table 3 Calculation results of landscape pattern index of landscape level of the Yellow River Basin

2.1.3 分类型景观特点

由表4可以看出，2000—2018年黄河流域各类型的斑块数量（NP）和斑块密度（PD）均有不同程度的增减，其中，耕地、林地、水域和建设用地的斑块数量和斑块密度持续增加，草地和未利用地减少，表明各地类的破碎化程度及变化趋势不同。其中，林地、耕地破碎化程度较高，主要是因为建设用地的不断扩张和交通的发展，大片的耕地、林地被占用、分割。各类型的景观形状指数（LSI）与斑块数量（NP）的变化趋势一致，但变化程度较小，其中，草地和耕地的形状最复杂。耕地、林地和未利用地的分离度指数（SPLIT）呈上升趋势，草地、水域和建设用地呈下降趋势，且建设用地的下降幅度最大，下降速度最快，说明建设用地斑块不断趋于聚集，但水域和建设用地的指数值仍较大，说明其斑块间的分离程度与其他类型相比仍较高。

2.2 黄河流域生态风险分异

2.2.1 生态风险时空演变

根据1.3.2中的生态风险评价模型，测算得到每个风险单元的生态风险值，在ArcGIS中将其赋值给各单元的中心点，运用半方差变异函数优化数据，通过克里金插值，得到2000—2018年黄河流域在网格尺度和县域尺度的生态风险空间分布图（图4），可以看出黄河流域生态风险的时空分异特征十分明显。

基于自然断点法，将生态风险分为5级，对应5个等级的风险区：低风险区（ERI≤0.46）、较低风险区（0.46＜ERI≤0.50）、中等风险区（0.50＜ERI≤0.55）、较高风险区（0.55＜ERI≤0.60）、高风险区（ERI ＞0.60），如图4所示。时间序列上，2000—2018年黄河流域生态风险值的最高值呈上升趋势，高风险区所占比例逐渐增加，网格尺度生态风险最高值由0.62增长至0.74；县域尺度生态风险最高值由0.62增长至0.73。空间分布上，2000年高风险区主要集中在下游的山东段、河南段，上游的四川段（网格尺度）；2010年，中游的陕西段南部（西安、咸阳、渭南、商洛），山西段南部（晋城、运城、长治）等地区也呈现高风险；2018年网格尺度上，内蒙古的鄂尔多斯市等地区也呈现高风险。这些地区生态用地分散，建设用地面积不断扩张，造成生态环境压力，导致生态风险持续升高。

统计得到各风险等级的面积占比，如表5所示。网格尺度上，2000年以低风险和较低风险区为主，2010年和2018年以中等风险和高风险区为主，整个研究期间，中等风险、较高风险和高风险区面积不断增加，其余等级区面积减少，其中高风险区面积占比提高了22.87%，低风险区面积占比下降了32.03%。县域尺度上，2000和2010年以较低风险区和中等风险为主，2018年以中等风险和较高风险区为主，整个研究期间，较高风险和高风险面积不断增加，其余等级区面积减少，其中高风险区面积占比提高了22.22%，低风险区面积占比下降了24.55%。

表5 生态风险各等级面积占比Table 5 Area proportion of each ecological risk level %

两种尺度对比可以看出，网格尺度各等级面积的变化幅度较大，且高风险区和低风险区面积占比均大于县域尺度，说明网格尺度下的生态风险变化敏感，较易产生两端值。

2.2.2 生态风险转移特征

通过叠加2000年和2018年生态风险分布图，得到2000—2018年黄河流域各生态风险等级转移情况。由表 6可以看出，在网格尺度上，面积转移前3名依次为：低风险转为中等风险（转移面积119 437.00 km2），低风险转为较高风险（转移面积98 816.10 km2），较低风险转为中等风险（转移面积77 808.70 km2）。在县域尺度上，面积转移前3名依次为：中等风险转为较高风险（转移面积123 095.00 km2），较低风险转为中等风险（转移面积108 254.00 km2），中等风险转为高风险（转移面积101 085.00 km2）。由此可见，黄河流域的生态风险大都是由低等级向高等级转移，生态风险加剧，进而导致生态问题严峻。

2.2.3 生态风险地类分布

运用ArcGIS统计工具得到2000—2018年黄河流域不同地类在各等级生态风险中的面积分布情况，如表7所示。2000年耕地主要分布在低风险区，占耕地总面积的47.97%，2010年和2018年主要分布在高风险区，占比分别为33.70%、37.79%，说明耕地所在区域的生态风险加剧，主要是人类活动干扰的结果，由于开发建设的需求，城镇和农村周边建设用地增加，导致耕地趋于破碎化，生态风险值增加。2000年林地主要分布在低风险区，占林地总面积的43.24%，2010年和2018年主要分布在高风险区，占比分别为40.30%、31.07%，说明林地所在区域的生态风险先增加后降低，生态风险的降低主要与国家实施的退耕还林及林地保护等政策相关。2000年草地主要分布在低风险区，占比43.24%，2010年主要分布在较低风险区，占比35.63%，2018年主要分布在中等风险区，占比37.00%，整体来看草地所在区域的生态风险增加，但分布在高风险等级区的面积比例较小。2000年和2010年水域主要分布在中等风险区，2018年分布在高风险区，占比46.18%，说明水域在2010—2018年生态风险剧增，相关水利设施的建设以及建设用地的开发对水域生态造成影响。2000—2018年，建设用地和未利用地的分布区域逐渐以高风险为主，这与建设用地规模增加，未利用地开发利用和存量土地更新相关。

表6 2000—2018年黄河流域生态风险转移情况Table 6 Ecological risks transfer of the Yellow River Basin from 2000 to 2018 km2

表7 2000—2018年黄河流域各地类的生态风险面积Table 7 Ecological risk area of various types in the Yellow River Basin from 2000 to 2018

2.2.4 生态风险空间关联

运用GeoDa软件，根据黄河流域2000、2010、2018年的生态风险情况，分析其全局自相关性，结果如表8所示。在整个研究期间，黄河流域生态风险的Moran’s I大于0，呈现空间正相关，即生态风险在空间上相互影响，具有空间相似性。网格尺度上，三期的Moran’s I先上升后下降，县域尺度上，三期的指数值持续下降，且网格尺度的Moran’s I大于县域尺度，表明网格尺度生态风险的空间正相关性更强。总体来看，黄河流域生态风险的Moran’s I值在近20年间呈下降趋势，这表明随着社会经济的发展和土地利用的变化，生态风险的空间聚集程度和空间分异性减弱。

同理，运用GeoDa软件分析生态风险的局部自相关性。局部自相关可以反映生态风险的聚集类型和空间位置，如图5 所示。网格尺度的“高-高”聚集单元（“热点”区域）主要集中在下游的河南和山东，中游的山西南部和陕西南部，上游的四川等区域，2000—2010年“低-低”聚集单元（“冷点”区域）主要集中在上游的甘肃中部（白银、兰州、庆阳、定西）、中游的陕西北部、山西北部（朔州、沂州）、内蒙的鄂尔多斯、呼和浩特以及青海的玉树、果洛等区域。2018年网格尺度的“冷点”区域减少，主要集中在陕西北部、山西北部以及内蒙的呼和浩特。县域尺度2000年的“热点”区域主要集中在下游的河南和山东以及上游的四川，2010年“热点”区域集中在中游的山西、陕西南部以及下游的大部分区域，2018年“热点”区域仅剩21个，集中在河南的洛阳、三门峡以及山东东营区域。县域尺度的“冷点”区域与网格尺度基本一致。

表8 黄河流域生态风险的Moran’s I值Table 8 Moran’s I of ecological risk in the Yellow River Basin

整体来看，“热点”区域和“冷点”区域的数量先增加后减少，总体上呈减少趋势。空间聚集不显著单元数量增加，说明生态风险局部空间聚集程度降低。

2.3 黄河流域生态风险管控区划分

黄河流域生态风险的空间差异较大，本文综合网格尺度和县域尺度生态风险诊断结果（图4），根据生态风险防控的需要，将两尺度空间上风险等级不一致区域以较高的风险等级为主进行调整，且为方便区域管理，按照不跨越市级行政区的原则，将黄河流域划分为生态风险重点管控区、严格管控区和一般管控区，如表9所示。

表9 黄河流域生态风险空间管控区Table 9 Ecological risk spatial control area of the Yellow River Basin

1）生态风险重点管控区：该区域生态用地分散，建设用地所占比例较大，且面积不断扩张，造成生态环境压力，生态风险也持续升高，应减少对建设用地的增量开发，转向存量更新和集约高效利用，并全面保护林地、草地、水域等生态用地和高质量耕地。

2）生态风险严格管控区：该区域各类用地分布较为均衡，应发挥国土空间规划的引领作用，建设用地的增量开发与存量更新同行，在促进经济发展的同时，尽量减少对生态的负面影响，维持当前生态风险水平。

3）生态风险一般管控区：该区域主要以农牧业为主，生态风险较低。应合理开发布局建设用地，加强基础设施建设，同时维持农业发展优势，推进农业现代化。

3 讨 论

2000、2010和2018年黄河流域（九省区）GDP分别为2.41×104、1.18×105和2.39×105亿元，年均增长率14.43%，城镇化率分别为32.11%、45.16%、56.20%，虽然城镇化率不断提升，但均低于全国平均水平（36.22%、49.68%、59.58%），说明黄河流域现阶段仍处在城镇化快速发展时期。与此同时，由于发展方式粗放、产业结构不合理以及建设用地无序扩张等问题，导致黄河流域的整体景观格局破碎、生态风险加剧，生态环境压力加重[31]。在黄河流域内部，区域发展不平衡现象突出，山东、河南等下游城市的经济发展和城镇化水平高于大部分中游、上游的城市，土地利用变化较为剧烈[20-21,31]，这也导致了生态风险的空间分布差异：高风险区主要集中在下游的山东、河南区域，低风险区主要集中在上游的宁夏、甘肃区域。总体而言，黄河流域的生态环境保护水平滞后于城镇化发展，且快速城镇化加剧了区域生态风险。

划分生态风险评价单元是诊断区域生态风险的基础。基于网格的划分方式有利于生态风险空间异质性的表达，是现阶段应用较广的划分方法。多数学者在研究中、小流域尺度或市、县域尺度时（如喀斯特山区小流域、南四湖流域），都会通过划分网格单元研究生态风险[3,14-15]。也有部分学者通过设置不同大小的网格单元，研究生态风险变化的尺度特征。基于行政区的划分方式可辅助决策者制定更具有适宜性的风险防控及管理政策[32]，其应用较少且多应用于大尺度的研究。本文从网格尺度和县域尺度分别诊断黄河流域生态风险的时空分异特征，且由上述分析结果可知，两尺度生态风险的空间分布情况和变化趋势相近，可以互相印证，提升生态风险评价结果的科学性。综合两尺度的结果，可以为流域高质量发展提供更为合理的建议。但两尺度的结果又有不同，网格尺度更微观，对生态风险变化的反映更敏感，可为县域尺度或更高一级的行政区内部微观生态风险管理提供参考，这也是论文需要进一步探索的方向。

4 结 论

1）土地利用结构分析表明，黄河流域的主要土地利用类型为草地和耕地。2000—2018年，耕地面积减少14 243 km2，未利用地面积减少9 410.3 km2，建设用地面积增加12 179.8 km2。现阶段黄河流域仍处在快速城镇化建设时期，对建设用地的刚性需求强劲。

2）景观格局计算结果显示，2000—2018年，黄河流域整体斑块数、斑块密度、最大斑块指数和香浓多样性指数的增长幅度分别为5.83%、6.12%、7.14%和1.49%，蔓延度指数的下降幅度为3.80%。黄河流域整体景观格局趋于复杂化、破碎化和分散化。其中，林地、耕地破碎化程度最高，草地、耕地的景观形状和结构最复杂。

3）生态风险评价结果显示，2000—2018年，黄河流域网格尺度和县城尺度生态风险的最高值由0.62分别增长至0.74和0.73，网格和县域尺度高风险占比分别提高了22.87%、22.22%。高风险主要集中在下游的山东段、河南段，上游的四川段，中游的陕西段南部和山西段南部。且生态风险大多由低等级向高等级转移。

4）流域空间关联分析看出，2000—2018年，黄河流域生态风险的Moran’s I均大于0，呈现空间正相关，但Moran’s I呈下降趋势，说明生态风险的空间聚集程度和空间分异性减弱。网格尺度的Moran’s I大于县域尺度，表明网格尺度生态风险的空间正相关性更强。