郭超凡,陈雯璟,牛明艳,张志高
(1.衢州学院商学院,浙江 衢州 324000;2.安阳师范学院资源环境与旅游学院,河南 安阳 455000)
地上生物量(Aboveground Biomass,AGB)是草地生态系统的物质基础,不仅能够客观反映草地的生长状况和草场载畜量,还是草地碳库的重要组成部分[1]。及时、准确地掌握草地地上生物量的累计状态、分布及动态变化对于评估草地生态系统、研究全球碳循环和确保草地生态安全等均具有重要意义。
目前,遥感技术被认为是最有效、最具潜力的草地生物量估算方法。遥感估算技术是以植被指数为主要输入变量,基于像元数理统计的回归分析方法,通过不同尺度数据之间建立函数关系来完成由点及面的转换[2]。植被指数的本质是多波段反射率的数学变换,使变换后的数据在增强植被信息的同时最小化非植被特征[3],因此植被指数可以较好地反映植物的生长状况及空间分布。Tucker等[4]发现通过对近红外波段和红光波段反射率进行归一化处理获取的NDVI(Normalized difference vegetation index)对于生物量的变化十分敏感。此后,为了突出植被特征、消除噪声因子影响,国内外学者相继提出了一系列的衍生植被指数(至今已有40余种),如反映植被衰老变化的NDVIre(Red-edge NDVI)[5],消除大气干扰的ARVI (Atmospherically resistant vegetation index )[6]。虽然大多数衍生指数与生物量存在较好的拟合关系,但研究发现不同指数对于生物量的敏感性存在最优适用范围,且一些指数模型的适用范围存在互补性。如郑阳等[7]的研究结果表明,NDVI、NDVIre等归一化差值指数在生物量低于1 000 g·m-2保持较高的敏感性,而SR (Simple ratio)、CIre(Red-edge chlorophyll index)等比值指数在生物量高于1 000 g·m-2时敏感性增强。近年来,植被指数适用性研究逐渐受到相关学者关注[7-9],但基于适用范围互补性的多个指数模型的综合运用研究尚不多见。因此,评估不同植被指数在草地生物量估算中的适用性并根据其适用范围构建多指数模型的协同估算方案,对于提高草地资源的定量监测精度具有重要的意义。
Sentinel-2相较于Landsat系列数据,具有更高的空间分辨率(10、20 m)和更多的光谱波段(13个波段),其独有的红边波段对于植被多种理化参量的变化敏感,可以衍生出更多具有不同生态学意义的植被指数[10-12],因此具有更大的草地生物量估算潜力。本研究以青海省金银滩草场为研究区,以Sentinel-2数据提取的18种植被指数为数据源,建立植被指数与草地生物量的拟合模型,并对其进行敏感性分析;对比不同指数模型的适用范围,探索通过构建多指数模型协同反演方案来提高草地生物量空间分布制图精度的可能性。研究结果以期为草地生物量遥感监测提供理论依据,为青海草原“智慧畜牧业”发展提供借鉴。
金银滩草原位于青海省海北藏族自治州海晏县境内,地理位置36°53′30″~37°5′30″N, 100°47′30″~100°59′10″E(图1);年日照时数2 980 h,年平均温度1.7℃,年均降水量499 mm,夏秋降水多,春冬降水少,属于高原内陆型气候,全县草地面积占总面积49.35%,草场资源丰富,草种类型多样,是全国草地生态畜牧业试验区。草地类型以高寒草甸类、高寒草甸草原类和温性草原三大类型草地为主。
地上生物量(AGB)数据采集于2017年8月5—6日进行。实地调查发现,研究区内草场实行轮牧制度,即草场只选择在非生长季节放牧(11月—次年4月,称为冬季牧场)或生长季节放牧(5—10月,称为夏季牧场)。根据研究区内牧草生长状况和轮牧类型,设置了3个采样区,包括1个夏秋草场(采样区Ⅰ)和2个冬季草场(采样区Ⅱ,Ⅲ)。3个采样区随机采集97个样方(图1),样方尽可能代表整个研究区域牧草的生长状况。样方设置在地势平坦、草地优势种单一、面积大于10 m×10 m的样地内,且样方应位于样地的中心。样方规格为0.5 m×0.5 m,齐地刈割,挑出石子和动物粪便等不可食部分并称重。最后中心样方的称重结果作为10 m×10 m样地的草地AGB。
图1 研究区位置及样方分布Fig.1 Location of the study areas and distribution of sampling points
研究所选用的遥感数据为Sentinel-2影像,包含13个波段,其中490、560、665 nm和842 nm波段空间分辨率为10 m,705、740、783、865、1 610 nm和2 190 nm波段分辨率为20 m,443、945 nm和1 375 nm波段分辨率为60 m。影像过境时间为2017年8月4日。Sentinel-2数据下载网址:https://scihub.copernicus.eu/dhus/#/home。使用SNAP对影像进行预处理,经过辐射定标,大气校正后得到反射率数据。Sentinel-2数据各个波段的空间分辨率有所不同,本文使用最近邻插值法,将处理后的各波段重采样至10 m。
选取常用于草地生长状况监测的18种植被指数,不同指数的计算公式见表1。由于Sentinel-2包含3个红边波段,而相关研究表明705 nm和740 nm处的反射率与叶绿素含量均具有较高的相关性[10],因此本文分别选用中心波长在705 nm和740 nm的波段作为计算中的红边波段。同时,由于中心波长在2 190 nm处的反射率与植被水分相关性优于1 375 nm和1 610 nm处的反射率,因此本文选取中心波长在2 190 nm的波段作为计算中的短波红外波段。
表1 植被指数公式Table 1 Vegetation index formula
Nguy-Robertson等[27]研究表明,在单变量模型拟合过程中,当估算方程呈非线性时,R2和RMSE的评估结果会存在偏差,需采用噪声等效 (Noise equivalent,NE) 进行模型的敏感性分析,对模型的适用性给出科学合理的分析评价。NE值越小,则说明植被指数对生物量的敏感度和适用性就越强。地上生物量估算的噪声等效NEΔAGE(Noise equivalent AGE)计算公式[28]:
式中,RMSE(VIvsAGB)表示植被指数关于生物量最优拟合函数的均方根误差,d(VI)/d(AGB)表示植被指数关于生物量的最优拟合函数对生物量的一阶导数。
对获取的样本点生物量进行统计分析,结果如表2所示。97个样方的生物量波动区间为0.040~1.744 kg·m-2,均值为0.606 kg·m-2。其中采样区Ⅰ生物量波动区间为0.040~0.484 kg·m-2,均值为0.247 kg·m-2;采样区Ⅱ生物量波动区间为0.300~1.444 kg·m-2,均值为0.763 kg·m-2;采样区Ⅲ生物量波动区间为0.364~1.744 kg·m-2,均值为0.810 kg·m-2。结果表明,受到轮牧制度、放牧强度、草地类型以及人类活动等因素影响,青海省金银滩草原牧草生物累积量存在空间差异。
表2 不同样区生物量统计结果/(kg·m-2)Table 2 Biomass statistics results in different sampling areas
为了直观地展示植被指数与草地生物量之间的拟合关系,绘制了18种指数与生物量之间的1∶1关系图,并进行了拟合模型的构建与最优模型的筛选,结果如图2所示。从图中能够看出所选植被指数与草地生物量之间均呈现显著的非线性拟合关系,R2为0.52~0.73,RMSE为0.230~0.283 kg·m-2,说明选择的指数均能够较好地反映草地生长状况,但效果存在差异。其中SR、mSR、SRre、mSRre、CIgreen、CIre、MTCI和MTVI2与生物量呈幂函数关系;而NDVI、mNDVI、NDVIre、mNDVIre、WDRVI、NDII、EVI、OSAVI、NDWI和GVMI与生物量呈指数函数关系,且NDII模型所对应的R2=0.73、RMSE=0.234 kg·m-2,在所有指数模型中精度最高。同时研究结果表明,18种指数拟合模型均受到不同程度“过饱和”问题的影响,且不同指数的饱和点存在差异。植被指数的“过饱和”问题已经成为定量遥感研究中的一个重要制约[7]。
图2 调查点植被指数和地上生物量的最优拟合散点图Fig.2 Best-fit scatter plot of vegetation index and aboveground biomass at survey site
留一交叉验证的结果表明(表3),不同植被指数最优模型评价结果和交叉验证结果基本一致。其中,典型植被指数NDVI和SR模型均具有较高的精度,而改进的mNDVI和mSR模型精度均低于NDVI和SR,这与郑阳等[7]的研究结果不一致,可能原因是枯黄期的干枯牧草影响了该类指数的灵敏性[29]。与植被叶面/冠层水分密切相关的NDWI、GVMI和NDII均具有较高的拟合精度,由于该类型指数包含了短波红外波段,研究表明短波红外对于植被水分含量变化十分敏感[23]。叶绿素作为反映植被生长状况的重要指标,对于植被生物量的累计同样十分敏感。由表3可知,能够反映植被叶绿素的CIgreen和Clre模型同样具有较好的估算结果。而其他反映植被叶面积指数(MTVI2、WDRVI)、缓解“过饱和”问题(EVI、WDRVI)、降低背景噪声(OSAVI)的指数并未在草地生物量估算中表现出优于NDVI和SR的结果。综上所述,影响草地生物量估算精度的主要因素是水分和叶绿素。
表3 不同指数最优模型精度评价及交叉验证Table 3 Accuracy evaluation and cross validation of the optimal models with the selected vegetation indices
NEΔAGE能够反映植被指数对于生物量变化的响应能力。RMSE越小,一阶微分绝对值越大,NEΔAGE越小,就证明植被指数对生物量的敏感度和适用性越强。以NDVI和SR建模精度和验证结果为标准,选择NDVI、SR、CIgreen、CIre、NDII、NDWI和GVMI 7个植被指数进行NEΔAGE分析,结果如图3所示, 7种植被指数的NEΔAGE均随着生物量的增加呈上升趋势。其中,NDVI、NDII、NDWI和GVMI所对应的NEΔAGE与生物量的变化呈经过(0,0)点的线性关系,而CIgreen、CIre和SR所对应的NEΔAGE与生物量的变化呈幂函数关系。两种不同NEΔAGE的变化趋势在0.6 kg·m-2附近相交。
图3 植被指数的敏感性分析Fig.3 Sensitivity analysis of vegetation index
综合考虑不同植被指数所对应的生物量最优拟合模型精度和NEΔAGE,幂函数关系模型中CIgreen的模型精度和验证精度高于CIre和SR,而对应的NEΔAGE值均小于CIre和SR,说明在幂函数关系模型中CIgreen指数适用性最佳;指数函数关系模型中NDII的模型精度和验证精度均高于GVMI、NDVI和NDWI,而对应的NEΔAGE小于GVMI、NDVI和NDWI,说明在指数函数关系模型中NDII适用性最佳。同时,NDII和CIgreen所对应的NEΔAGE在0.65 kg·m-2附近相交,那么在生物量低于0.65 kg·m-2时,NDII在所有选择的指数中具有最佳的适用性,而CIgreen在生物量高于0.65 kg·m-2时适用性最优。
上述研究结果表明NDII和CIgreen模型适用范围具有互补优势,因此本文尝试采用NDII和CIgreen模型进行研究区草地生物量协同反演制图。具体流程是:分别利用NDII和CIgreen最优拟合模型对研究区观测时间段的生物量空间分布进行制图。在NDII和CIgreen制图结果均大于0.65 kg·m-2的区域,选用CIgreen的制图结果;在NDII和CIgreen制图结果均小于0.65 kg·m-2的区域,选用NDII的制图结果;剩余区域为了消除模型间的误差采用NDII和CIgreen制图结果的平均值(即均值模型)。剔除城区、道路和水域等非植被区域获取的研究区生物量反演制图结果如图4所示。结果反映出研究区草地生物量分布具有明显的空间差异性。生物量最高的区域是优质草地种植基地,同时远离城区的草地生物量较高,而城区周边的草地生物量明显较低,可能是由于城区周围多为夏季牧场,牛羊放牧制约了草地生物量的累积,此外旅游开发以及人为活动也会在一定程度影响草地的生长。
图4 研究区草地生物量估算结果Fig.4 Estimation results of biomass in study area
图5 生物量预测值与实测值的关系Fig.5 Relationship between predicated biomass and measured biomass
牧草的水分和叶绿素含量是表征草地地上湿生物量的主要指标[30]。本文通过模型拟合发现,能够反映植被叶面(冠层)水分含量的NDWI、GVMI和NDII指数和能够反映植被叶绿素含量的CIgreen和Clre指数均对于草地生物量变化十分敏感。这个结果与董建军等[31]利用多源卫星遥感数据得出的结果一致,说明这类植被指数与牧草生物量不仅存在统计意义上的相关性,同时可以从植被理化特征角度解释两者的关联,具有很好的稳定性。
近年来,遥感估算研究开始关注不同植被指数的适用性[7-9],但基于适用范围的多指数模型综合运用研究则相对薄弱。本文利用NEΔAG对比分析了不同植被指数在牧草生物量估算研究中的适用范围,并在此基础上提出了基于NDII和CIgreen的多指数模型牧草生物量协同估算方案。结果表明,上述协同方案较常规的单变量模型具有更高的精度,能够应用于整个研究区生物量的反演制图。理论上,该方法也可应用于其他植被理化参量如叶绿素、水分含量的估测。本研究的方案针对其他生理生化参量反演的适用性有待今后进一步研究。
采用Sentinel-2遥感影像结合地面实测数据进行牧草生物量估算研究,综合分析、评价了18种具有不同理化特征含义植被指数在牧草生物量估算中的适用性,并在此基础上利用提出的多指数模型协同估算方法完成对研究区生物量空间分布的高精度制图,结论如下: