基于5种植被指数的荒漠区植被生物量提取研究

梁博明，刘 新，郝媛媛，楚 彬，唐庄生

（甘肃农业大学草业学院，草业生态系统教育部重点实验室，国家林业草原高寒草地鼠害防控工程技术研究中心，甘肃 兰州 730070）

土地荒漠化严重威胁着人类的生存环境，制约着全球社会和经济的稳定发展，受到国内外学者的高度重视。中国是世界上荒漠分布广、面积大、危害严重的国家之一［1］。第五次全国荒漠化和沙化土地监测结果显示，中国荒漠化土地面积2.61×106km2，占国土面积的27.20%［2］，近4×108人口受到荒漠化的影响［3］。中、美、加国际合作项目研究表明，我国因荒漠化造成的直接经济损失约为每年5.41×1010元人民币［4］。

在荒漠区，地上生物量是反映植被生长发育状况的重要指标，对植被准确估产和品质测定具有重要作用。因此，准确获取植被地上生物量信息对植被长势监测、产量估算、植被区维护与治理等具有重要意义。传统的植被地上生物量测量多是破坏性取样，费时费力，难以满足生态环境中对大面积地上生物量准确及时估测的需求。多光谱遥感监测具有快速、准确、无破坏性等特点，已广泛应用于植被地上生物量的快速获取［5］。植被指数是不同遥感光谱波段间的线性或非线性组合，主要反映植被在可见光、近红外波段反射与土壤背景之间差异的指标，各个植被指数在一定条件下能用来定量说明植被的生长状况［6］，而建立植被指数与植被地上生物量的估算模型已广泛应用于现实研究中。叶静芸等［7］利用样方调查数据和Quick Bird 影像数据对乌兰布和沙漠东北缘荒漠-绿洲过渡带植被地上生物量的研究表明，比值植被指数（Ratio Vegetation Index，RVI）线性模型估算荒漠植被地上生物量的效果最好（R2=0.82，RMSEP=15.07）；陈琪等［8］为了分析内蒙古阿拉善盟草地地上生物量的时空分布规律，构建了1981—2000 年NOAA（National Oceanic and Atmospheric Administration，美国国家海洋和大气管理局）NDVI（Normalized Difference Vegetation Index，归一化植被指数）、2001—2018年MODIS（Moderateresolution Imaging Spectroradiometer，中分辨率成像光谱仪）NDVI 遥感数据产品和实测地上生物量的关系模型，发现幂函数的相关性最优，模型精度可达78%；Wu等［9］基于天宫一号高光谱数据的荒漠化地区稀疏植被生物量估测中，SAVI（Soil-Adjusted Vegetation Index，土壤调节植被指数）的一元线性模型R2高于NDVI；杨丽萍等［10］在内蒙古额济纳的研究表明，基于OSAVI（Optimization Soil-adjusted Vegetation Index，优化型土壤调节植被指数）的TVDI（Temperature Vegetation Dryness Index，温度植被干旱指数）是研究土壤水分反演的最佳模型，可以间接反应研究区地上生物量的分布状况；张锦丽等［11］研究表明，利用Lanndsat-8 OLI 影像提取的植被指数DVI（Difference Vegetation Index，差值植被指数）与研究区地上生物量实测数据具有相关性。由此可知，数据源、研究区、植被指数以及拟合关系等都是决定植被地上生物量信息提取模型的因素。

哨兵（Sentinel）数据作为目前分辨率最高（多光谱10 m）的全球性开源数据［12］，在植被生物量提取方面相较于MODIS、Landsat等开源数据无疑具有更高的精度，但已有研究中利用其提取荒漠区植被生物量的研究相对较少。位于巴丹吉林沙漠和腾格里沙漠之间的民勤县是典型的荒漠-绿洲生态系统，是荒漠植被生物量研究的理想场所。因此，本研究以民勤县为研究区，基于Sentinel-2 数据构建并筛选不同植被指数与实测地上生物量之间的最优估算模型［13］，并评估民勤县的植被生物量分布状况，以期为荒漠区植被地上生物量快速提取研究提供参考。

1 数据与方法

1.1 研究区概况

以甘肃省民勤县为研究区，地理坐标（101°47′39″～104°57′3″E，38°0′39″ ～39°27′49″N，图1），总面积1.58×106hm2，海拔1298～1936 m，地处河西走廊东北部，石羊河流域下游［14］。东西北三面被腾格里沙漠和巴丹吉林沙漠包围，多年（1953—2021年）平均气温8.5 ℃，平均降水量114 mm，平均蒸发量2412 mm，属于典型的温带大陆性干旱气候。区域植被特征明显、群落片层简单，以“灌草丛”为主［15］。由荒漠（沙漠）、荒漠-绿洲过渡带和绿洲3 种西北干旱区典型的地表形态组成［16］，分布着绿洲灌区、荒漠平原区、沙漠区和山地4 种景观，是典型的荒漠-绿洲生态系统［17-18］。

图1 研究区地理位置及采样点分布Fig.1 Location and distribution of sampling points in the study area

1.2 数据来源及预处理

1.2.1 野外实测数据 野外实测点主要位于荒漠-绿洲过渡带［19］，为生态缓冲区［20］，兼顾两地植被分布差异，尽可能涵盖全县植被类型［21-22］。采样时间为植被茂盛期的2020年7月18—31日；同时为了与陆地卫星遥感数据空间分辨率保持一致［23］，样方大小设置为10 m×10 m，共64个采样点（图1）。在样方内按照五点采样法设置1 m×1 m小样方，记录样方内植被信息（涉及18科、33属、33种，分类地位见表1）后齐地面剪取5个小样方内植物并编号分类［14］，带回实验室烘干后取干重平均值，乘100 得到百平方米面积上的地上生物量，因此试验中地上生物量统一单位为kg·(100m2)-1。

表1 样方内植物物种及其分类地位Tab.1 Species and taxonomic status in quadrat

1.2.2 遥感数据 选取2020年7—8月（以匹配地面实测数据）经过预处理后的哨兵2号（Sentinel-2，表2）L2A级数据产品［23］（https:/scihub.copernicus.eu/）的6景影像（表3），并进行影像镶嵌、裁剪等预处理。

表2 哨兵2号数据波段特征Tab.2 Band characteristics of Sentinel-2 data

表3 研究区遥感影像数据列表Tab.3 List of remote sensing image data in the study area

1.3 研究方法

1.3.1 植被指数选取及相关性分析 选取估算荒漠植被生物量较好的比值植被指数RVI［7］、应用最为广泛［24］的归一化植被指数NDVI、对环境变化极为敏感且算法简单［25］的DVI、加入土壤调节系数L后可以减弱土壤背景影响［26］的土壤调节植被指数SAVI 以及适合于监测耕地人工林［27］的优化型土壤调节植被指数OSAVI 5 种植被指数，计算公式见表4，与实测地上生物量之间进行Pearson 相关性分析后用于植被地上生物量的反演。

表4 植被指数及其计算公式Tab.4 Vegetation index and its calculation formula

1.3.2 最优模型筛选及精度验证 为防止过度拟合［28］，设置64 个样点数据中训练集与验证集比例为2:1［29-30］。

基于经验统计法，构建5种植被指数与实测地上生物量（训练集：43个样点）之间的一元线性、对数、指数和二项式估算模型，并根据决定系数［31］R2（Coefficient of Determination）评估并筛选最优估算模型。

再依据估算值与实测值（验证集：21个样点）之间线性拟合的R2和均方根误差RMSE（Root Mean Square Error）对最优模型进行精度评价。计算公式如下：

式中：n、i分别为数据个数第几个数；Xˉ表示地上生物量实测的均值；Mi和Oi分别为地上生物量的估算值和观测值。

1.3.3 研究区植被生物量反演 研究区群落整体植株矮小、植被稀疏，生物量低［32-33］，利用植被类型实地占比分布状况［34］，结合植被间的共生关系［35］，按生物量大小划分民勤县植被分布状况［36-38］（表5）。

2 结果与分析

2.1 植被指数与地上生物量的相关性

RVI、NDVI、DVI、SAVI 和OSAVI 5 种植被指数与植被地上生物量均呈显著正相关关系（P＜0.05），相关系数为0.67～0.79（表6），均可用于植被地上生物量估算。

表6 5种植被指数与地上生物量的相关性Tab.6 Correlation between five vegetation indexes and aboveground biomass

2.2 最优模型筛选

基于5种植被指数分别构建地上生物量的一元线性、对数，指数和二项式估算模型（表7）。基于SAVI 的地上生物量模型效果表现最佳，4 种拟合方法的决定系数R2分别为0.62、0.21、0.73和0.76；NDVI 次之，R2分别为0.46、0.15、0.60 和0.63；RVI 4 种拟合方法的平均R2较大（0.53），仅次于SAVI（0.58），且对数估算模型是5种植被指数对数估算模型中最好的（R2=0.44），但指数和二项式2 种拟合方式均不及其他4 种指数；DVI 和OSAVI 最差。因此，基于SAVI 的二项式模型为估算民勤荒漠区植被地上生物量的最优模型。

表7 基于遥感影像植被指数的地上生物量估算模型Tab.7 Estimation model of aboveground biomass based on vegetation index from remote sensing image

2.3 估算精度评价

由地上生物量估算值与验证集之间的线性拟合关系可知（图2），决定系数R2为0.73，通过P＜0.05显著性水平检验，均方根误差RMSE 为0.12，模型精度符合要求，可用于民勤荒漠植被地上生物量信息的提取。

2.4 植被生物量反演

基于SAVI 二项式模型的民勤植被地上生物量空间分布（图3）结果表明，民勤县近2/3 区域为沙漠、戈壁等无植被区（66%）；低植被区面积次之，占研究区面积的1/5，主要分布在民红公路沿线以及绿洲外围的沙漠和荒漠中；中植被区面积最小，仅为县域面积的1/20，多分布在人工植被和高植被区域中天然植被的外围；高植被区主要分布在红崖山、环河、昌宁和南湖4大灌区（人工植被）以及青土湖、红沙岗镇西北部等区域（天然植被），占比不足10%。

图3 基于SAVI二项式模型的研究区地上生物量空间分布Fig.3 Spatial distribution of aboveground biomass in the study area based on SAVI index binomial model

3 讨论

民勤地处腾格里沙漠与巴丹吉林沙漠的交接处，其土地干旱、植被稀疏［39］。植被地上生物量为体现植被生长状况最重要的指标，其值的大小直接反应了该地区植被生长的优劣［40］。为了能更好的反演地上生物量的分布状况，研究采用了空间分辨率高、云量小的哨兵2 号遥感影像与植被地上生物量相结合的方式，估算荒漠区植被地上生物量。

通过对植被稀疏的荒漠区进行地上生物量提取分析，结果表明SAVI指数与地上生物量的相关性最好，模拟效果均优于其他植被指数，与马中刚［41］对河北康保县和张锦丽等［42］对准噶尔盆地的荒漠区植被地上生物量反演研究结果一致。这是由于SAVI指数通过向NDVI指数的分母中引入土壤调节系数L，进一步降低了土壤背景变化对植被指数的影响。L是SAVI指数提取植被信息的关键，取值范围为0～1，L=0 时，表示植被覆盖为0，L=1 时，表示植被覆盖度非常高。这与Huete等［43］提出的对于调整不同的L，SAVI 指数几乎可以排除由土壤背景等因素所造成的植被指数变化的研究结果基本保持一致。

利用经验统计法，构建基于不同植被指数的地上生物量估算模型（一元线性、对数，指数和二项式模型），并利用最优模型进一步对地上生物量的空间分布展开研究。该模型方法参数输入单一，计算简单，经常应用于遥感图像数据与野外实测数据相结合的模型构建，因此，在本次研究中，取得的结果较为良好，同时也为荒漠地带植被的地上生物量研究提供依据。

4 结论

通过对民勤地区植被地上生物量与植被指数的研究，得出如下结论：

（1）RVI、NDVI、DVI、SAVI 和OSAVI 5 种植被指数与植被地上生物量均呈显著正相关关系（P＜0.05），且SAVI的相关性最高（r=0.79）。

（2）基于SAVI 指数的二项式模型是研究区地上生物量估算的最优模型（R2=0.76），且精度较高（R2=0.73，RMSE=0.12）。

（3）民勤县植被相对密集区主要分布于四大灌区（红崖山、环河、昌宁、南湖）、青土湖周边以及红沙岗镇西北区域，其他地域植被较为稀疏，无植被区［＜0.005 kg·(100m2)-1］、低植被区［0.005～0.2 kg·(100m2)-1］、中植被区［0.2～0.5 kg·(100m2)-1］和高植被区［＞0.5 kg·(100m2)-1］的占比分别为66%、21%、5%和8%。