近35 a新疆天山巴音布鲁克草原退化程度评价
2023-05-29邓成军李文利公延明李凯辉
赵 剑,邓成军,李文利,赵 金,公延明,2,李凯辉,2
(1.中国科学院新疆生态与地理研究所,新疆 乌鲁木齐 830011;2.中国科学院新疆生态与地理研究所巴音布鲁克草原生态系统研究站,新疆 巴音布鲁克 841314;3.中国科学院大学,北京 100049;4.新疆巴音郭楞蒙古自治州草原工作站,新疆 库尔勒 841000)
草原是我国陆地生态系统的重要组成部分,占国土面积的41.67%[1],其中各类天然草地约有4×108hm2,可利用面积为3.9×108hm2,集中分布在中西部地区[2]。新疆是我国五大草原畜牧业省份之一,天然草地毛面积5.6×107hm2,可利用面积4.8×107hm2,位居全国第三位,草地面积约占新疆总面积的34.44%[3]。然而,受自然和人为因素影响,我国草原正发生不同程度的退化[1]。新疆天然草地已有85%出现不同程度的退化[4],与20 世纪60 年代相比,巴音布鲁克草原草地植被覆盖度由原来89.4%下降到30%~50%[5],巴音布鲁克草原退化面积达到30.48×104hm2,占草地总面积的19.61%[5]。2000 年以来,我国实施了包括退牧还草和草原生态奖补机制在内的一系列重大生态工程[6-7],重点进行退化草原修复治理工作。
草原幅员辽阔,仅采用野外调查手段不能满足大尺度评价要求,且监测周期长工作量大,因此3S技术被广泛应用到草地监测方面[8]。草原生态系统研究中,学者们大多选用植被指数、植被覆盖度、植被生产力等指标来评价草地状况。吴晓全等[9]利用CASA 模型估算了新疆天山地区的植被净初级生产力的时空分布并分析其驱动因素,结果表明天山地区植被-净初级生产力空间上呈西高东低、由北向南递减的特征;针对草原退化或健康状况评价问题,杜自强等[10]选取植被覆盖度、牧草地上生物量和牧草可食率3个指标加权平均构建了草地退化指数,通过TM 影像遥感监测黑河中上游草地1986—2003 年的退化情况,结果表明草地退化加剧、局部改善、整体恶化的变化格局;赵玉婷等[11]基于“压力-状态-响应”(PSR)模型,从压力、状态、响应3个层面选取17 个指标构建了甘南藏族自治州高寒牧区草地生态系统健康评价体系,并表明13 a来该地区高寒草地健康状况负向变化趋势明显;陆均等[12]通过CVOR 指数研究表明,巴音布鲁克高寒草原生态系统健康状况由2004 年的一般病态逐渐恢复为2012年健康状态,且退化状况由强度退化恢复为中度退化;吴志丰等[13]探究了一种基于遥感的参照覆盖度提取及草地退化评价研究的方法,在一定程度上解决了目前大范围草地退化评估过程中面临的缺乏参照系统及遥感数据应用误区等问题;李霞等[14]选取盖度、总产量和可食草产量3个指标并结合GIS中的空间插值法反演甘肃省草地退化状况,结果显示甘肃省退化草原面积1.79×107hm2,占全省草原总面积的69.65%。草原退化程度反映在植被覆盖度、地上生物量、可食草比例、草地退化指示物种以及土壤有机质含量等指标的变化上,诸多研究对不同指标进行了分析,但退化表征指标复杂且受限于各种调查技术,各个分指标相互独立,综合评估指标的建立是面临的问题。
目前,对草原退化程度的分级方法种类繁多,各退化评价指标差异大且难统一,缺乏综合评价指标体系,通过计算分指标的权重值,实现不同指标的有机综合是提高草原监测精度的关键[15]。因此,本研究选取植被覆盖度、总产草量以及草层平均高度3 个指标,并引入了Min-Max 标准化方法对分指标进行处理,进一步构建草地退化指数(Grassland Degradation Index,GDI),以1986年野外调查数据作为未退化参照,通过中高分辨率的Landsat影像反演并计算得出GDI 的变化率来表征草原退化情况,旨在建立草原退化遥感综合评估指标,解决草原退化大尺度监测效率低的问题,评估20 世纪80 年代以来巴音布鲁克草原的退化情况。
1 材料与方法
1.1 研究区概况
巴音布鲁克草原位于新疆维吾尔族自治区巴音郭楞蒙古自治州和静县西北天山南麓,东西长292 km,南北宽108 km,主要由大尤尔都斯盆地和小尤尔都斯盆地以及天山山脉中段组成,海拔为1624~4606 m。其气候为典型高寒气候类型,雨雪天气多,冬季漫长夏季短暂,年均气温-4.8 ℃,全年积雪日约137 d,积雪深度45 cm[16]。草地类型以沼泽类、高寒草甸类和高寒草原类为主,地表优势植物有线叶嵩草(Kobresia capillifolia)、珠芽蓼(Polygonum viviparum)、紫花针茅(Stipa purpurea)、寒生羊茅(Festuca kryloviana)等。
地面调查样地点数据共52个,其中建模集样点35个,验证集样点17个(图1)。
图1 研究区位置及样地点位Fig.1 Location and sample points of the study area
1.2 数据来源与处理
本文采用与研究相关的遥感数据和地面调查数据。其中20世纪80年代草原普查数据和2021年地面样地调查数据均来自于中国科学院巴音布鲁克草原生态系统研究站;选取的遥感影像为2000年7 月6 日Landsat ETM+、2009 年7 月31 日Landsat ETM+和2021 年7 月24 日Landsat OLI 共3 期影像,云量低于20%且与实地采样时间最为接近,数据间隔控制在10 a 左右,保持一定的逻辑;数字高程DEM 数据来自于地理空间数据云(https://www.gscloud.cn)。
预处理:
(1)遥感影像(ENVI 5.3软件):Radiometric Calibration、QUAC、Seamless Mosaic、边界裁剪。
(2)地面调查数据(ArcGIS 10.6):Add Data 工具添加坐标点,Extraction Values to Points模型构建。
(3)数字高程数据(ArcGIS 10.6):拼接、边界裁剪和格式转换。
1.3 统计学方法
1.3.1 回归分析 回归分析(Regression)是一种试图以一个或多个自变量来解释另一个因变量的统计方法,在非线性回归模型中,一种类型是可以通过变量变换,化为线性模型进行分析,即曲线回归(曲线拟合),一般形式为:
式中:yi为因变量(被解释变量);f(x,θ)为参数估计方程;ei为误差的函数,拟合后取最小值。
1.3.2 Min-Max 标准化 由于不同的数据有不同的量纲,而进行综合指标构建时需要有统一的标准,因此需要对数据首先进行标准化处理。数据标准化可以消除不同数据之间不同量纲的影响,解决数据不匹配的问题,常用的标准化方法有Min-Max 标准化(Min-Max Normalization)[17]、z-score标准化方法等,考虑数据量和数据处理的便捷性,本研究采用的方法为Min-Max 标准化方法,标准化处理之后数据在[-1,1]范围内。
式中:X*为归一化处理后的数据;X为原始数据;Xmin为数据的最小值;Xmax为数据的最大值。
1.3.3 主成分分析法 指标权重常用的确定方法有专家打分法[14]、主成分分析法[18]、聚类分析法[19]等。但专家打分法掺杂了太多的主观意识,主成分分析(Principal Component Analysis,PCA),是一种通过正交变换将一组存在相关性的变量转换为一组线性不相关变量的数学统计方法,将数据的大部分信息以贡献程度靠前的变量集中表示。已有多项学者的研究表明,主成分分析法应用于草地退化分指标的权重确定具有一定的科学性和准确性[20-21]。
1.3.4 精度评价 GDI 模型的精度采用验证点位的遥感拟合值与地面点位实际计算值的均方根误差结果来进行衡量,均方根误差(RMSE)是用来衡量观测值同真值之间的偏差,RMSE 对相差较大的数据相关性更高,所以其针对异常值更为敏感,计算公式为:
式中:RMSE 表示均方根误差;n表示样本的个数;GDIg表示地面调查点位实际计算值;GDIrs表示遥感拟合值。在用于评价反演模型的精度时,RMSE值越小表示模型的精度越好。
1.4 草原退化分级方法
1.4.1 植被指数 归一化植被指数(Normalized Difference Vegetation Index,NDVI)是应用最广泛的植被指数,可以有效反映植被状况,且计算简单,研究也最深入[22]。此外,本文又引入了其他植被指数,差值植被指数(Difference Vegetation Index,DVI)、比值植被指数(Ratio Vegetation Index,RVI)、增强型植被指数(Enhanced vegetation index,EVI)、土壤植被指数(Soil-Adjusted Vegetation Index,SAVI)(表1),并给出其计算公式。
表1 植被指数及计算公式Tab.1 Vegetation index and calculation formula
1.4.2 草地退化指数(GDI)草地退化指数(Grassland Degradation Index,GDI)是大尺度上用于评价草原退化程度的常用指标[15]。本文草地退化指数是植被覆盖度、总产草量和草层平均高度3 个因子的加权综合,是草原植被状况的数字化表示,取值范围为0~1,计算公式为:
式中:GDIg表示地面调查点位草地退化指数;n表示指标的个数;vi表示i指标的标准化值;wi表示i指标的分权重。
1.4.3 草原退化分级标准 以1986 年草原样方调查资料作为未退化参照标准,根据《天然草地退化、沙化、盐渍化的分级指标》(GB 19377-2003)[29]和前人的研究[21],参考了植被覆盖度减少率(未退化<10%、轻度退化11%~20%、中度退化21%~30%、重度退化>30%),草层平均高度减少率(未退化<10%、轻度退化11%~20%、中度退化21%~50%、重度退化>50%)和总产草量减少率(未退化<10%、轻度退化11%~20%、中度退化21%~50%、重度退化>50%)的单指标分级阈值,加权平均得到了综合指标草地退化指数的分级阈值,将研究区草原退化程度分为4 个等级:未退化(GDI 的减少率<10%)、轻度退化(GDI的减少率为10%~20%)、中度退化(GDI的减少率为20%~40.8%)、重度退化(GDI 的减少率>40.8%)。
2 结果与分析
2.1 20世纪80年代草地退化指数的计算结果
草原退化的参考基准值根据20世纪80年代草原普查数据确定,新疆天山巴音布鲁克草原分为高寒草甸类、高寒草原类、沼泽草地类3 种草地类,地表植被包含19种植物(表2)。
表2 主要植物种类Tab.2 Main plant species
结合遥感影像和野外实测样地精细划分29 个草地类型,计算20 世纪80 年代草原普查数据形成了退化背景的GDIg参考值(表3)。GDIg最大值出现在高寒草甸类的线叶嵩草、杂类草类型,表明该类型植被状况最好,GDIg最小值出现在高寒草原类寒生羊茅、冰草,该类型植被状态最差,整体来看,高寒草甸类和沼泽类植被要好于高寒草原类。
表3 20世纪80年代草原普查数据集Tab.3 Grassland census data set in 1980s
2.2 2021年草地退化指数的计算结果
2021年地面调查点位的70%作为建模集样本,共35 个样地点位,即高寒草甸类样本点6 个,高寒草原类样本点24个,沼泽类样本点5个(表4),采样时间为7月10—31日,植被覆盖度主要通过目视估测法得到,每个样地包含3个1 m×1 m的小样方,通过刈割法获取总产草量鲜重,草层平均高度是样地内多个物种的高度平均值。
表4 2021年建模样地数据集Tab.4 Data set of construction points in 2021
2.3 草地退化指数的遥感模型
将归一化处理后的植被覆盖度、草层平均高度、总产草量3个变量通过SPSS提取出1个主成分,通过成分矩阵和分指标的方差贡献率计算得到综合得分模型中的系数(表5),最后对其进行归一化处理和求平均值。
表5 指标权重系数Tab.5 Index weight coefficient
植被覆盖度、草层平均高度、总产草量的权重分别为:0.46、0.30 和0.24(表5)。草地退化指数GDIg计算公式如下:
式中:GDIg表示地面调查点位的草地退化指数;Ci表示i群落标准化的植被覆盖度;Hi表示i群落标准化的草层平均高度;Yi表示i群落标准化的总产草量。
建模样地GDIg与6 种植被指数相关性最高的为NDVI,皮尔逊相关性系数为0.743,sig<0.01,在0.01 水平上相关性显著,优于DVI 的0.694**、RVI的0.568**、EVI 的0.642**、SAVI 的0.637**及TGDVI的-0.435**,其中,*代表在0.05 级别(双尾)相关性显著,**在0.01 级别(双尾)相关性显著。在GDIg为因变量NDVI 为自变量的5 种回归模型(线性、指数、对数、二次、幂函数)中,线性函数模型显示出了最佳拟合(R2=0.553)。草地退化指数GDIrs的计算方程如下所示:
式中:GDIrs表示草地退化指数的遥感拟合值;NDVI2021表示2021年的NDVI值。
为了测试GDIrs回归模型遥感拟合值的准确性,采用剩余30%的样地点位(17个样地点位)用于模型精度验证,GDIg与GDIrs之间RMSE 为0.140,RMSE值较小,说明反演模型精度较高,可用于草原退化程度评价。
2.4 巴音布鲁克草原的退化现状
相比于20 世纪80 年代,2021 年巴音布鲁克草原退化分级以未退化为主,未退化面积占总面积的比例为60.51%,退化面积占总面积的比例为39.49%,其中重度退化面积、中度退化面积和轻度退化面积分别占总面积的比例为14.71%、12.55%和12.23%。
在3 种草地类型中,高寒草甸类草地面积最大。高寒草甸类重度退化面积占该类草地总面积的比例为21.19%,大于高寒草原类的0.67%和沼泽草地类的14.39%,高寒草原类未退化面积占该类草地总面积的91.54%,大于高寒草甸类的49.06%和沼泽草地类的41.61%。高寒草甸类草地退化现状比高寒草原类和沼泽类更为严重,高寒草原类草地植被退化状况最轻。
由于巴音布鲁克草原山间盆地地形特殊,海拔相差大,草地主要集中分布在2347~2663 m、2663~2980 m、2980~3322 m 的3 个海拔梯度内(图2)。1624~2347 m、3698~4606 m海拔梯度内重度退化草地面积占该梯度草地总面积的比例分别为63.03%、57.32%,山地和平原地区草地退化现状严重;2347~2663 m、2663~2980 m 海拔梯度内未退化草地面积占该梯度草地总面积的比例为75.03%、73.40%,中高海拔地区尤其是盆地内草地退化状况较轻。
图2 不同海拔梯度退化面积Fig.2 Degradation area map of different altitude gradients
2.5 GDI模型的草原退化评价应用
GDIrs模型可以应用于前期年份的草原退化程度评价,由于缺乏历史年份的地面调查数据,又不能将2021 年GDIrs模型直接应用于其他年份,故需要通过2021 年NDVI 值对2000 年、2009 年的NDVI值做辐射配准,提取2021 年52 个样地点位对应的2000 年、2009 年 和2021 年 的NDVI 值,通 过SPSS 25.0完成回归分析获得辐射配准模型。
根据公式(7)、公式(8)对2000年和2009年NDVI 值进行波段运算,将计算结果代入公式(6)计算2 a 的GDIrs值,最后通过GDIrs的变化率分级确定2000年和2009年草原退化程度,评价近20 a的草原退化变化特征。
通过2000年、2009年和2021年3期退化结果明显看出(表6),2000年以来,巴音布鲁克退化情况有了一定的改善,2000—2009 年重度退化、中度退化和轻度退化面积占总面积的比例分别下降了5.58%、15.84%、4.4%,未退化面积占总面积的比例增加了25.82%,退化情况得到了明显的改善;2009—2021 年变化幅度不大,重度退化、轻度退化面积占总面积的比例提高了3.82%、0.45%,中度退化的面积占总面积的比例下降了0.27%,未退化面积占总面积的比例降低了4%,草原退化程度有轻微波动,但未退化面积仍大于退化面积。
表6 2000年、2009年和2021年草原退化情况Tab.6 Vegetation degradation in 2000,2009 and 2021
3 讨论
巴音布鲁克草原海拔较高,通过遥感手段监测时难以找到无云干扰的多期影像,长时间序列研究比较困难,本文尽量保证了每10 a 评估1 次的逻辑性。与其他学者草地退化遥感评估的区别在于本文形成了20世纪80年代草原背景值,以GDIrs的变化率来表征了退化程度。目前存在3点问题需要深入讨论:
3.1 草原退化遥感评估的误差分析
本研究中GDIrs模型拟合度为0.553,相关性一般,主要有3个方面的原因导致:(1)由于巴音布鲁克遥感影像云量较大,选择和地面采样点时间最为接近的Landsat 数据,二者之间仍然存在时间异质性,时间上不能做到完全匹配;(2)Landsat 影像分辨率为30 m,但样方为1 m×1 m,空间异质性会导致草原退化评估结果产生误差;(3)样方点位过少,过拟合现象会导致评估结果产生误差。针对草原退化的本底参照系统,以评价指标的最大值[30]或平均值作为理想样地参数也是一种重要思路,但会有退化评估程度过低的问题,还有学者选取起始年份作为背景值[31],但起始年份草原的退化情况难以准确把握。为了更接近草原未退化实际情况,本研究的巴音布鲁克参照系统按照20世纪80年代草地类型大类赋值,将退化背景值控制在合适范围内。
3.2 GDI模型草原退化评估的适用性
目前,应用于草地生态系统的遥感应用研究中,存在定量反演精度不高的问题[32],本文构建的GDI 模型可以解决草原退化定量评估的科学问题;并且通过辐射配准可以完成GDI 模型在前期年份的应用,填补由于历史数据缺失导致连续时间序列研究的空白,时间尺度上模型更具优势,可以用于指导草原生态系统保护的数据基础,对于长时间序列草原退化评估该模型普适性更高。并且通过引入权重可以实现多指标的综合,提高遥感草原退化评估的效率,满足草原生态系统管理的应用需要。
2000年未退化面积占总面积的比例为38.69%,2009年未退化面积占总面积的比例为64.51%,2021年未退化面积占总面积的比例为60.51%,本文结果与焦阿永等[33]通过NDVI值计算得出的2000—2019年北疆地区草地覆盖显著增加的结论相符;但其中2000—2009年草原退化程度相差较大,新疆草地植被-净初级生产力与降水呈正相关关系[34],通过辐射配准方法完成前期年份的草原退化评估存在气象背景值不统一的误差,优先完成研究区不同降水年份的分类,可以提高模型的精度[35],后期工作需要把巴音布鲁克草原分为降水平年、涝年和旱年,比较相同降水背景下的草原变化情况。
3.3 草原退化的影响因素
植被长势好坏与气候因素密切相关,Fu等[36]研究显示,1961—2009年新疆天山开都河流域气温以0.167 ℃·(10a)-1的趋势增加,降水以0.167 ℃·(10a)-1的趋势增加;刘濛濛等[37]对巴音布鲁克气候的研究结果表明,20 世纪90 年代中期以后,巴音布鲁克气温升高、降水增多,因此本文研究时间段内的1986—2015年巴音布鲁克草原气候总体呈现“暖湿化”;Gang 等[34]研究表明,新疆1981—2010 年的草地NPP 同降水和气温呈正相关关系,巴音布鲁克草原温度回暖和降水增加趋势使得当地草地植被得以恢复,这与赵鹏等[38]研究发现,1982—2015 年气候要素导致新疆草地NPP增加的结果相似,与本文研究得出2000—2009 年草原退化程度明显改善的结果相符,而对于2009—2021年草原退化程度轻微波动的结果,姚俊强等[39]研究结果指出,20 世纪90年代之后新疆气候多呈暖湿配置,但在2022年时段内有所减缓,“暖湿化”趋势减缓导致了2021年草原植被状态的非线性增加结果。此外,人类活动对草原的影响复杂多变,曹永香等[40]通过主成分分析法研究发现,人类活动是策勒绿洲-沙漠过渡带植被变化的主导影响因素。除了牲畜量的增加导致草原退化外,退牧还草工程下的围栏封育、季节性休牧、划区轮牧、牧草补播等[41]是目前人类保护草原生态的主要措施,人类活动在草原退化程度改善中起到了多大的作用需要选取合适的科学方法进一步分析。
4 结论
本文通过构建GDIrs模型和Landsat影像对巴音布鲁克草原多年退化程度进行评价,提高了草原退化遥感监测效率,主要结论有:
(1)GDIg综合指标为植被覆盖度、草层平均高度、总产草量的加权平均,引入该指标可以提高草原退化遥感评价的效率;GDIg与NDVI 的相关性最好,且GDIrs模型方程的拟合度为0.553,RMSE 为0.140,适用于巴音布鲁克草原退化程度遥感评估。
(2)2021年巴音布鲁克草原未退化面积占总面积的比例为60.51%,高寒草甸类草地退化最严重,高寒草原类草地退化状况最轻,空间分布上表现为由盆地向山地退化加重的趋势。
(3)利用辐射配准方法对前期年份草原退化程度进行评估,研究的时间尺度更长,结果表明2000—2009年巴音布鲁克草原退化程度明显改善,2009—2021年草原退化程度轻微波动。