生态脆弱性是在自然或人为因素作用下, 生态系统的正常结构被扰动且超过自我调节的"阈值"而导致的生态功能大幅降低、生态恢复能力减弱甚至完全丧失恢复能力的现象^[1]。随着经济、社会的不断发展, 城镇化进程加快, 特别是不合理的人为活动, 导致生态脆弱性问题不断涌现。开展生态系统脆弱性研究, 能够在明确脆弱性的基础上, 整体掌握其分布特征, 为生态环境保护和可持续发展提供科学依据^[2]。

国内外学者针对生态脆弱性研究, 大多从生态系统的状态、压力和恢复能力等角度考虑^[3], 采用层次分析法(AHP)^[4]、模糊评价法^[5-6]、主成分分析法^[7-8]、"P-S-R"模型^[9-10]和支持向量机模型^[11]等方法, 选取相应的评价指标体系进行分析^[12]。近年来, KANG等^[13]利用AHP模型, 依据自然、环境和社会经济因素, 运用GIS对研究区生态脆弱性动态变化进行分析; ZHANG等^[14]提出"压力-支持-状态-响应(PSSR)"模型, 对黄河三角洲湿地生态系统脆弱性进行评价; HONG等^[15]针对高度城市化地区的生态脆弱性, 建立了生态敏感性、生态压力和自我恢复力等9个要素12个指标的评价指标体系, 开展城市区域生态脆弱性评价。随着研究的深入, 部分学者发现景观格局对生态过程起着重要作用, 并从景观格局角度对生态脆弱性进行研究。SATI^[16]研究了Garhwal地区与水电工程相关的景观脆弱性和生态恢复问题, 认为人类活动是导致景观脆弱性变化的主要原因; 任志远等^[17]认为景观脆弱性能够反映生态系统的不稳定性和对外界干扰的敏感性, 有助于了解生态环境的状态和变化趋势; 封建民等^[18]基于景观格局指数, 从生态系统压力度、敏感度和恢复力3个方面构建评价指标体系, 评价榆阳区生态脆弱性的空间差异。张月等^[19]利用景观指数构建景观脆弱度指数模型, 将新疆艾比湖流域典型区域划分为低、较低、中、较高和高5个脆弱度等级进行分析。

目前, 针对生态脆弱性的评价方法多样, 但指标体系尚未统一。国内外学者常以省域、市域, 甚至国家、地区等较大区域开展研究, 导致范围过大而缺乏局部指导意义^[20]。胡宝方等^[21]利用3S技术和生态系统健康理论以及景观生态学理论, 采用"P-S-R"模型研究草海湿地生态脆弱性动态变化, 但未对其时空分异特征进行分析。鉴于此, 笔者利用3S技术, 结合景观生态学理论构建景观生态脆弱性评价模型, 基于景观格局对研究区景观类型及景观脆弱度变化特征进行分析, 继而以研究区53个行政村作为评价单元, 从较小尺度和可控单元出发, 研究草海区域景观生态脆弱性变化及时空分异特征, 以期为草海湿地区域生态环境保护、修复和重建, 以及研究区土地利用规划提供科学依据与参考。

1 研究区概况与数据来源 1.1 研究区概况

草海位于贵州西部威宁彝族回族苗族自治县, 地处云贵高原中部, 处于滇东高原向贵州高原过渡的顶点区域, 乌蒙山脉的腹心部分^[22]。草海是我国3大高原淡水湖泊(草海、滇池、青海湖)之一, 同时也是贵州最大的高原天然淡水湖, 正常蓄水面积为1 980 hm², 平均水深1.35 m, 最深处约2.8 m, 受季节性降水影响, 丰水期水位可达2 172.0 m, 枯水期水位降至2 171.2 m^[23]。由于紧邻威宁县城, 草海受人类活动影响较大, 特别是自1958年开始"围湖造田"以来, 到1972年草海湖区水面仅存约5 km², 造成了较为严重的生态环境问题。为有效保护草海及其周边生态环境, 1982年草海开始恢复蓄水, 1985年经贵州省政府批准建立草海自然保护区, 1992年经国务院批准建立国家级自然保护区^[21]。2017年8月中央第七环境保护督察组反馈督察情况时指出, 威宁县县城发展与草海国家级自然保护区保护严重冲突, "城进湖退"问题突出。

由于地处喀斯特生态脆弱区, 草海周边经济、社会发展在受到自然条件制约的同时, 也影响着景观格局及生态脆弱性。该研究以草海为中心, 研究区包括威宁县小海镇黑岩洞至草海一带的集水区域, 涉及草海镇、小海镇和双龙乡部分区域, 共计53个行政村。

1.2 数据来源与处理

遥感影像来源于地理空间数据云平台(http://www.gscloud.cn)和美国地质调查局(http://earthexplorer.usgs.gov)网站提供的Landsat系列的TM、ETM+和OLI传感器类型的影像数据, 轨道号为129/41, 影像以夏季时相为主。基础数据还包括中国科学院资源环境科学数据中心(http://www.resdc.cn)提供的1: 10万比例尺土地利用现状遥感监测数据, 威宁县国土资源局提供的1: 1万地形图和村级行政界线。

为提高遥感影像可视化效果及解译精度, 对所获取的遥感影像进行诸如图像增强、影像数据修复、图像配准、研究区裁剪和数据重采样等预处理。参考GB/T 21010—2017《土地利用现状分类》, 结合威宁县土地利用调查数据库, 将研究区土地利用类型划分为耕地、林地、草地、水域和建设用地, 遥感解译结果经精度验证评价, Kappa系数均大于0.80, 能够满足研究需求。

2 研究方法 2.1 基于景观格局的景观脆弱度构建

景观格局是在自然条件和人类活动的共同作用下形成的, 其形状、大小及排列的不同, 不但可以表现出地物类别的异质性, 而且是生态过程在不同尺度下受到干扰的结果。为研究草海湿地区域景观生态脆弱性及空间分异时空变化特征, 基于景观生态学理论和方法, 以景观格局为基础, 利用景观敏感度指数(LSI, I_LS)、景观适应度指数(LAI, I_LA)构建景观脆弱度指数(LVI, I_LV)。I_LV用于描述景观在受到外界干扰时所表现的敏感性及景观格局的结构、功能和特性发生改变的程度, 其值越大, 表示生态系统越脆弱, 计算公式为

$ {I_{{\rm{LV}}}} = {I_{{\rm{LS}}}} \times \left( {1 - {I_{{\rm{LA}}}}} \right)。$

(1)

I_LA与景观脆弱度密切相关, 景观系统的多样性越大、结构越复杂且分布越均匀, 系统越稳定, 抗干扰适应恢复能力就越强, 计算公式为

$ {I_{{\rm{LA}}}} = {D_{{\rm{PR}}}} \times {I_{{\rm{SHD}}}} \times {I_{{\rm{SHE}}}}。$

(2)

式(2)中, D_PR为斑块丰度密度指数; I_SHD为香农多样性指数; I_SHE为香农均匀性指数。

I_LS表示遭遇干扰时各类景观所受到的生态损失差别, 即其自然属性损失的程度, 是某一类型景观结构指数和脆弱度指数的综合, 计算公式为

$ {I_{{\rm{LS}}}} = \sum\limits_{I = 1}^N {{U_i}} \times {V_i}。$

(3)

式(3)中, n为景观类型数; i为景观类型; U_i为景观类型i的景观干扰度指数; V_i为景观类型i的易损度^[19]。

2.2 空间统计学分析方法 2.2.1 空间自相关分析

将空间自相关分析应用到景观生态脆弱性评价中, 能够发现景观生态脆弱性的空间分布特征^[24]。全局Moran's I指数用于反映空间邻接或邻近区域(行政村)单元属性值(景观生态脆弱性)的相似度, 用于检验某一要素(行政村的景观生态脆弱性)的空间分布模式^[25]。采用全局Moran's I指数表示相邻行政村生态脆弱性的相关程度, 以揭示空间结构性要素对脆弱性的影响程度, 计算公式为

$ \begin{array}{l} I = \frac{{n\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{w_{ij}}} \left( {{x_i} - \bar x} \right)\left( {{x_j} - \bar x} \right)} }}{{\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{w_{ij}}} {{\left( {{x_i} - \bar x} \right)}^2}} }}\\ \;\; = \frac{{\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{w_{ij}}} \left( {{x_i} - \bar x} \right)\left( {{x_j} - \bar x} \right)} }}{{{s^2}\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{w_{ij}}} } }}。\end{array} $

(4)

式(4)中, $ {s^2} = \frac{1}{n}\sum\limits_{i = 1}^n {{{\left( {{x_i} - \bar x} \right)}^2}} ;\bar x = \frac{1}{n}\sum\limits_{i = 1}^n {{x_i}} $; n为研究单元(行政村)个数; x_i、x_j分别为景观生态脆弱性在行政村i和j的观测值; w_ij为空间权重; Moran's I指数的数值范围一般介于[-1, 1]之间, 通常将其解释为相关系数, 即小于0表示呈负相关, 等于0表示不相关, 大于0表示呈正相关。

2.2.2 局部自相关分析

为揭示研究区各行政村单元与相邻行政村之间的关系, 采用局部空间自相关(LISA)分析方法研究53个行政村的景观生态脆弱性空间自相关水平的差异性。LISA(I_LISA)指数可有效反映空间变量值的高值或低值的局部空间聚集, 从而反演局部空间异常特征, 计算公式为

$ {I_{{\rm{LISA}}}} = \frac{{{x_i} - \bar x}}{m}\sum\limits_{j = 1}^n {{w_{ij}}} \left( {{x_i} - \bar x} \right)。$

(5)

式(5)中, $ m = \left( {\sum\limits_{j = 1}^n {x_j^2} } \right)/\left( {n - 1} \right) - \bar x $, I_LISA值为正, 表示该空间单元(行政村)周围相似(高值或低值空间集聚), I_LISA值为负值, 表示非相似值之间的空间集聚^[24]。

3 研究区景观类型及景观脆弱度指数变化特征分析 3.1 研究区景观类型变化特征分析

利用ArcGIS 10.2对遥感解译所获取的研究区1995、2005和2015年3期土地利用数据进行可视化制图(图 1), 并对土地利用景观类型面积进行统计汇总(表 1)。

由图 1和表 1可知, 从研究区整体变化情况来看, 20 a间, 研究区林地、草地景观类型面积整体呈下降趋势。其中, 林地主要分布在研究区周边山区, 林地面积减少113.94 hm²; 草地分布较分散, 主要分布在草海湖区周边, 草地面积减少1 864.35 hm², 草地退化较严重。水域景观主要是草海湖区和杨家湾水库, 20 a间研究区水域面积持续增加, 增加734.13 hm²。随着经济、社会的不断发展, 城镇化水平的加快, 建设用地面积持续增加, 20 a间增加648.72 hm²; 建设用地主要集中在草海保护区毗邻的县城区域, 从空间变化来看, 县城建设用地向保护区方向发展, "城进湖退"问题依然存在。

3.2 研究区景观生态脆弱度指数变化特征分析

采用景观软件Fragstats 4.2对研究区各年份土地利用景观数据进行分析, 获取研究区斑块数量(n_i)、景观类型面积(A_i)、景观类型周长(P_i)、斑块丰度密度指数(D_PR)、香农多样性指数(I_SHD)和香农均匀性指数(I_SHE)等景观格局指数。结合式(1)~(3)和表 1统计结果, 参照文献[26-27]得出研究区景观脆弱度指数(表 2)。

由表 2可知, 1995年研究区景观类型脆弱度指数范围为0.135~0.627, 各景观类型脆弱度指数由大到小依次为草地、林地、耕地、建设用地和水域; 2005年研究区景观类型脆弱度指数范围为0.122~0.953, 各景观类型脆弱度指数由大到小依次为草地、林地、耕地、建设用地和水域; 2015年研究区景观类型脆弱度指数范围为0.121~0.957, 各景观类型脆弱度指数由大到小依次为草地、林地、耕地、建设用地和水域。1995—2015年研究区景观脆弱度指数范围逐渐变大, 各景观类型脆弱度指数排序不变; 其中, 林地、草地景观脆弱度指数逐渐增加, 耕地、建设用地景观脆弱度指数先减小后增加, 水域景观脆弱度指数呈减小趋势。

3.3 土地利用景观类型对研究区景观脆弱性的影响分析

从不同土地利用类型来看, 20 a间草地景观脆弱度大于其他各景观类型, 表明研究区草地景观的稳定性较弱, 对外界干扰较敏感。草地景观多分布在草海湖区周边, 以及研究区周边坡度较大的地带, 易受水土流失、水淹等自然因素和人为因素的影响, 生态被破坏后不易恢复, 导致草地景观脆弱度最大。

林地景观脆弱度次之, 这是由于林地景观类型主要分布在研究区周边的高山地带, 而以草海湖区为中心的周边区域林地景观类型受到人为因素干扰较多。由于受到当地喀斯特自然条件的限制, 为获取更多的耕地, 人为因素导致森林被砍伐, 部分较平缓地区的林地转为耕地, 特别是草海湖区周边区域的低山缓坡地带, 这些地区生态环境极易受损, 导致林地景观脆弱度较高, 其生态问题需进一步得到关注。

研究区内耕地景观以坡耕地为主, 特殊地形以及威宁县城的扩张, 导致部分耕地转为建设用地。耕地区域由于自然条件较好, 人口密集, 极易受到人类活动的影响, 同时人类长时间的改造自然条件, 使得这些区域的生态恢复能力较低, 受到的影响因素也较多。因此, 耕地景观生态脆弱的趋势比较明显, 如果不及时加以治理, 耕地区域的生态环境问题将进一步恶化。

研究区建设用地以威宁县城区域为主, 由于受到当地经济、社会发展等因素影响, 县城建设用地主要在原有基础上向外扩张, 同时受到草海湖区和山区地形限制; 其他建设用地主要分布于山脚下的道路两侧, 具有一定的随意性, 受到干扰较少。水域以草海湖区和杨家湾水库为主, 近年来草海保护区管理委员会开展了一系列以草海湖区水面面积恢复为主的生态保护措施, 使得水域景观较为稳定, 受到人为因素的干扰最小, 脆弱度也最小。

4 研究区村级景观脆弱性变化及时空分异特征分析

通过以上研究, 从整体上分析了研究区景观类型及景观脆弱度变化特征。为具体了解研究区景观生态脆弱性时空变化特征, 使生态保护政策更具针对性和可操作性, 以研究区53个行政村作为评价单元, 从较小尺度和可控单元出发, 利用景观脆弱度评价模型, 运用ArcGIS 10.2融合工具汇总得出各行政村的景观脆弱度指数。

4.1 研究区村级单元景观脆弱性变化特征分析

利用ArcGIS 10.2空间叠加分析功能, 将研究区土地利用数据叠加到行政村多边形, 采用景观分析软件Fragstats 4.2获取各村级景观格局指数, 结合式(1)~(3), 计算得出村级单元景观生态脆弱度。为使不同单位、量级的指标能够进行对比分析, 将景观脆弱度指数进行归一化处理, 其值越接近于0, 表示景观生态脆弱性越差, 越接近于1则越好^[20]。利用自然断点聚类法(natural breaks)将研究区景观生态脆弱性划分为5个等级^[28](表 3), 并进行可视化处理, 得出研究区景观生态脆弱性分布图(图 2)。

由图 2可知, 从空间分布来看, 20 a间处于重度脆弱区、中度脆弱区的行政村主要集中在研究区四周, 例如:1995—2015年, 开华村、凉山村一直处于重度脆弱区, 寒洞村一直处于中度脆弱区; 2005—2015年陕桥村、鸭子塘村一直处于重度脆弱区, 板仓村、石龙村、塔山村一直处于中度脆弱区。进一步分析可知, 处于重度脆弱区的行政村以林地、草地景观类型为主, 而草地、林地景观脆弱度较高, 从而导致这些行政村的脆弱性程度较明显, 易受到人类活动及自然因素的影响, 需引起人们的注意; 研究区中部的行政村以耕地、建设用地和水域为主, 这些区域由于耕地景观分布较广, 面积较大, 其景观脆弱性程度较低, 以潜在脆弱区和轻度脆弱区为主。

从数量上来看, 处于重度脆弱区、中度脆弱区行政村的数量呈减少趋势, 20 a间由9个减少到4个; 处于一般脆弱区的行政村有所增加, 由13个增加到17个; 处于潜在脆弱区的行政村增加较多, 由7个增加到11个; 其中, 位于草海湖区周边的燕山村、阳关村、东山村、西海村和海边村一直处于潜在脆弱区。

结合表 1和图 1进一步分析可知, 处于重度脆弱区、中度脆弱区的行政村主要分布于研究区周边以林地景观类型为主的山区, 近20 a来林地景观遭受损失较严重, 部分林地转为耕地。由表 2可知, 耕地景观脆弱性要小于林地景观, 因此导致处于重度脆弱区、中度脆弱区的行政村数量减少, 部分行政村转为一般脆弱区。增加的潜在脆弱区主要分布于草海湖区周边, 草海湖区周边区域以耕地、草地和建设用地为主, 该部分区域受到当地经济、社会发展等因素影响, 建设用地的增加受到草海湖区和山区地形限制, 受到干扰较少; 近年来草海湖区水面面积得以恢复, 使得水域景观较为稳定。

4.2 研究区村级单元景观生态脆弱性时空分异特征 4.2.1 景观生态脆弱性的全局空间自相关分析

利用ArcGIS 10.2空间统计工具中的空间自相关分析(Moran's I)工具对各行政村的景观生态脆弱性进行统计分析, 得出1995、2005和2015年研究区Moran's I指数分别为0.340 011、0.153 091和0.144 734, 标准化Z得分分别为4.376 702、2.233 296和2.119 562, 标准化Z值均大于2。由此可知, 研究区在不同时间段内各行政村的景观生态脆弱性水平呈现显著自相关, 在空间上存在空间集聚效应, 即对于景观生态脆弱性高的行政村, 其周边相邻行政村的景观生态脆弱性也高, 对于景观生态脆弱性低的行政村, 其周边相邻行政村的景观生态脆弱性也低。20 a间, 研究区Moran's I指数整体呈下降趋势, 表明1995—2015年研究区景观生态脆弱性接近的行政村在空间上表现出的聚集性有所下降。

景观生态脆弱性低的行政村以耕地、建设用地和水域景观为主, 易受到人类活动的影响。随着经济、社会的发展和城镇化进程的加快, 人类干扰活动进一步增加, 建设用地向草海核心区域扩张, 造成草海湖区周边行政村的生态环境受干扰程度加剧, 生态脆弱性维持低水平并呈现集聚趋势; 而景观生态脆弱性高的行政村多处于研究区周边的山区, 林地、草地景观较多, 人类活动干扰相对较弱, 加之近年来退耕还林还草政策的有效实施, 使得景观生态脆弱性高的行政村明显减少。

4.2.2 景观脆弱性的局部空间自相关分析

为揭示研究区53个行政村景观生态脆弱性空间分布的异质性, 利用ArcGIS 10.2空间统计分析工具中的聚类和异常值分析(Anselin Local Moran's I)工具对各行政村景观生态脆弱性进行局部空间自相关分析, 得出各行政村景观生态脆弱性局部空间自相关分布图(图 3)。

经分析可知, 研究区各行政村20 a间的景观生态脆弱性空间分布由低-低聚集类型转为高-高聚集类型; 其中, 1995年海边村、西海村、富民村和前进村为低-低聚集, 2005—2015年鸭子塘村、陕桥村为高-高聚集。20 a间发生的变化以草海中部的行政村为主, 且由低-低聚集转为高-高聚集, 该区域景观生态脆弱性空间分布异质性呈增加趋势。研究区的行政村大部分空间自相关属于"不显著"类型, 并没有出现低-高聚集区。

5 结论

利用3S技术, 通过遥感影像解译获取研究区多时段土地利用景观数据, 基于景观生态学理论, 参考前人研究构建了景观生态脆弱性评价模型。通过对草海湿地区域景观类型及脆弱度指数变化进行分析, 整体了解研究区生态脆弱性变化情况。为有效解决研究区范围过大而缺乏局部指导意义的问题, 从较小尺度和可控单元出发, 以研究区53个行政村作为评价单元, 研究了村级单元景观生态脆弱性时空变化特征, 并进一步采用空间统计学分析了村级单元景观生态脆弱性的空间关联关系。得出以下结论:

(1) 从整体来看, 随着经济社会的发展以及城镇化进程的加快, 草海湿地区域林地、草地面积减少较多, 建设用地面积持续增加; 经草海保护区管理委员会一系列诸如草海湖区水面恢复、退耕还湖、生态保护与修复等政策的实施, 研究区水域面积持续增加, 草海湖区得以恢复, 但建设用地有向草海湖区方向发展趋势, "城进湖退"问题依然存在。

(2) 1995—2015年研究区景观脆弱度指数范围逐渐变大, 但各景观类型脆弱度指数排序不变; 20 a间草海湿地区域景观生态脆弱性空间分布变化较明显, 景观生态脆弱性整体呈下降趋势。

(3) 从村级单元的角度分析可知, 研究区各行政村景观生态脆弱性水平呈现显著自相关, 在空间上存在空间集聚效应, 20 a间景观生态脆弱性接近的行政村在空间上表现出的聚集性有所下降; 各行政村景观生态脆弱性空间分布由低-低聚集类型转为高-高聚集类型, 大部分行政村空间自相关属于"不显著"类型, 并没有出现低-高聚集区。

景观格局能够反映资源和环境的分布形式, 对自然与人类活动所造成的干扰、生物多样性及生态系统稳定性具有一定影响。从景观格局角度出发分析草海湿地景观生态脆弱性时空变化特征及空间关联性和异质性, 具有一定的现实意义。