鄂尔多斯—榆林地区景观生态风险评价及其驱动因子分析

马 骏, 裴燕如, 王慧媛, 于 强, 牛 腾, 岳德鹏

(北京林业大学 精准林业北京市重点实验室, 北京 100083)

随着社会经济发展，人类城镇化加快。一方面土地利用类型不合理利用，人地矛盾尖锐；另一方面大多数国内的矿产资源不能够合理地开发利用，环境代价巨大，使绝大多数的自然生态系统都受到人类活动的影响，而这些人类活动影响会造成诸多的生态风险，因此科学管理这些生态风险成为人类与自然和谐发展的重要前提[1]。生态风险评价作为现阶段生态环境修复、治理以及生态安全格局构建和生态风险预警等诸多工作的重要参考依据，引起了全球学者的高度重视。

生态风险指生态系统及其组分在受到外界干扰后所承受的风险，其会对生态系统的结构和功能会造成不利影响[1-2]，生态风险起源于小尺度的人体健康评价，后来演变为衡量一种或多种胁迫因素对生态系统造成正在发生或潜在发生的具有负生态效益的可能性的研究手段[2-4]。随着生态风险评价理论和方法的完善，区域生态风险评价的提出，将空间异质性和风险因子等级关系纳入到生态风险评价体系。区域生态风险主要应用在土壤、大气以及水环境的生态评估。赵晓光等人[5]对矿区不同土地利用类型的土壤重金属污染程度和综合潜在生态风险危害进行评价,何瑞东等人[6]分析郑州市某生活区大气重金属元素的污染特征，对郑州市某生活区进行潜在生态风险和居民健康风险的评估。李捷等学者[7]探究了北方六湖沉积物中重金属污染源的差异性。为了探究景观格局对生态风险的定量影响，基于区域生态风险评价提出景观生态风险评价，景观生态风险评价更加注重尺度效应、空间异质性以及景观格局对的生态功能、过程的影响。目前国内外对景观生态风险的研究已经相对成熟，从研究对象上来看，流域[8-10]、矿区[11-12]、城市[13]、自然保护区[14]等已经成为研究的热点，从研究尺度来看，综合考虑研究区面积和景观斑块的平均面积来划分网格是主流的研究手段[15-17]，也有部分学者从县域、市域等行政单位为研究单元进行生态风险研究[18-19]；从评价方法上来看，大多基于风险源汇和景观格局两种方法，有的学者还基于生态系统服务功能进行景观生态风险评估的研究[20]，当前关于景观生态风险评价的研究主要将土地利用与生态风险相结合并探究其时空分布，鲜有涉及到不同驱动力影响生态风险演变过程的定量分析，在本研究中使用地理探测器探测鄂尔多斯—榆林地区(以下简称“鄂榆地区”)生态风险演变过程中的关键驱动因子。

鄂榆地区矿产资源丰富，被誉为中国的“能源走廊”。进入21世纪，西部大开发战略实施，西部资源开发产业化进程加快，煤炭能源开发成为鄂榆地区的支柱产业[21]。随着人为开发活动的进行，对整个鄂榆地区的草地、林地、未利用地等景观要素都造成了不同程度的影响，在很大程度上改变了景观结构的布局，并对生态安全产生威胁。基于此，本文以鄂榆地区为研究区，基于土地利用变化信息和景观格局构建景观生态风险评价模型，测定各评价单元的景观生态风险值，通过地理探测器对引起生态风险的驱动力因素进行探测，确定生态风险时空演变的主要驱动因子，为生态环境治理、修复以及生态安全格局的构建提供科学依据。

1 数据来源与研究方法

1.1 研究区概况

研究区位于我国第二大沉积盆地——鄂尔多斯盆地腹地，地处内蒙古自治区鄂尔多斯市和陕西省榆林市两个地级市，并以此为边界向外扩长5 km,将对鄂榆地区生态环境有重要影响的黄河及黄河周边的地区也能纳入研究区范围内。

鄂尔多斯市位于内蒙古自治区西南部，介于北纬37°35′24″—40°51′40″，东经106°42′40″—111°27′20″，总面积为86 752 km2。榆林市位于陕西省最北部，介于北纬36°57′—39°35′，东经107°28′—111°15′，总面积为43 578 km2。鄂榆地区的东部、北部、西部三面被黄河环绕，北侧与库布齐沙漠相接，毛乌素沙地则地处陕西省榆林地区和内蒙古自治区鄂尔多斯市之间，属于北温带半干旱大陆性气候区，地势由西北向东南倾斜，土地利用类型主要以草地和未利用地为主，冬夏寒暑变化大，多年平均气温7.2 ℃，日最高温度39 ℃，日最低气温-31.4 ℃。多年平均降水约357 mm，降水主要集中在7—9月，占全年总降水的70%左右。鄂榆地区人口结构主要以汉族占大多数，根据第7次人口普查数据，鄂尔多斯市常住人口为2 153 638人，榆林市常住人口为3 624 750人。鄂榆地区的产业结构主要以第二产业为主，第二产业的比例大致保持在50%左右。2020年，鄂尔多斯市实现生产总值(GDP)为3 533.66亿元，榆林市实现生产总值(GDP)为4 089.66亿元。

1.2 数据来源及处理

本文数据主要选取鄂榆地区2000和2010年30 m空间分辨率的Landsat 7 TM/ETM+遥感影像和2020年30 m空间分辨率的Landsat 8OLI遥感影像为数据源，将其在ENVI 5.3中对影像进行辐射校正、大气校正、几何校正、图像增强以及镶嵌裁剪等预处理，得到鄂榆地区3期的遥感影像数据。利用极大似然法进行监督分类，根据目视解译并结合研究区实际情况，采取LUCC分类体系将土地利用类型分为耕地、林地、草地、水域、城乡，工矿，居民用地(以下称建设用地)及未利用地6类，在解译结果随机抽取验证样地与Google Earth和实地进行选取的1 217个验证点进行验证，总体分类精度均达88%以上，满足研究需求。驱动因子分析中使用的人为干扰度数据是通过于立忠等[22]提出的人为干扰度赋值表，参考谷东起等[23]建立的人为干扰度计算模型计算并将其空间化所得；NDVI数据为NASA提供的MOD13A1产品(http:∥ladsweb.nascom.nasa.gov/)；年降水量数据、年均气温数据来自于国家地球系统科学数据共享服务平台(http:∥www.geodata.cn/);DEM数据来自地理空间数据云平台(www.gscloud.cn/)；地表蒸散发数据来源于NASA官方遥感影像数据(https:∥ladsweb.modaps.eosdis.nasa.gov/)。

1.3 研究方法

1.3.1 土地利用变化指标构建

(1) 土地利用动态度。在景观生态风险的研究中，土地利用类型通常被看作是景观嵌体的类型。土地利用动态度可以描述景观更替变换的速度与强度，通常定义为指定时间范围内区域土地利用类型变化的快慢[24]，以表现区域内土地利用的剧烈程度。通常分为以下两类：

①单一土地利用类型动态变化度。单一土地利用类型动态变化度指的是某研究区一定时间范围内，某一土地利用类型的面积变化情况[24]，计算公式为：

(1)

式中：K为研究时段内某一土地利用类型动态度；Ua，Ub分别问研究初、研究末期某一土地利用类型面积；T为研究时段(a)。

②综合土地利用类型动态变化度。综合土地利用类型动态变化度表示所有土地利用类型变化的整体情况，计算公式为：

(2)

式中：LC为研究时段内土地利用类型年变化率； LUi为研究初期第i类土地利用类型面积； ΔLUi-j为研究初期至研究末期时段内第i类土地利用类型转为非i类土地利用类型面积的绝对值。

(2) 土地利用类型转移矩阵。土地利用作为景观最直观的表达形式，当强烈人类活动或自然灾害发生时，会对土地利用类型造成变更，景观要素的空间结构以及生态功能也随之发生改变，从而使区域的景观生态风险发生变化。土地利用类型转移矩阵描述研究时段内土地利用类型的变更过程，能够充分反映研究时段内土地利用类型的变化方向，计算公式为：

(3)

式中：Aij为土地利用类型面积；n为土地利用类型数量；i，j为转移前后的土地利用类型。

1.3.2 景观生态风险指数构建 土地利用速度、强度以及方向的变化会引起景观格局的改变，不同的景观格局对生态风险的抵御能力是不相同的。为了将景观生态风险在空间上以定量化的形式呈现，并改善现状中不利的环境要素，在未来的设计规划中提出发展性意见，通常通过计算多种景观指数，按照不同的权重构建景观生态风险指数。

(1) 划分评价单元。为了将生态风险指数空间化，结合苏海民等研究[25]，评价单元大小的选择应为斑块平均面积的2～5倍，利用ArcGIS 10.2渔网构建工具对研究区处理，分割为5 km×5 km大小相同的生态风险单元，共划分5 864个，作为评价单元。

(2) 景观生态风险模型构建。景观格局是在多种驱动力在景观结构综合作用形成的，反映了人类活动对景观结构作用的强度、频率以及范围，同时也影响着生态系统的格局和过程。根据相关文献[17,23]利用景观干扰度和景观脆弱度构建景观生态风险评价模型，计算每个评价单元的景观生态风险指数，将其作为评价单元的中心值。

景观干扰度是景观受到外界因素干扰程度的一种定量表征，将景观破碎度、景观分离度、景观优势度按一定的权重叠加计算所得[26-30]，结果详见表1。

表1 景观指数计算公式及意义

景观脆弱度(Fi)表示不同类型景观对外界干扰的敏感程度和易损程度，景观脆弱度在不同阶段的景观自然演替过程中是不同的[31]。根据前人的研究成果[32]，将研究区的景观分别进行如下赋值：未利用地6，水域5，耕地4，草地3，林地2，建设用地1，归一化得到景观类型脆弱度指数。

通过景观生态风险指数可以建立起景观结构与生态风险的关系，能够定量地反映景观格局与生态风险的相关性，计算方法如下为：

(4)

式中：Aki表示第k个样区i类景观的面积；Ak表示第k类样区的总面积。

1.3.3 空间自相关分析 景观生态风险在空间上量化后，通常会在空间上按一定的规律分布，而空间相关性分析是对景观生态风险空间分布模式描述的一种有效方法。空间相关性是检测空间中某一要素的属性值是否与相邻空间点属性值相关联的指标。正相关表示某空间单元的属性值与其相邻单元的属性值有相同的变化趋势，负相关则相反。常用的指标有Moran’sI指数和LISA指数。

(1) Moran’sI指数(全局自相关)是描述空间要素及属性值在区域内的空间依赖程度，其计算公式为：

(5)

(2) LISA指数(局部自相关)是描述空间单元与其相邻单元的相似程度，其计算公式为：

(6)

1.3.4 地理探测器 景观生态风险是由多种风险因子对区域内景观的综合影响所产生，是由多种驱动力，包括自然环境、社会经济环境以及人类活动等多种条件综合作用所致。地理探测器是探测空间分异性并揭示其背后驱动因子的一种统计学方法。其核心思想是：如果自变量对某个因变量有很强的关系，那么自变量和因变量在空间分布上具有一定的相似性。其中地理探测器中的主要用来因子探测器用来探测自变量对因变量的解释力,其大小用q值来衡量[33]。因此本文选取因子探测器来分析影响鄂榆地区景观生态风险的各个驱动因子,计算公式为：

(7)

(8)

2 结果与分析

2.1 土地利用类型时空变化分析

2.1.1 土地利用结构变化 根据图1和表2可得，在鄂榆地区，草地、耕地和未利用地所占的面积比例最大，共占90%以上，在2000—2010年，林地、建设用地、未利用地呈现增加的趋势，其中建设用地的面积增加最多，相比2000年增加了0.29%；耕地、草地及水域分别呈现减少的趋势，耕地的面积减少最多，相比2000年减少了0.21%，其他土地利用类型的变化幅度比较小，总的来说2000—2010年的趋于相对稳定的状态；在2010—2020年，耕地、水域和建设用地的面积有较明显的增长，其中耕地面积增加比例最大，相比2010年增加2.01%，其次为建设用地，相比2010年增加1.32%。相反草地的面积减少比例最大，相比2010年减少2.17%，其他土地利用类型变化幅度较小，处于相对稳定的状态。

图1 鄂榆地区2000，2010和2020年土地利用变化

表2 鄂榆地区2000-2020年土地利用类型面积变化

2.1.2 土地利用速度变化 鄂榆地区的单一土地利用动态度和综合土地利用动态度详见表3，在2000—2010年和2010—2020年2个时间段内，建设用地的土地利用动态度均为最高，分别为5.75和16.47，表明在2000—2020年建设用地的增加速度始终正向增长，且2010—2020年建设用地的扩张增加速度比2000—2010年要快，增加的建设用地主要分布在乌海市、达拉特旗、鄂尔多斯市、准格尔旗伊金霍洛旗、府谷县、神木县、榆阳区、靖边县以及定边县等地区。在研究期间，2000—2010年与2010—2020年的综合土地利用动态变化度分别为0.19与0.69，表明在2010—2020年研究区整体的土地利用类型间变化比较活跃。

表3 鄂榆地区2000-2020年单一动态度和综合动态度

2.1.3 土地利用方向变化 2000—2020年土地利用类型发生了明显的变化(表4)，其中在2000—2010年，草地是主要的转出类型也是主要的转入类型，草地转出为耕地的面积为987.18 km2，占转出类型总比重的45.52%，转出为未利用地的面积为504.17 km2，占转出类型总比重的23.25%，转出为建设用地的面积为211.26 km2，占转出比例的9.74%;转入草地的来源主要是耕地，转入面积为1 152.12 km2，占总转入比例的57.93%，其次为未利用地，转入面积为441.33 km2，占总转入比例的22.19%。其他土地利用类型的转入和转出面积较少。在2010—2020年，草地仍然为主要的转出类型和主要的转入类型，转出耕地面积为5 154.39 km2，占总转出比例的51.30%，转出为未利用地的面积为2 156.52 km2，占总转出比例的21.46%，转出为建设用地的面积为1 375.14 km2,占总转出比例的13.69%；转入草地的主要来源是未利用地和耕地，其中未利用地转入面积为3 131.30 km2，占转入比例的44.71%，耕地转入面积为2 483.85 km2，占转入面积的35.47%，其他土地利用类型的转入转出变化很小。草地景观作为一种生态脆弱性比较低的景观，且在鄂榆地区分布相当广泛，当受到自然和人为因素等作用容易产生性变，而草地作为鄂榆地区主要的景观类型，其数量的增减会导致景观生态风险的直接变化。

表4 鄂榆地区2000-2020年土地利用类型面积转移矩阵 km2

2.2 景观生态风险时空变化分析

在ArcGIS 10.2中采用克里金插值分析方法，对研究区生态风险进行空间插值，然后根据自然断点法，将生态风险区划分为低生态风险、较低生态风险、中生态风险、较高生态风险和高生态风险5类风险区域。

由图2可知，总体上鄂榆地区的中生态风险、较高生态风险和高生态风险地区主要分布在北部的库布齐沙漠和中部的毛乌素沙漠地区，较低生态风险和低生态风险地区主要分布在东部、北部边缘地区、西南部分地区以及库布齐沙漠和毛乌素沙漠之间的地区。由表5可知，2000—2010年，鄂榆地区的中生态风险地区、较高生态风险地区以及高生态风险地区面积呈现扩张趋势，扩大比例分别2.02%，0.22%和0.57%，与较高生态风险地区和高生态风险地区相比，中生态风险地区的面积增幅较大，主要由于在研究期间鄂榆地区实施退耕还林(草)、沙漠农业试验田以及围封禁牧等大量生态保护政策，对鄂榆地区库布齐沙漠和毛乌素沙漠地区沙漠化逆转起到一定的作用，沙漠地区出现草地、林地和耕地等景观类型，植被作物从无到有，区域景观类型多样性增多，导致斑块破碎度增大，脆弱性增强，生态风险呈现升高趋势。

表5 鄂榆地区景观生态风险等级面积及比例

在2010—2020年，鄂榆地区的中生态风险地区、较高生态风险地区及高生态风险地区的面积呈现大幅降低趋势，减少比例分别为4.28%，1.09%和2.06%，同时低生态风险地区和较低生态风险地区面积大幅增加，增加比例分别为4.48%和2.95%。由图2可知，在西北部的库布齐沙漠地区较高生态风险区和高生态风险区面积减少显著，主要原因是随着多年生态保护政策的实施，鄂榆地区生态环境演变具有明显的生态累积效应，生态功能的可恢复性发挥了作用，草地、耕地景观类型在空间分布上逐渐呈现均匀趋势，景观优势度增加且破碎度减少，脆弱性降低，抵抗外界干扰的能力增强，最终使得生态风险降低。

图2 鄂榆地区2000-2020年景观生态风险空间分布

2.3 景观生态风险空间聚集特征分析

鄂榆地区2000，2010和2020年的全局Moran’sI指数均在0.5以上，分别为0.823，0.826和0.845，说明研究区景观生态风险的空间分布呈显著的正相关关系。散点主要集中在一、三象限，相似值呈现聚集状态，大多数的风险单元在空间上的分布规律呈现高—高和低—低，高生态风险区周围分布着相应的高生态风险区，低生态风险地区分布着相应的低生态风险区。全局自相关Moran’sI指数在2000—2020年呈现上升趋势，表明鄂榆地区风险聚集趋势呈现上升趋势。

鄂榆地区2000，2010和2020年的局部自相关LISA指数分布图如3所示。从空间分布上看，3个时期的高—高生态风险地区主要分布在杭锦旗西北部和东北部、鄂托克前旗的西部、乌审旗及榆阳区，少量高—高生态风险地区分布在靖边县和准格尔旗周围，这些区域往往景观类型多样性高，景观优势度低且景观分布破碎化，人类活动极易引起以林草退化为主等生态问题。低—低生态风险地区主要分布在鄂榆地区杭锦旗、达拉特旗北部边缘地区、准格尔旗、府谷县、神木县、佳县、吴堡县和清涧县东部边缘地区，以及正定县、乌海市、鄂托克旗中部、杭锦旗东南部、鄂尔多斯市、伊金霍洛旗西北部和神木县北部地区，这些地区景观类型以耕地为主，景观优势度高，且斑块的破碎度较低。极少数的高—低和低—高生态风险区主要在分布在高—高和低—低生态风险聚集区周围，说明研究区的生态风险比较稳定，出现急剧变化的可能性较低。

图3 鄂榆地区2000-2020年景观生态风险LISA指数分布

2.4 景观生态风险驱动力分析

鄂榆地区景观生态风险的时空演变受多种外界干扰因素综合影响。例如人为、自然及社会经济因素等。本文采用地理探测器探测各个驱动因子对2000，2010及2020年景观生态风险的解释力(图4)。本文选取人为干扰度、人造地表作为人为因素，利用NDVI值、降水量、年平均气温、高程以及地表蒸散发作为自然因素综合分析鄂榆地区景观生态风险影响机制的分析指标。基于ArcGIS平台对7个因子使用自然断点法进行分级，使用地理探测器测算每个驱动因子的q值，q值越大说明对景观生态风险的解释力越强。

图4 鄂榆地区2000-2020年鄂榆地区景观生态风险因子解释力

由图4可知，在2000—2020年不同时期影响鄂榆地区景观生态风险的驱动因子的解释力变化相对显著。在人为因素中，人为干扰度的q值均大于0.7，且逐年增加，反映了人类活动强度对景观格局的扰动呈现增强趋势，人为干扰度对鄂榆地区的景观生态风险的解释力越来越增强，成为影响鄂榆地区景观生态风险的最主要驱动因素；而人造地表的q值小于0.1，对鄂榆地区的景观生态风险解释力并不强；在自然因素中，NDVI的q值均在0.6左右，对景观生态风险有很强的解释力，年降水量的q值在2000和2010年大于0.3，2020年大于0.1，年平均气温的q值保持在0.1附近，这两种驱动因子在自然因素中仅次于NDVI，但作为植被生长必不可少的因素，在一定程度上也影响着生态风险的变化。高程和地表蒸散发的q值均在0.1以下，其对生态风险的影响作用微乎其微。

在自然因素中，NDVI、年降水量和年均气温是主要的驱动因子，其中景观生态风险与植被覆盖、降水和气温因子的起伏变化密切相关，结合鄂榆地区土地利用数据以及NDVI数据，当地的植被覆盖范围在2000—2020年期间有明显的扩大，库布奇沙漠和毛乌素沙地的治理中，以生态保护为主的防沙治沙取得到了明显的成效，沙区的生态环境从不稳定逐渐趋于稳定，由人工干预的植被群落结构也趋于稳定，景观生态风险逐渐变低；根据鄂榆地区的气象资料，2000—2010年，鄂榆地区的气温呈现降低的趋势，降水呈现波动上升的趋势。2010—2020年，气温呈现平稳状态，降水仍呈现波动上升趋势。与全球变暖和全球降水呈现极端的趋势不同，最可能原因就是人类活动的影响，改变了地表景观类型布局，从而改变了当地的气候条件，进而导致景观生态风险也随之发生变化。

在人为因素中，以生态保护为主的人类活动是造成鄂榆地区生态风险变化的主要驱动因子。为切实推动生态建设发展，从21世纪初西部大开发中央政府政策的实施，众多生态保护措施接连落实。这些年来，鄂榆地区普遍实施退耕还林、退牧还草，禁牧、休牧、轮牧、生态移民等保护措施的推进。鄂榆地区的植被在恢复的过程中，景观的结构发生剧烈变化，景观的数量从单一向综合发展，景观的布局也由简单向复杂转变，从而导致景观格局在空间的布局发生显著变化，景观生态风险也随之出现先增减后降低的变化。

3 结 论

(1) 鄂榆地区主要的土地利用类型为草地、耕地和未利用地，约占整个研究区的90%以上，在2010—2020年，土地利用强度和速度最为活跃。在研究期间，草地作为主要的转入转出类型，转出为耕地的比例增加5.68%，转出为建设用地比例增加3.95%，转出为未利用地的比例减少1.79%，耕地转入的比例减少22.46%，未利用地转入的比例增加22.52%。

(2) 鄂榆地区在2000—2010年，中生态风险、较高生态风险和高生态风险地区面积呈现扩张趋势，低生态风险、较低生态风险地区面积呈现收缩趋势，其主要原因是人为生态修复活动的影响，沙漠地区逐渐出现耕地、草地、林地等土地利用类型，景观类型多样性增强，破碎度变大，景观优势度降低，导致生态风险变高。在2010—2020年，景观中生态风险、较高生态风险和高生态风险低面积呈现收缩趋势，而低生态风险和较低生态风险地区呈现扩张趋势，主要由于多年生态政策的实施，景观类型逐渐均一化，破碎度降低且优势度增强，导致生态风险降低。

(3) 鄂榆地区2000，2010，2020年的全局自相关分析Moran’sI指数分别为0.823，0.826和0.845，说明生态风险在空间分布上呈显著正相关，大多数生态风险单元呈现高—高和低—低分布，少数的生态风险单元在高—高和低—低风险单元周围分布，说明研究区的生态风险比较稳定，出现急剧变化的可能性较低。

(4) 影响鄂榆地区景观生态风险的主要驱动因子是人为因素中的人为干扰度，自然因素中的主要驱动是NDVI，年降水量和年均气温属于次要驱动因子，人造地表、高程以及地表蒸散发对鄂榆地区景观生态风险影响不大。