草原植被长势遥感监测方法适宜性研究*

饶新宇, 李红军, 张圣微, 雒 萌, 刘志强, 张静文

(1.内蒙古农业大学水利与土木建筑工程学院 呼和浩特 010018; 2.中国科学院遗传与发育生物学研究所农业资源研究中心/中国科学院农业水资源重点实验室/河北省节水农业重点实验室 石家庄 050022; 3.内蒙古自治区水资源保护与利用重点实验室/内蒙古自治区农牧业大数据研究与应用重点实验室 呼和浩特 010018; 4.衡水市园林中心 衡水 053000)

草地资源是陆地生态系统的重要组成部分, 在人类活动、自然生态平衡及陆地生态系统的能量流动和物质循环中发挥着重要作用[1-2]。我国是草地资源大国, 拥有各类天然草地约3.9×108hm2[3], 牧草种类繁多, 保障了我国畜牧业的发展[4]。但由于长期的掠夺式放牧、气候变化、鼠虫害以及大面积草地的开垦, 草地资源大规模退化, 恶化的草原生态环境威胁到畜牧业的可持续发展[5-7]。草地资源作为畜牧业发展的基础, 保护和建设草地资源也是我国西部大开发战略的重要内容, 对促进农业经济结构的战略性调整具有重要意义[8]。2021年我国中央一号文件指出, 要科学开展大规模国土绿化行动, 稳步恢复草原生态环境。草地作为可更新资源, 地域辽阔, 要实现可持续利用离不开对草原的快速监测。因此开展草原植被长势监测方法的研究和应用, 对于草地资源的管理、保护以及生态环境的改善, 具有重要的科学意义和应用价值[9]。

草原植被长势是指其植被的生长状况和趋势。传统的草原植被长势监测多采用地面人工调查的方法, 主要测量和记录牧草的植被覆盖度、高度、牧草发育期和产草量等植被长势信息。人工监测方法耗时费力, 且只能得到若干点的数据, 覆盖范围有限,代表性较差, 无法快速反映辽阔草原长势的时空变异, 影响到草原管理措施的时效性[2,10-12]。遥感具有高空间、高时间分辨率, 以及敏感波段与地表植被状况存在密切相关的特征[13-14], 为植被长势的大范围快速监测提供了支持[15-16]。遥感技术首先在农作物长势监测中得到发展和广泛应用。作物长势遥感监测的实质是通过遥感反演获得可反映其生长特征的生物学指标。因为植物的生长都需要叶片的光合作用, 因此叶面积指数(LAI)就是一个与作物长势特征密切相关的综合指数。基于近红外和红光波段对植被特征的敏感反映, 将这2个波段组合后形成的植被指数则可实现对LAI的估算[17]。常用的植被指数包括归一化植被指数(NDVI)、垂直植被指数(PVI)、比值植被指数(RVI)和增强型植被指数(EVI)等, 其中NDVI最为常用[18-20], 可直接用于对作物长势的监测。由于农作物和草原植被长势监测对象均是地表植被, 因而草原植被长势遥感监测的方法多是借鉴于农作物长势的遥感监测。李楠等[21]基于长时间序列的MODIS NDVI数据分析了若尔盖高原植被生长季NDVI的变化趋势、波动程度及未来变化趋势。田海静等[5]利用MODIS NDVI对我国北方草原草地植被长势动态变化进行动态监测, 分析发现近20年我国北方草原整体植被呈现恢复状态, 但局部地区植被存在退化现象。

农作物长势的遥感监测方法主要包括: 直接监测法、植被生长过程曲线法、同期对比法, 其中同期对比法的应用最为广泛[22]。考虑到同期对比法中不同时期土地利用变化以及参照年份作物长势不确定性对监测评价的影响, Li等[23]提出了利用作物多年同期NDVI大数据的百分位数评估法, 从而实现了作物长势的定量遥感监测评价。上述遥感监测方法都是利用了植被的NDVI, 从不同角度对其长势进行评估。这些方法在草原植被长势监测中适宜性如何? 其评估结果对草地长势的反映有何差异? 与农作物长势遥感监测相比, 草原植被遥感监测的研究相对较弱, 对于上述问题的研究, 将有助于遥感技术在草原植被长势监测中应用和发展。基于此, 本研究以内蒙古自治区西乌珠穆沁旗(简称西乌旗)草原植被监测为例, 对目前常用的长势遥感监测方法进行适宜性分析, 比较不同方法对草地植被长势监测的差异及存在的问题, 以期为不同遥感监测方法的选择或综合应用提供依据, 促进草原植被遥感监测技术的研究和应用。

1 研究区概况与数据来源

1.1 研究区概况

西乌旗位于锡林郭勒盟东部(116°31′～119°46′E,43°93′～45°37′N), 下辖2个苏木, 5个镇, 土地面积22 434.5 km2, 总人口9万。全旗可利用草场面积20 290 km2, 占土地总面积的90%, 草原类型主要包括草甸草原和典型草原(图1)。西乌旗属于温带干旱半干旱大陆性气候, 冬春寒冷漫长, 4月中旬−10月中旬日平均气温高于0 ℃。近5年(2016−2020年)年平均气温为2.82 ℃, 年平均降水量为402 mm。

图1 西乌旗草原类型及采样区分布Fig.1 Distribution of grassland types and sampling areas of the study area

1.2 数据来源与处理

1.2.1 遥感数据与数据预处理

本研究使用MODIS的植被产品MODIS13Q1数据 (http://ladsweb.modaps.eosdis.nasa.gov), 数据格式为hdf, 空间分辨率为250 m, 时间间隔为16 d, 时间跨度为2016−2020年。MODIS NDVI产品提高了空间分辨率及其对地表植被叶绿素的相关性, 能够很好地反映植被长势的变化, 是目前植被长势评估中应用较为广泛的数据产品。利用MRT (MODIS Reprojection Tools)对MODIS13Q1数据进行NDVI数据提取、拼接、投影及格式转换等预处理。为了去除噪声的影响, 利用TIMESAT软件里的Savitzky-Golay滤波对MODIS NDVI时间序列进行平滑处理。

获得内蒙古草原NDVI植被指数序列数据后, 参考http://www.globallandcover.com官网提供的土地覆盖类型数据, 对研究区草原用地类型进行提取。利用决策树法通过NDVI阈值调整, 获得研究区典型草原与草甸草原的分布, 对分类结果进行去斑处理, 结合地面调查数据对草原类型分类精度进行检验, 总体精度达81% (图1)。

1.2.2 地面实测数据

研究区草地生物量地面监测调查时间为2017年7月31日至8月1日。在西乌旗选择长势均一、具有代表性的24个样区(图1), 每个样区设置3个1 m×1 m的样方, 记录每个样区的地理位置, 采集全部地上生物量并称重, 在实验室内105 ℃恒温下烘干至恒重, 称取地上生物量干重, 用于指示其草原长势。

2 研究方法

2.1 直接监测法

直接监测法利用遥感反演的植被指数通过不同等级阈值的划分来直接区分草原长势的差异。本研究利用MODIS NDVI数据, 通过对NDVI的等级划分来实现对草原植被长势的差异判读。

NDVI计算公式:

式中: NIR为近红外波段反射率,R为红光波段反射率。

2.2 植被生长过程曲线法

草原植被生长是一个随着时间渐变的过程, 按时间序列统计监测区域内草原NDVI的平均值, 即构建监测地区草原植被NDVI随时间的变化过程曲线,通过当年草原生长过程与参考年份的生长曲线的比较, 来评价当年长势好于、持平于或差于参考年份。如果生长曲线高于参考年份, 则表示其长势要好于同期, 反之亦然。

2.3 同期对比法

利用当年草原植被NDVI与去年或某一参考年份同期的草原植被NDVI相减, 得到相比参考年份草原植被长势的变化情况, 按照差值结果设定相应的阈值, 将草原植被长势分为: 好于参考年, 与参考年持平, 差于参考年。计算公式为:

式中: Δ NDVI 为当年与参考年草原植被NDVI的差值, NDVI1为 当年NDVI值, NDVI2为参考年份草原的NDVI值。

2.4 基于NDVI的百分位数法

基于NDVI的百分位数法利用了多年同期同类地物NDVI的大数据, 分析获得任一NDVI数值在这一大数据中对应的百分位数, 进而实现对其长势的定量评价。该方法仍属于同期对比, 对研究区内本年度及过去4年的同期草原NDVI大数据进行统计,计算不同NDVI值对应的像素数量, 按照NDVI值的大小进行排序, 统计不同NDVI值对应的累计百分位, 建立近5年同期的NDVI百分位查询表。计算公式如下:

草地NDVI取值范围为0到1, 为了便于计算,将NDVI值乘以100后取整(NDVIint), 统计不同NDVI值对应的像素数量, 再进行百分位数的计算。

式中:Pn为NDVIint=n时的百分位数(%),n的取值范围为0到100, NUMndvi,j为NDVIint=j时的像素数量,NUMall为所有像素的数量。

建立NDVI百分位查询表后, 查询任一草原像素NDVIint对应的百分位数, 以其百分位数代表其当前的长势评分, 生成对应的百分位数图像, 从而评价当前草原植被长势在这几年中的长势水平, 如果百分位数较高, 则反映其长势在近5年中处于较好的水平, 反之亦然。本研究考虑到草甸草原与典型草原整体长势的自然差异, 对两种草原类型分别建立NDVI百分位查询表进行长势评价。

3 结果与分析

3.1 直接监测法结果

利用西乌旗2020年7月26日草原的NDVI数据, 按照0.1的间隔进行长势等级划分, 实现对其长势差异的直接监测, 其分布如图2所示。该时期西乌旗草原植被处于盛草期, 草地植被覆盖度普遍较高, 植被指数NDVI趋于全年最大值。在草地长势分布上, 西乌旗东南部大部分地区草原植被的NDVI>0.6, 西北部大部分地区NDVI<0.6, 长势存在明显的区域差异。西乌旗草甸草原集中在其东南部, 而西北部多为典型草原, 草地长势与草原类型的空间分布基本一致。直接监测法能够区分草地长势的空间差异, 但由于NDVI等级划分阈值的限制, 对同一类型草地内部的长势差异信息反映不够清晰, 即使典型草原长势较好, 草甸草原长势偏差, 其二者在自然条件上的差异仍会造成其整体长势的地区差别。图2中草甸草原NDVI>0.7的比例为60%, 典型草原NDVI<0.6的比例为73%。

图2 基于NDVI数值的西乌旗草原长势等级划分(2020年7月26日)Fig.2 Division of grassland growth based on NDVI for the study area (July 26, 2020)

3.2 植被生长过程曲线法监测结果

从11月中旬到翌年1月底, 西乌旗大部分地区被雪覆盖, 长势监测没有意义。本研究选取雪地消融后的时段, 计算各个时期该区域全部草原的NDVI平均值, 按照时间序列构建其生长过程曲线。按照同样的时间步长构建2016−2020年平均NDVI的变化曲线, 并对二者进行比较, 如图3所示。

图3 2020年草原植被生长过程曲线与近5年(2016-2020年)平均生长过程曲线比较Fig.3 Comparison of grassland vegetation growth process in 2020 and the average one for 2016−2020

2020年西乌旗草原生长曲线从2月初至4月中旬低于同期近5年的平均NDVI值, 2020年4月中旬之后的NDVI曲线明显高于平均曲线, 说明本年度草原长势在前期偏差, 进入5月份后逐渐好转并明显好于近5年平均长势。2020年与2016−2020年平均NDVI曲线最大值都出现在7月底左右, NDVI最大值分别为0.60和0.57。8月中旬之后, 草地逐渐干枯, NDVI值开始减小, 2020年与近5年NDVI曲线接近重合。就全年生长过程曲线监测点的NDVI累计值而言, 2020年为6.27, 近5年平均累计值为6.12, 2020年草原整体长势略好于近5年平均。

3.3 同期对比法监测结果

对2020年7月26日和2019年7月26日同期草地的NDVI做差值比较, 将NDVI差值>−0.05且<0.05设定其长势为“与去年持平”, 相应地, NDVI差值<−0.05和>0.05分别设定为“差于去年”和“好于去年”, 分级结果见图4。相较于2019年同期, 西乌旗草原大部分地区好于去年或与去年持平, 其中好于去年的比例为40.5%, 多集中于西乌旗西部和北部地区的典型草原。与去年持平的比例为50.9%, 集中在西乌旗东南部地区的草甸草原。草地植被长势差于去年的比例为8.6%, 零星分布于西乌旗草原中部地区。整体来看, 2020年7月26日盛草期的草地植被长势明显好于2019年的同期。

图4 2020年与2019年7月同期草原植被长势比较Fig.4 Comparison of grassland vegetation growth in July of 2020 and 2019

3.4 百分位数法的长势监测结果

利用2016−2020年7月26日同期NDVI数据,基于西乌旗草甸草原和典型草原的分布, 分别统计近5年NDVI的数据分布, 建立2类草原类型的百分位查询表, 对2020年的草甸草原和典型草原长势分别进行评价然后将其合并成图, 结果如图5所示。2020年草甸草原长势评价得分较高的像素较多, 高于50分的草原比例为63%, 得分高于80的地区多集中在西乌旗东南部的边界地区, 得分低于20的像素很少, 零散分布于中西部。典型草原评价得分高于50的比例为66%, 得分高于80的像素不多, 分布较为分散。整体来看, 东部地区草地长势得分值要高于西部地区, 2020年西乌旗草原植被长势在近5年中处于中等偏上水平。

图5 基于NDVI百分位数法的草原植被长势评分Fig.5 Evaluation of grassland vegetation growth based on NDVI percentile method

3.5 4种草原植被长势遥感监测方法适宜性分析

由于遥感信息的瞬时性, 直接监测法只能通过草原植被NDVI的分级, 直观地反映某一时刻其长势的空间分布。要通过草原NDVI数值实现准确地进行长势分级, 需要根据地面调查确定各长势等级的NDVI阈值。本研究中简单地对草地NDVI进行了等级划分, 用于展示其长势的空间异质性, 属于半定量的监测方法。西乌旗草原NDVI分布显示东南部地区草地植被长势好于西北部地区, 这是因为西乌旗东南部地区大部分是草甸草原, 而典型草原集中分布于西北部地区。由于自然环境自身的差异, 草甸草原更适合牧草的生长, 类似于农业生产中的高产区, 植被生长明显好于典型草原, 仅单从NDVI值大小而不考虑草原类型来评估植被长势容易造成草原长势评价的偏差, 即使典型草原长势再好, 与草甸草原相比总是相对较差, 这类似于农作物遥感长势监测中高产农田与低产农田的比较, 低产区长势的变化很难被监测出来。

植被生长过程曲线法监测的是整个区域草原平均长势与参照年份的差异, 是对整体长势的相对评价。该方法通过整体长势(NDVI平均值)随时间的变化曲线, 显示不同时期该区域在整体长势上好于、持平于、差于某一年份的同期, 适用于监测区域植被长势空间相对均一的情形。本研究中, 西乌旗草原类型包括草甸草原和典型草原, 草原类型自身的差异往往会形成不同的生长曲线, 而集总式的监测结果无法提供不同草原类型的差异。图6显示了2020年草甸草原(a)与典型草原(b)不同的生长曲线及其与近5年各自平均生长曲线的比较。2种类型草原生长曲线在走势上虽然相似, 但在变化幅度上存在显著差异。草甸草原NDVI的峰值可达0.73, 而典型草原的峰值为0.55。

图6 2020年草甸草原(a)与典型草原(b)植被生长过程曲线与近5年(2016-2020年)平均生长曲线比较Fig.6 Comparison of vegetation growth process for meadow grassland (a) and typical grassland (b) in 2020 and the average one for 2016−2020

相对于植被生长过程曲线法, 同期对比法通过与不同年份草原同期NDVI的差值来半定量地评价与参考年份长势的差异。该方法的评价结果容易受参考年份的影响, 如果参考年份草原受灾或者风调雨顺而导致其长势较差或较好, 将会造成本年度长势评价过于乐观或悲观。同时, 如果土地利用类型发生变化, 如退耕还草等, 则会因为参照地物类型的不一致造成不同用地类型长势的比较。在本研究中,如果将参考年份设定为2018年或2017年, 其评价结果则如图7所示。图7a中2020年好于2018年的比例为36.97%, 与2018年持平的比例为50.66%, 差与2018年的比例为12.37%。图7b中2020年好于2017年的比例为53.27%, 与2017年持平的比例为39.29%, 差与2017年的比例为7.44%。可见对比法的评价结果易受参考年份的影响。

图7 2020年与2018年(a)和2017年(b) 7月同期草原植被长势比较Fig.7 Comparison of grassland vegetation growth in July of 2020 with 2018 (a) and 2017 (b)

百分位数法的长势监测方法利用相同草原类型多年的NDVI大数据, 建立相应的百分位数查询表,分别对典型草原和草甸草原进行评价, 解决了草地类型变化、参考年份草地长势偏差过大等原因造成的长势评价误差。百分位数法将监测区域的所有草原像素NDVI值纳入大数据中进行评比分析, 通过其所处的百分位来定量评价其长势。因为采用多年同期NDVI数据, 该方法纳入了同类草原类型可能存在的从最差(最小NDVI)到最优(最大NDVI)长势信息, 提高了长势监测评分的合理性。

3.6 基于草原采样的长势监测方法验证

因为植被生长过程曲线法与同期对比法提供的是基于比较的半定量评价结果, 不适合结合地面监测进行精度的验证。对于直接监测法(草原NDVI数值)和基于NDVI百分位数法, 本研究利用2017年7月31日至8月1日在西乌旗24个草原样区的地面监测对其长势评价结果进行精度验证。图8为采样区单位面积产草量干重与其NDVI数值、百分位数法长势评分的相关性分析。草原样区产草量干重与其NDVI数值相关性的决定系数为0.5502, 略高于其与百分位数法长势评分的相关性(0.5047), 二者均达到显著相关(P<0.01), 可见这2种方法适合于对草原长势的定量监测。

图8 2017年草原植被NDVI及百分位数长势评分与地面植被生物量相关性Fig.8 Correlation between grassland NDVI, vegetation growth score by percentile method and ground vegetation biomass in 2017

4 讨论与结论

农作物长势遥感监测关注的是其营养生长和产量形成过程, 到作物生长后期长势信息(NDVI)则不再是重点, 而应利用长势的过程信息评价其对最终产量形成的影响。草原长势监测更注重全年中生物量的增减变化, 因此植被的NDVI数值可一直作为长势监测的数据源。基于NDVI数值提出的直接监测法、植被生长过程曲线法、同期对比法都属于半定量的监测方法[24]。直接监测法利用的植被NDVI数值, 也可定量地评估植被的长势, 如本研究中植被NDVI数值与其产草量干重的显著相关性, 但要实现对草原长势的等级分类, 则需要实时的地面监测辅助信息来确定不同长势等级的NDVI阈值。随着基于物联网的草原生态监测系统的建设, 以及无人机搭载多光谱仪快速监测技术的发展, 天-空-地一体的草情长势监测平台可为这些半定量监测方法的准确分类提供实时的辅助信息[25]。

草原NDVI数值常被用于长势的评价, 而NDVI在高植被覆盖区(如草甸草原)存在饱和现象[26], 在低植被覆盖区又易受到土壤背景的影响(如荒漠草原), 虽然针对不同的草原状况可以选择适宜的植被指数以提高监测精度, 如增强植被指数(EVI)、土壤调节植被指数(SAVI)等, 但不同植被指数间的可比性及长势监测的持续性需要研究。本研究中草原样区的选点基本属于典型草原(图1), 因而其百分位数法长势评分与其NDVI数值密切相关(R2=0.55)。如果样区分别来自典型草原和草甸草原, 因为不同草原类型百分位数评分的查询曲线不同, 即不同草原类型下相同NDVI数值对应的评分标准不同, 则会降低其评分与采样区单位面积产草量干重的相关性。因此, 如果利用百分位数评分监测草原的产草量, 必须基于相同的草原类型。

及时掌握草原长势信息, 对牧民生产、政府禁放牧管理、草地资源的保护都具有重要的应用价值。遥感技术高时效、低成本、大覆盖的优势使得其特别适用于辽阔草地的长势监测。本研究以内蒙古西乌旗为研究区域, 对常用的草地长势遥感监测方法:直接监测法、植被生长过程曲线法、同期对比法、基于NDVI百分位法进行了应用和适宜性分析, 以期明确各类监测方法的技术优势、监测信息的差异及其局限性。研究结果表明: 直接监测法技术简单, 能够直观反映草地植被长势的空间差异, 且与草原单位面积的产草量显著性关, 但要进行长势等级的划分则需要地面监测数据支持, 其监测结果对不同草原类型内部长势的变化信息辨识度不高。生长过程曲线法提供了草原整体长势随时间的变化, 因其监测信息过于集总, 适用于长势较为均一的地区, 若要反映不同草原类型长势过程的差异, 则需分类型进行监测。同期对比法反映了当前与参考年份同期草地长势差异的空间异质性, 作为比较结果, 容易受参考年份草原长势变化的影响。基于NDVI的百分位法利用了监测区多年草原NDVI的大数据, 包含了从最好到最差的长势信息, 因而能合理地定量评价草地的长势, 但该方法需要多年草地NDVI数据的支持。在草原长势遥感监测工作中, 需要关注这些方法的局限性, 根据监测目的选择适宜的方法组合, 从不同角度展示草原长势信息, 为草原科学管理提供辅助依据。