基于MODIS数据的黄河源区植被覆盖度反演研究

刘 晶，党星海,2，陈丽丽，边雁君，黄立鑫，李志红，周兆叶*

（1.兰州理工大学 土木工程学院，甘肃 兰州 730050；2.甘肃省应急测绘工程研究中心，甘肃 兰州 730050；3.中铁二十局集团市政工程有限公司，甘肃 兰州 730030）

高寒草地是维持高原地区生态平衡的关键物种，高寒草地的退化将导致不可逆转的生态问题，因此对于高寒草地的监测至关重要。高寒草地的监测方法主要包括传统的地面调查和利用植被指数进行遥感估 测[1]，地面调查方法费时费力，样本数量少，空间代表性差，得出的结论不足以很好地反映高寒草地覆盖度的真实空间分布状态[2-3]；遥感估测的常用手段包括卫星遥感、无人机遥感等[4]。利用无人机能低空飞行、悬停等特点进行取样[5-6]，并配合少量地面取样，可成为样方尺度高寒草地覆盖度、高度、生物量调查的主要手段，在很大程度上解决传统观测方法的不足[7-8]。 无人机操作灵活，可通过设置其飞行高度和飞行模式获取与遥感像元尺度匹配的图像，能架起样地实测数据与遥感影像之间的桥梁。目前基于植被指数的遥感估测方法包括回归方程法、像元二分法和机器算 法等[1,7]。

本文以黄河源区的高寒草地为研究对象，以无人机航拍为获取数据的主要手段，结合地面样方调查数据和MODIS数据提取了植被指数；再利用回归模型法和像元二分法反演了源区高寒草地覆盖度；最后筛选最优指数和算法进行黄河源区高寒草地覆盖度估算以及变化趋势分析，以期为高寒草地的合理利用与生态治理提供数据支撑。

1 数据来源与研究方法

1.1 数据来源

1.1.1 研究区概况与数据获取

本文所指的黄河源区为唐乃亥水文站以上区域，面积为1.2×105km2，主要植被类型为高寒草地、高寒灌丛化草地和高寒草甸等[8-9]。该区域植被生长季节性变化明显，对气候变化极其敏感，属于生态脆弱带和气候变化敏感区[10-12]。

本文采用大疆精灵Pro四旋翼无人机，航高设置为20 m，样地大小布设为250 m×250 m（以便与MODIS数据像元大小匹配），样地内以矩形的方式获取照片，每块样地共获取13张照片，单张照片覆盖范围约为 35 m×35 m[3]。2016－2018年7月底8月初在黄河源区共拍摄326块样地，4 238张照片。样地位置分布如图1所示。

图1 研究区航拍样地位置分布

1.1.2 MODIS数据下载

根据黄河源区的经纬度范围，确定行列号为h25v05、h26v05的两景数据可覆盖全区域。该数据的时间分辨率为16 d，空间分辨率为500 m。本文下载数据时间序列为2001－2018年7月11日（以便与野外航拍时间一致）的MOD11A1和MOD13A1数据，共72景影像。

1.2 研究方法

无人机样地植被覆盖度采用项目团队自主开发的覆盖度提取软件进行计算。首先利用MODIS Tool软件对下载的源区MOD13A1数据进行预处理，包括定义投影、裁切等，并提取、计算得到归一化植被指数（NDVI）、增强型植被指数（EVI）、土壤调节植被指数（SAVI）和修正的土壤调节植被指数（MSAVI）[3]； 然后利用回归模型法和像元二分模型，分别基于4种植被指数对黄河源区植被覆盖度进行反演，并比较各植被指数采用上述两种方法的反演精度；最后比较采用同种反演方法，哪种植被指数的反演精度更高就采用该植被指数进行逐年反演，并分析黄河源区 2001－2018年植被覆盖度的变化情况。

2 黄河源区植被覆盖度反演

2.1 黄河源区样地植被盖度提取

本文采用研究团队自主开发的Pixel Based Manual Classifier V1软件提取航拍照片的植被覆盖度和非植被等信息。该软件能利用阈值法快速、高效、准确地区分植被和非植被，并根据分类结果进行植被覆盖度的计算。

2.2 植被覆盖度反演

2.2.1 植被指数数据获取

首先从经过预处理的MOD13A1数据中提取NDVI、EVI，再提取红光波段和近红外波段，用以获取SAVI和MSAVI[4]。

2.2.2 基于植被指数的回归模型法

基于植被指数的回归模型法首先建立植被指数与实测植被覆盖度的关系，再结合卫星遥感数据将该关系反推到整个黄河源区，从而获得黄河源区的植被覆盖度数据。

本文首先分别计算2001－2018年各种植被指数的均值，再采用“提点法”得到各植被指数的像元值，并利用其中75%的样本点数据与植被覆盖度建立拟合关系，然后通过拟合关系进行反演，最后利用25%的样本点进行验证，比较其反演精度。

1）利用75%的植被覆盖度对各植被指数均值进行采样，并建立植被指数与实测植被覆盖度之间的散点图，得到拟合关系式，如图2所示。

图2 植被指数与实测FVC建立的拟合关系

2）在ArcGIS中，分别利用各植被指数的拟合关系式反演整个黄河源区的植被覆盖度，结果如图3所示。

图3 回归模型法各植被指数的反演结果

3）利用剩余的25%的样本点数据进行精度验证，结果如图4所示。

图4 回归模型法精度验证

2.2.3 像元二分模型法

像元二分模型法是将像元中的地物组成分为植被与非植被两部分，通过获取纯植被像元与非纯植被像元的光谱信息来达到分解混合像元的目的。通过分解可计算得到该混合像元中植被因子所占的比例。以NDVI为例，植被覆盖度的计算公式为：

式中，NDVIsoil为纯土壤像元对应的NDVI值；NDVIveg为纯植被像元对应的NDVI值；NDVI为需要计算植被覆盖度的混合像元的NDVI值。

具体取值方法为：若有实测数据，则分别取实测植被覆盖度的最大值和最小值作为FVCmax和FVCmin，遥感影像中的NDVI值作为NDVImax和NDVImin，与实测的FVCmax和FVCmin数据相对应[6]。

与回归模型法一样，采用75%的样本点进行计算，25%的样本点进行精度验证，拟合结果如图5所示。结合样地所在像元拟合得到的植被覆盖度与样地实测植被覆盖度进行精度验证，结果如图6所示。

图5 像元二分模型法各植被指数反演结果

图6 像元二分模型法精度验证

2.3 逐年反演植被覆盖度

两种方法精度验证的相关系数如表1所示，可以看出，对于4种植被指数来说，回归模型法的模拟精度均优于像元二分模型法；对于回归模型法来说，4种植被指数中MSAVI的反演精度最高，R2达到了0.815 6。 因此，本文采用回归模型结合MSAVI建立的关系：y=146.4x+0.853 3，作为反演2001－2018年黄河源区植被覆盖度的依据。

表1 两种方法精度验证的相关系数R2

3 研究结果分析

3.1 2001－2018年植被覆盖度变化分析

本文采用最小二乘线性回归模型计算黄河源区2001－2018年植被覆盖度的变化情况。该方法常用于分析长时间序列的植被变化趋势，将时间变量t作为独立变量，与各年的NDVI时序数据进行最小二乘回归分析，从而得到一个线性方程；再利用该方程，根据其斜率数值，获得源区2001－2018年植被覆盖度的变化趋势。若斜率小于零，则植被生长呈下降趋势；反之，则呈增长趋势。时间变化公式为：

式中，y为各年的平均植被覆盖度；a为拟合参数；b为线性回归方程的斜率[6]。

2001－2018年黄河源区植被覆盖度变化趋势如图7所示，可以看出，2001－2018年源区的植被覆盖度整体呈增长趋势，但增长不是很明显，且自2001年 以来，黄河源区植被覆盖度在50%～70%之间。

图7 2001－2018年黄河源区植被覆盖度变化趋势分析

3.2 黄河源区多年植被覆盖度均值统计分析

本文根据植被覆盖度划分标准（表2）对黄河源区植被覆盖度多年均值进行分级[13-14]，结果如图8所示，并统计了不同覆盖度等级所占的比例（表3）。

表2 植被覆盖度划分标准

结合图8和表3可知，整个黄河源区内植被覆盖度＜10%的裸地，占比为1.41%；植被覆盖度为10%～30%的低覆盖区，占比为11.69%；植被覆盖度为30%～45%的中低覆盖区，占比为17.66%；植被覆盖度为45%～60%的中覆盖区，占比为19.59%；植被覆盖度＞60%的高覆盖区，占比为49.65%；因此黄河源区约50%的地区其植被覆盖度在60%以上。

图8 黄河源区植被覆盖度分级结果

表3 黄河源区植被覆盖度不同等级的面积统计

4 结 语

1）2001－2018年黄河源区的植被覆盖度整体呈增长趋势，但该趋势不是很明显，且自2001年以来，黄河源区植被覆盖度在50%～70%之间。

2）整个黄河源区内植被覆盖度＜10%的裸地，占比为1.41%；植被覆盖度为10%～30%的低覆盖区，占比为11.69%；植被覆盖度为30%～45%的中低覆盖区，占比为17.66%；植被覆盖度为45%～60%的中覆盖区，占比为19.59%；植被覆盖度＞60%的高覆盖区，占比为49.65%；因此，黄河源区约50%的地区其植被覆盖度在60%以上。

本文只选用了MODIS数据进行研究，没有结合其他较长时间序列数据（如GIMMS）进行长时间序列的黄河源区植被覆盖度的变化趋势分析；计算得到了黄河源区2001－2018年的植被覆盖度变化趋势，但未进行影响因素分析，这也是今后需继续研究的内容。