何国兴， 柳小妮*， 张德罡， 李 强， 蒲小鹏， 刘志刚，关文昊， 杨军银， 韩天虎， 孙 斌， 潘冬荣

(1. 甘肃农业大学草业学院，甘肃 兰州 730070；2. 草业生态系统教育部重点实验室(甘肃农业大学)，甘肃 兰州 730070；3. 甘肃省草原技术推广总站，甘肃 兰州 730000)

作为我国重要的国土资源和可再生性的自然资源，草地具有重要的生态价值和经济价值[1]。FVC不仅可以直观反映草原地上植被生长状况，同时也是指示草地健康状况及监测环境变化的常用评估指标之一[2]。地面实际测量和遥感反演是获取植被覆盖度数据的两种主要方法[3-4]。地面测量通常采用针刺法和目测法[5]，但地面测量具有较强的主观性、效率较低、成本较高等缺点，导致其无法大面积推广[6]。随着遥感技术的不断发展，利用遥感影像反演植被覆盖度已成为较成熟的技术，其主要方法有像元二分法、经验回归模型法和混合像元分解法[7]。由于遥感反演具有效率高、成本低、客观性较强等优点，已成为在大空间尺度和长时间序列范围内获取草地覆盖度数据和监测评估草地健康状况的重要手段[8]。因此，构建遥感反演模型提取区域植被覆盖度的方法不仅是当前研究区域生态模型的基础工作，也是今后开展生态工程评价的主要方向。

中分辨率成像光谱仪(moderate resolution imaging spectroradiometer,MODIS)是地球观测系统(earth observation system,EOS)的核心传感器，它的任务主要是每日获取4次包括大气、海洋、陆地等在内的地球系统相关要素变化数据[9]。游浩妍等[10]基于MODIS数据对研究区的植被覆盖度利用像元二分模型进行估算，发现基于NDVI反演模型拟合效果最佳;周伟等[11]基于MODIS NDVI数据反演了中国草地覆盖度的空间格局和变化趋势，揭示了草地覆盖度变化机理;宋清洁等[12]基于MODIS增强型植被指数(enhanced vegetation index,EVI)和NDVI建立了研究区草地植被覆盖度的反演模型，结果表明基于EVI构建的对数模型为最优反演模型;Wiesmair等[13]使用随机森林回归建立了估算格鲁吉亚高加索草地FVC变化的遥感模型。

甘肃省作为全国重点牧区之一，草原面积2.57×108hm2[14]，其中高寒草甸是一重要类型，具有调节气候、涵养水源、保持水土、净化空气等重要的生态功能[15]。高寒草甸作为全球公认的气候变化敏感区和生态脆弱区，它的形成与青藏高原独特的地形和气候条件有着密不可分的关系[16]。同时，高寒草甸也是调节高寒区草地生态平衡和气候变化的基质，不仅有助于畜牧业的发展，而且对维持生态系统稳定性具有十分重要的作用[17]。目前利用遥感技术获取植被覆盖度数据、监测草地生长状况的研究中大多使用的是基于像元二分法或直接利用NDVI去反映植被覆盖度，但对于某一特定草地类型覆盖度反演的研究较少，不同草地类型对植被指数的敏感程度和模型的精度又不尽相同[2]。因此利用经验回归模型法将地表实测数据和VI相结合构建甘肃省高寒草甸FVC遥感反演模型，可以提高甘肃省高寒草甸FVC的遥感反演精度，对高寒草甸草地生态可持续发展和恢复重建具有十分重要的意义。

本研究基于地面实测样点数据和MODIS产品数据，在研究高寒草甸不同植被覆盖度遥感监测模型基础上，模拟分析2000—2019年甘肃省高寒草甸生长季植被覆盖度的时空变化特征，以期为高寒草甸草地资源的动态监测和健康评价提供理论依据和技术支撑。

1 材料与方法

1.1 研究区概况

以甘肃省高寒草甸草原为研究区(图1)，行政区域包括甘南藏族自治州、天祝藏族自治县、肃南裕固族自治县和山丹县，位于97°23′～104°45′ E，33°06′～39°43′ N。境内地形十分复杂，山地、高原占总面积的95%以上，平均海拔3 000 m，从西北至东南地形逐渐降低。研究区属于高原大陆性气候，年降水量400～800 mm，年均温1～3℃，年日照时数2 200～2 400 h，自西北至东南逐渐降低[12]。研究区高寒草甸面积40 510 km2，占总面积的36%。

高寒草甸是指以冷中生多年生草本植物为优势形成的植物群落，该类草地主要以嵩草属植物为主，草层低矮、丛生，结构简单，层次分化不明显，叶片缩小，被绒毛或呈垫状；主要分布在甘南高原和祁连山山地并向东延伸到海拔3 000 m以上的西秦岭山地。

1.2 数据来源及预处理

1.2.1草地类型数据 将来源于中国草业与生态网络(http://www.ecograss.com.cn/)的《中国1∶100万草地资源图》经过校正后得到本研究中的草地类型数据[18]；甘肃省1∶100万行政区划图来源于地球科学数据共享平台(http://www.geodata.cn)。

图1 研究区及采样点示意图Fig.1 Study area and distribution of sampling points

1.2.2地面实测数据 2014年研究区草原调查数据，来源于甘肃省草原总站，共370个样点(图1)，包括经纬度、海拔高度、草地类型、草地覆盖度等。所有样点的调查时间集中在2014年8—9月，该时期是草地植被生长旺盛期。将样点数据在R 3.6.1中随机抽样70%和30%，分别作为建模数据组和检验数据组。

1.2.3遥感数据 本研究选用美国航空航天局(national aeronautics and space administration,NASA) 对地观测系统数据共享平台(https://lpdaacsvc.cr.usgs.gov/appeears/explore)提供的2000—2019年MOD13Q1产品，包括NDVI和EVI，空间分辨率250 m，时间分辨率16 d。红光波段和近红外波段数据为MOD09 GQ. 006产品，空间分辨率250 m，时间分辨率8 d。

1.3 研究方法

1.3.1植被指数选择 差值植被指数(difference vegetation index,DVI)虽与植被覆盖度相关性较高，但土壤背景较强的影响该指数。通常在较低覆盖度时，RVI随着植被的增加而增加，当植被覆盖度大于80%后，RVI的敏感程度呈现下降趋势[2]。NDVI具有植被探测范围较宽、灵敏度较高的优点，是一种被广泛使用的植被指数，该植被指数可以在不同程度上克服地势、太阳高度角、大气以及群落结构等因素的影响[19]。EVI能分离红光和蓝光通过的气溶胶，补偿气体残留单体对红光的吸收，能够反映植被的季节性[20]。土壤调整植被指数(soil adjusted vegetation index,SAVI)由于加入土壤调节系数L，可以减弱土壤背景的影响。修正土壤调整植被指数(modified soil adjusted vegetation index,MSAVI)也能够有效减弱土壤背景的影响[21]。基于此，选择以上植被指数(表1)，进行研究区高寒草甸FVC反演。

表1 植被指数计算公式Table 1 Vegetation index formula

运用MRT (modis re-projection tool)读取数据，同时转换遥感影像的数据格式和地图投影；利用ENVI 5.3中的Band math工具将MOD13Q1数据和MOD09 GQ.006数据的每个像元值缩小10 000倍，得到真实的NDVI、EVI值和红光波段、近红外波段值[22]。根据计算公式对DVI，RVI，NDVI，EVI，SAVI和MSAVI进行波段组合计算。将植被指数数据采用最大值合成法(maximum value composite,MVC)[23]生成2000—2019年8月份月植被指数，以减少云、积雪、大气和太阳高度角等外界因素的影响。

1.3.2FVC反演模型构建

(1)经验回归模型法 利用ArcGIS 10.2中Spatial Analyst工具提取各经纬度点植被指数值，将FVC实测值与提取的各植被指数在R 3.6.1中进行Pearson相关性分析。选择线性模型、对数模型、指数模型、二次多项式模型和幂函数模型等具有代表性的常用回归模型[16]，构建FVC与各植被指数间的模型，再根据精度结果筛选最佳模型。

式中：FVC为植被覆盖度；NDVIsoil为裸地土壤；NDVIveg为有植被完全覆盖区域的NDVI值。

(3) 模型精度检验采用回归方程的决定系数R2、均方根误差(root mean square error,RMSE)来评价模型精度。RMSE常用来衡量模型的估测能力、描述模型的不确定性以及量化模型精度，RMSE数值越低，表明回归模型精确度越高。

以30%草原调查中的草地FVC为真实值，对所建立的回归模型进行精度评价，分别计算各反演模型估测的草地FVC与真实值之间的RMSE(公式如下)，筛选出研究区不同草地类型草地FVC监测的最佳遥感反演模型。

综合模型的R2，RMSE、实测值与预测值回归方程的R2和斜率等因素，确认最优模型，并应用其反演并分析2000—2019年高寒草甸FVC动态变化。

1.3.3FVC分级 依据全国《草场资源调查技术规程》及植被划分标准[25-27]，并结合甘肃省高寒草甸植被特点，将FVC分为5个等级：低FVC(0～20%)、较低FVC(20～40%)、中FVC(40～60%)、较高FVC(60～80%)、高FVC(80～100%)。

1.3.4FVC趋势分析 采用一元线性回归分析法和趋势法(Slope)进行像元尺度的计算，分析20年间高寒草甸FVC的变化趋势和变化速率。计算公式如下[28]：

式中：θSlope为草地FVC变化趋势；FVCi为第i年的FVC；i为年变量。当θSlope<0时，FVC减少；θSlope>0时，FVC增加。

采用F检验法进一步分析草地FVC随时间发生的变化趋势。计算公式如下[29]：

式中：n为时间序列长度20；i为2000—2019年的年份序号1～20；i为年份序号的平均值；FVCi为第i年的FVC值；FVC为2000—2019年FVC的平均值。

参照已有的研究[29-31]，根据θslope值的大小将研究区高寒草甸FVC变化划分为7个等级：显著增加、中度增加、轻度增加、基本不变、轻度减少、中度减少、显著减少，并统计各区间的面积百分比。

1.3.5植被覆盖度变异系数 变异系数(coefficient of variation,Cv)可以反映高寒草甸FVC多年的变化波动情况，Cv≤0.15说明高寒草甸FVC变化较小，草地生长状况相对稳定；Cv>0.15说明FVC变化较大，草地生长状况不够稳定[32]。本研究通过计算像元尺度的变异系数，进而了解研究区高寒草甸FVC的波动情况。计算公式如下[33]：

参照已有的研究[34]，根据Cv值的大小将高寒草甸FVC稳定性分为5个等级：十分稳定(Cv≤0.10)、稳定(0.100.30)，并统计各区间的面积百分比。

2 结果与分析

2.1 FVC遥感反演模型的构建

2.1.1FVC与各植被指数间相关性 由6种植被指数DVI，RVI，NDVI，EVI，SAVI和MSAVI与FVC实测值的相关性分析(图2)结果可知，各植被指数与FVC均显著正相关(P<0.05)。

其中，NDVI与FVC间的相关系数最高，达到0.90；EVI和SAVI与FVC间的相关性高于DVI和RVI。6种植被指数中RVI与FVC间的相关系数最低(0.59)，其余植被指数与FVC的相关性系数均高于0.65，因此选择DVI，NDVI，EVI，SAVI和MSAVI与FVC进行经验模型的构建。

2.1.2FVC与植被指数间反演模型的建立 基于FVC与DVI，NDVI，EVI，SAVI和MSAVI构建反演模型，选择R2>0.60的模型[35]和像元二分模型计算RMSE(表2)。

5个植被指数模型中，NDVI模型优于其他模型。其中NDVI的指数模型、二项式模型和幂函数模型的R2均高于0.80，而NDVI的二项式模型、幂函数模型和线性模型的RMSE均小于其他模型。EVI的二项式模型(R2=0.64)优于幂函数模型(R2=0.63)。SAVI的对数模型(R2=0.64)和幂函数模型(R2=0.68)优于MSAVI的二项式模型(R2=0.61)。

2.1.3FVC与植被指数间反演模型精度检验 NDVI的线性模型、指数模型、对数模型、二项式模型和幂函数模型，以及EVI的二项式模型，SAVI的对数模型和幂函数模型、像元二分模型共9个模型的检验结果如图3所示。

图2 Pearson相关性矩阵Fig.2 Pearson correlation matrix注：相关性矩阵分为两部分，上三角为相关系数与显著性检验，下三角为散点图与拟合曲线，**与***代表在0.05和0.01水平上差异显著。各变量间相关关系相交于变量的行列处Note:The correlation matrix is divided into two parts,the upper triangle is the correlation coefficient and the significance test,the lower triangle is the scatter chart and the fitted curve,** and *** indicate significant difference at the 0.05 and 0.01 level,respectively. The correlations between variables can be found at the row and column intersection of the two variables

表2 植被覆盖度与各植被指数间的反演模型Table 2 Inversion model of vegetation fractional cover and vegetation indices

图3 FVC实测值与预测值比较Fig.3 Comparison of measured and predicted values of FVC注：(a)NDVI线性模型；(b)NDVI指数模型；(c) NDVI对数模型；(d) NDVI二项式模型；(e) NDVI幂函数模型；(f) EVI二项式模型；(g) SAVI对数模型；(h) SAVI幂函数模型；(i)像元二分法。图中虚线为预测值等于真实值，实线为预测值与真实值间的回归直线，实线离虚线越近说明预测结果越准确Note:(a) NDVI Linear model;(b) NDVI Exponential model;(c) NDVI Logarithmic model;(d) NDVI Quadratic polynomial model;(e) NDVI Power function model;(f) EVI Quadratic polynomial model;(g) SAVI Logarithmic model;(h) SAVI Power function model;(i) Pixel dichotomy. In the figure the dashed line is the predicted value equal to the measured value. The solid line is the regression line between the predicted value and the measured value. The closer the solid line is to the dotted line,more accurate the prediction

NDVI二项式模型的拟合R2=0.89，高于其他模型。像元二分模型线性拟合方程的斜率为0.90，高于其他模型。

综合模型R2，RMSE，实测值与预测值回归方程的R2和斜率等因素，确认NDVI的二项式模型：y=—0.65x2+1.97x—0.23为最优模型。

2.2 2000—2019年高寒草甸FVC的时空变化

2.2.1FVC的空间分布特征 应用NDVI二项式模型y=—0.65x2+1.97x—0.23反演2000—2019年高寒草甸FVC的动态变化。由20年高寒草甸均值空间分布(图4a)可知，高寒草甸FVC整体较高，其中FVC>80%的面积最大(47.44%)，60

从区域分布来看，研究区南北跨越较大，高寒草甸FVC呈现明显变化趋势。玛曲县、碌曲县、夏河县、合作市和天祝藏族自治县等海拔较高且人烟稀少的区域，高寒草甸FVC较高(>80%)；60

2.2.2FVC的时间变化规律 由高寒草甸FVC年际变化图(图4b)中可以看出，2000—2019年高寒草甸年均FVC呈现波动上升趋势(P<0.05)，分别在2001年、2007年和2012年FVC出现峰值，随后FVC下降。线性拟合发现，高寒草甸FVC以0.15%·a-1(R2=0.58)的速率呈显著上升趋势(P<0.05)。从整体变化来看，2002—2007年、2008—2013年和2015—2019年这3个时间段内FVC均呈现上升趋势。高寒草甸FVC大约每5年会出现明显波动。

将高寒草甸2000—2019年每年FVC等级进行统计(图4c)后发现，高寒草甸年均FVC大部分以中、高覆盖度为主，说明高寒草甸整体上植被覆盖度较高，生长状况良好。2000年极高覆盖度水平比例最小，极低和低覆盖度水平比例达到研究年份的最大。其余年份内低和极低覆盖度水平的比例基本趋于稳定，高和极高覆盖度水平的比例呈现此消彼长的变化趋势。

图4 2000—2019年高寒草甸FVC空间分布及年际变化趋势Fig.4 Alpine meadow FVC spatial distribution and annual change trends from 2000 to 2019

2.2.3FVC空间趋势变化分析 趋势分析结果表明(图5)：2000—2019年FVC呈上升趋势的面积为21 078 km2，占高寒草甸总面积的52.76%；稳定趋势的面积为11 017 km2，面积占27.58%；下降趋势的面积为7 852 km2，面积占19.66%。其中轻度增加的面积比例最高，占33.14%；显著减少所占的面积比例最低，占0.82%。

图5 2000—2019年高寒草甸FVC变化趋势分级Fig.5 Alpine meadow FVC change trends classification from 2000 to 2019

从空间分布来看，FVC显著增加的区域主要分布在肃南裕固族自治县；中度和轻度增加的区域主要分布在夏河县、合作市、天祝藏族自治县和肃南裕固族自治县；轻度和中度减少的区域主要分布在碌曲县的中部和玛曲县的东南部；显著减少的区域零散分布于玛曲县的边缘和碌曲县的中部。

2.2.4FVC稳定性分析 研究区2000—2019年高寒草甸Cv≤0.15的面积占总面积的86.48%，高寒草甸FVC这20年内波动较小；Cv<0.10的面积占49.59%；0.100.15)区域主要分布在肃南裕固族自治县、民乐县、天祝藏族自治县、夏河县、玛曲县和舟曲县的边缘地带，其面积占高寒草甸总面积的13.52%；高波动区域(Cv>0.30)主要分布在肃南裕固族自治县的西北部，总面积1 095 km2(图6)。

图6 2000—2019年高寒草甸FVC变异系数分级Fig.6 Alpine meadow FVC coefficient of variation classification from 2000 to 2019

3 讨论

3.1 高寒草甸FVC遥感反演模型

植被指数是将遥感数据相关的光谱信号综合在一起，通过不同波段的组合，增强不同的植被信息，降低干扰[36]。植被指数受到植被自身、大气、地物等因素的影响，能够反映绿色植物的相对丰度，与FVC和叶面积指数(leaf area index,LAI)等指标具有较强的相关性[37]。然而，植被指数不具备物理意义，只有将其转换为具体的生态参数(生物量、植被覆盖度、叶面积指数等)后，植被指数才能够定量分析生态环境的变化过程[38]。

本研究发现，NDVI和SAVI与FVC的相关性高于其他植被指数，而RVI与FVC间的相关性低于0.60，可能是由于高寒草甸地势复杂、日照强烈、紫外线强、潮湿寒冷、昼夜温差大、植被生长期短等因素导致[39]。利用像元二分法和经验回归模型法对高寒草甸植被覆盖度构建的反演模型，其中经验回归模型法NDVI和SAVI两个植被指数与植被覆盖度的相关性较强，NDVI二次多项式回归模型精度验证拟合效果达到0.89，能够较好地反演高寒草甸植被覆盖度，这与巴音布鲁克高寒草原遥感反演模型研究的结果相一致[2]。另外，陈建军等[38]研究也发现，NDVI二次多项式模型对黄河源区高寒草地FVC具有较高的反演精度。尽管部分研究表明NDVI会出现饱和现象，但由于高寒草甸是以嵩草为主的草甸草原，单位面积内的生物量有限、草层较低导致NDVI的饱和现象不明显[38]。像元二分模型的反演精度略低于NDVI二项式模型，这与李艺梦等[40]基于Landsat8影像植被覆盖度估算的研究结果有相似之处。

3.2 高寒草甸FVC趋势变化

本研究分析结果表明，2000-2019年高寒草甸年均FVC呈现波动上升趋势，高寒草甸FVC呈现西北低、东南高的分布格局。另外，20年平均FVC显示，高寒草甸FVC整体呈现较高水平，其中高FVC区域主要在甘南州地区，这与王伟军等[41]在甘南州草地NDVI动态变化的研究结果相一致。周伟等[42]利用MODIS数据研究中国草地覆盖度时空动态时也发现，甘南州整体FVC处于较高水平。趋势分析结果表明，20年间高寒草甸FVC呈增加趋势的面积高于降低趋势的面积。其中呈减少趋势的地区主要在玛曲县的东部和碌曲县的中部，这与杨瑞瑞[43]在黄河源植被覆盖度的研究结果相一致。其他地区均呈现增加趋势。变异系数结果表明，高寒草甸FVC相对稳定的比例占86.57%，说明20年间高寒草甸FVC整体呈现稳定趋势。

究其原因，是因为甘肃省高寒草甸大部分位于甘南州境内，而甘南州又地处青藏高原东北边缘，是青藏高原与黄土高原的交汇区，在气候上属于亚洲季风边缘区域和高寒气候的过渡区，属于华西秋雨区[44]。秋季连阴雨对甘南高寒草甸草原的返青和枯黄具有十分密切的关系[45]。秋季较高的降水量是该区域高寒草甸FVC高于其他区域FVC主要原因，适宜的水份直接影响植物的光合作用和呼吸作用，是决定植物生长状况的关键因素。由于甘南州高寒草甸面积大、分布广，其FVC的水平直接影响甘南高寒草甸的整体水平，是导致高寒草甸FVC较高的一个重要原因。加之上世纪60年代开始的全球气候变暖为高寒草甸植被生长提供了良好的温度条件，由于温度和水分均处于较适宜的条件，所以高寒草甸FVC呈现逐年增加的变化趋势。另外，高寒草甸FVC的增加与我国生态恢复工程的不断实施也有密切的关系。如1999年和2003年相继开始实施的退耕还林还草和退牧还草工程，对我国草地植被恢复和牧草生长具有重要的促进作用[46]，使得2002—2007年间甘肃省高寒草甸FVC呈现上升趋势。同样，2010和2015年相继实施的二轮生态补奖政策对草地生态保护起到了重要作用，也促进了草地FVC的稳定增长[46]。甘肃省高寒草甸区位于青藏高原东段，是两轮生态补奖中的重点补偿区域，2010年以后高寒草甸FVC出现拐点，呈上升趋势，尤其是2015—2019年高寒草甸FVC呈现明显的上升趋势，这与生态工程的实施有着密不可分的关系。

但对高寒草甸这一特殊地理环境下的独特植被研究较少，对于该草地类型的时空变化的具体原因，有待进一步研究。

4 结论

本研究基于MODIS数据和高寒草甸FVC实测数据，采用Person相关性分析后发现DVI，NDVI，EVI，SAVI和MSAVI与高寒草甸FVC具有较高相关性。NDVI的二次多项式模型：y=—0.65x2+1.97x—0.23能较好的反演高寒草甸FVC，模型精度较高。2000—2019年，甘肃省高寒草甸年均FVC时间上呈现波动上升的趋势，2002—2007年、2008—2013年和2015—2019年这3个时间段内FVC均呈现上升趋势，大约每5年会出现明显波动。空间上高寒草甸FVC由东南向西北逐渐递减，其中FVC较高区域主要分布在甘南藏族自治州和天祝藏族自治县，且具有较强的稳定性；FVC较低的主要集中于研究区的西北部且稳定性较差。