聂素云, 杨 彬, 夏 微, 张 远
(1. 华东师范大学 地理科学学院 地理信息科学教育部重点实验室, 上海 200241;2. 中共辽宁省委党校 现代科技教研部, 沈阳 110004)
小麦是世界三大粮食作物之一, 在我国广泛种植. 动态监测小麦的生长状况, 对于保障粮食安全具有重要意义. 小麦冠层含水量能够直接表征植株的水分胁迫状况, 是衡量小麦长势的重要参数之一,了解小麦冠层含水量状况十分必要. 植被含水量反演源自于Tanner[1]在1963年首次提出的以冠层温度指示植被水分亏缺, 然而温度受环境影响剧烈, 难以有效反演植被冠层含水量. 随着光谱仪等硬件技术的发展, 植被含水量反演更多聚集于近红外-短波红外波段.
Gao[2]于1996年利用860 nm和1240 nm的反射波谱建立归一化水分指数 (normalized difference water index, NDWI), 能够有效估算植被等效水厚度 (equivalent water tickness, EWT), 被广泛用于高光谱数据[3-5]和MODIS (moderate-resolution imaging spectroradiometer)数据[6]反演作物冠层含水量. 受限于高光谱数据难以获取和中低分辨率影像反演精度低等问题, 近些年, 研究人员基于中高分辨率的遥感影像, 采用物理模型和经验统计方法在作物冠层含水量反演方面进行了一系列研究[7]. 在物理模型方法方面, Boren等[8]和Pan等[9]利用Landsat-8和Sentinel-2数据, 基于查找表的辐射传输模型反演农田作物冠层含水量, 结果证实, 增加先验信息提高了农田作物冠层含水量的反演性能, 优于基于光谱指数建立的经验统计模型. 此外, 物理模型参数众多且难以有效获取, 模型复杂导致不易推广. 因此, 对于遥感反演应用而言, 经验模型仍不可或缺. Zhang等[10]基于Sentinel-2数据, 评估多种植被指数对于拔节期冬小麦冠层含水量的反演能力, 研究表明, NDWI的反演结果最优 (R2= 0.68,RMSE (root mean squared error)= 0.148 g·cm–2). Han等[11]对比基于Sentinel-1的水云模型和基于Sentinel-2的植被指数在小麦冠层含水量反演中的能力, 结果显示Sentinel-2光学影像的反演结果较优 (R2= 0.632, RMSE = 2.1%), 研究表明, 利用中分辨率影像构建的NDWI在反演作物含水量方面具有较大潜力. 程小娟等[12]利用Landsat TM5在北京顺义和通州进行了冬小麦冠层含水量反演, 采用了水分指数 (water index, WI) 和NDWI, 反演精度分别为R2= 0.5711、RMSE = 3.89%和R2= 0.5116、RMSE = 0.024 g·cm–2. 侯学会等[13]基于Landsat 8遥感数据, 利用垂直干旱指数和短波红外垂直失水指数的比值形式构建了一种新型植被水分指数, 反演了多个时期的冬小麦冠层含水量, 精度优于使用单一植被指数的结果,R2与RMSE分别为0.763和2.296%. 此外, Djamai等[14]使用 Sentinel-2 MSI 和Landsat-8 OLI 数据, 对用于绘制农田生物物理变量的 Sentinel 简化 2 级产品原型处理器 (simplified level 2 product prototype processor, SL2P) 进行了验证, 其中生物物理变量包括植被覆盖度、叶面积指数和冠层含水量, 结果表明, SL2P算法能够有效反演植被覆盖度和叶面积指数, 但对于冠层含水量反演结果欠佳.
综上所述, 基于光谱指数的小麦冠层含水量反演方法具有较强的普适性, 能够为小麦大范围精准灌溉提供有效的技术手段和信息支持[15-18], 但以上研究多为单个时期或短时期的反演, 未对小麦整个生长周期的冠层含水量进行探究. 动态观测小麦整个生长周期的含水量, 有利于实时监测小麦长势,及时指导农田灌溉. 此外, 在小麦生长的前/中期, 麦田像元多为混合像元, 不同时期的小麦覆盖度差异较大, 直接采用光谱指数进行反演, 易造成较大误差. 因此有必要排除土壤影响, 进一步提高冠层含水量反演精度. 基于此, 本研究以河北省无极县及周边地区的麦田为研究区, 利用Landsat-8和Sentinel-2数据时相互补的优势, 获得覆盖小麦整个生长期的遥感观测数据. 基于像元混合分解模型构建反演方程, 排除土壤影响, 对冬小麦8个物候期的冠层含水量进行反演并制图, 为区域作物水分遥感监测提供技术参考.
研究区位于河北省石家庄市无极县, 如图1所示. 地理位置为38°02′42″N ~ 38°26′26″N,114°43′58″E ~ 115°18′33″E. 无极县地处华北平原, 为冬小麦主产区, 境内地势平坦, 水源充足. 属温带季风气候, 四季分明, 光热充足. 该地冬小麦通常于10月上旬播种, 次年6月上旬收获.
图1 研究区示意图Fig. 1 Schematic diagram of the study area
选取Landsat-8和Sentinel-2多光谱遥感影像进行时相互补, 以提升遥感数据的时间分辨率. 影像获取时间与田间测量时间基本一致, 最大时间偏差不超3天, 具体信息见表1. Landsat-8 OLI和Sentinel-2 MSI数据的光谱范围分别为430 ~ 2300 nm和433 ~ 2280 nm, 分别有9个和13个波段,前者影像空间分辨率为30 m, 后者包括10 m、20 m和60 m共3个级别. 首先, 利用ENVI软件对Landsat-8原始影像进行辐射定标、大气校正及裁剪; 利用SNAP软件对Sentinel-2原始影像进行重采样和裁剪, 最终得到研究区遥感影像的地表反射率数据. 其次, 基于4月18日的Sentinel-2影像提取研究区的小麦种植信息. 结合地面实测数据进行目视判别, 将研究区地物分为七类: 小麦、建筑 (包括道路)、水体、裸地、林草地、温室大棚和未分类地物. 在影像上均匀地选取500个样块, 其中, 5%作为训练样本, 95%作为测试样本. 采用支持向量机方法进行分类, 总体分类精度为95%, kappa系数是0.913, 满足反演要求. 在分类基础上, 将小麦种植区域与其他地物进行分离, 见图2.
表1 遥感数据列表Tab. 1 List of remote sensing images
图2 冬小麦种植空间分布Fig. 2 Spatial distribution map of the winter wheat planting area
田间实测数据的采样时间为2017年1月22日—5月28日, 主要采集了149个样点的地理坐标、小麦冠层含水量实测值、土壤含水量实测值 (测至7 cm深度) 以及小麦生长状态实拍照片. 每个样方采集范围为10 m × 10 m, 在样方内随机选取4株小麦, 对植株样品进行称重, 以确定鲜重W1. 在实验室100℃的条件下对植株进行烘干, 获得干重W2. 本研究采用叶片含水量 (fuel moisture content,FMC)[19]表征冬小麦冠层含水量WFMC, 按照公式 (1) 进行计算. 小麦实拍照片用来获取实测小麦覆盖度, 采用阈值分割法[20-21]进行分类, 共计205张, 不同时期小麦实拍照片及分类结果如图3所示.
图3 实测覆盖度示意图Fig. 3 Schematic diagram of measured wheat coverage
植被覆盖度(fraction vegetation cover, FVC)是指植被 (包括枝、茎、叶) 在单位面积内垂直投影面积所占百分比[22-23]. 获取高精度和高时空分辨率的小麦覆盖度是衡量小麦长势的重要依据, 也是本次研究中影响小麦冠层含水量反演精度的重要参数. 像元二分模型是一种简单实用的混合像元分解模型, 在遥感估算植被覆盖度方面得到了广泛的使用[24]. 模型假定像元中只有植被和土壤组分, 其中植被占像元的百分比即为该像元的植被覆盖度, 采用公式 (2) 计算小麦覆盖度:
式(2)中:WFVC为小麦的植被覆盖度,S为某一植被指数值,Sveg和Ssoil分别为纯植被和裸土区的植被指数值.
常用于计算植被覆盖度的植被指数为NDVI (normalized difference vegetation index), 然而,NDVI在低植被覆盖区对土壤背景非常敏感, 并且, 近红外与红光波段的反射率比值形式容易在植被茂盛时期饱和, 造成估算误差. 因此, 利用NDVI计算得到的生长初期的小麦覆盖度偏高[25]. 研究表明,MSAVI(modified soil adjusted vegetation index)可以将土壤背景的影响减低, 增强对植被的敏感性[26].Huete等[27]提出增强植被指数EVI (enhanced vegetation index), 通过加入蓝光波段反射率, 采用非比值形式计算, 从而降低饱和. 处在不同物候期的小麦, 长势差异明显, 覆盖度变化显著. 为提高小麦整个生长期植被覆盖度的反演精度, 需寻求最佳的植被指数. 本研究选取NDVI、MSAVI和EVI共3种植被指数进行测试, 并以实测覆盖度验证估算精度, 得到最精确的反演结果, 用作小麦冠层含水量反演的关键参数. NDVI、MSAVI和EVI植被指数计算公式分别为
式(3)—(5)中:ρR、ρNIR、ρB分别为红光、近红外和蓝光波段反射率.
Sveg和Ssoil的取值是影响植被覆盖度反演精度的另一个重要因素. 在无实测数据的情况下, 最常用的方法是取给定置信度的置信区间内的最大值和最小值. 在本研究中, 涉及到小麦多个生长期. 在低植被覆盖期 (1—3月), 麦田稀疏, 直接取置信区间的最大值作为Sveg,Sveg偏低, 将导致覆盖度估算值偏高; 而高植被覆盖期 (4—5月), 麦田长势旺盛, 几乎达到郁闭, 取置信区间的最小值作为Ssoil,Ssoil偏高, 将导致覆盖度估算值偏低. 直接采用常规取值方法难以得到合理的估算结果. 因此, 本研究采用双时期法[28]分别确定Landsat-8和Sentinel-2数据的Sveg和Ssoil, 依据物候特征将小麦生长期划分为低植被覆盖期 (1—3月) 和高植被覆盖期 (4—5月). 1月22日和2月27日的小麦处在越冬期,影像中存在大量裸土区, 用于确定纯裸土像元的植被指数值比较理想, 根据频率累积表, 取频率为2%的植被指数值为Ssoil. 4月12和4月28日的小麦处在拔节期, 田间郁闭, 用于确定纯植被像元的植被指数值比较理想, 根据频率累积表, 取频率为98%的植被指数值为Sveg.
利用田间实测数据 (覆盖度、含水量), 建立植被水分指数与麦田实测整体含水量之间的回归关系. 构建反演方程, 通过遗传算法全局优化求解. 得到冬小麦不同时期的冠层含水量, 并验证其精度.最后, 将不同时期的冬小麦冠层含水量反演结果制图显示.
2.2.1 植被水指数的计算
Gao[2]提出基于近红外和短波红外的归一化水指数NDWI, 可有效估算植被含水量, 被广泛用于作物冠层含水量反演研究. 本文基于NDWI开展冬小麦的冠层含水量反演, 计算公式为
式(6)中:ρNIR、ρSWIR分别为近红外和短波红外波段反射率, 分别对应Landsat-8中的B5和B7,Sentinel-2中的B8和B12.
2.2.2 反演方程构建
小麦生长前/中期覆盖度较低, 在同一像元中, 土壤也有一定反射贡献, 为减少土壤对反演的影响,需要量化土壤贡献的遥感信息. 本研究将植被覆盖度这一参数引入到小麦冠层含水量的反演过程中.基于实测的小麦冠层含水量Wveg、土壤含水量Wsoil和麦田植被覆盖度FVC, 根据式 (7) 计算麦田的混合含水量W*, 并且建立NDWI与W*的回归关系. 由图4 (a) 可以发现, NDWI与实测小麦冠层含水量之间并无直接关系, 但与麦田的混合含水量W*具有较强的相关性,R2= 0.928. NDWI与W*的拟合方程可以表示为式 (8).
图4 NDWI与小麦冠层含水量散点图以及与麦田混合含水量的回归关系Fig. 4 Scatterplot of NDWI vs. wheat canopy water content, and regression relationship of NDWI vs. wheat canopy water content and mixed moisture content of wheat field
式(7)和(8)中:W*表示麦田实测值的混合含水量,Wveg为田间实测小麦冠层含水量(canopy water content, CWC),Wsoil为田间实测土壤含水量.
将式 (8) 代入式 (7) 中可得反演方程式 (9).Mveg和Msoil为需要求解的参数, 分别代表小麦冠层含水量和土壤含水量, FVC可由像元二分模型估算得到.
2.2.3 遗传优化求解反演方程
为求解冬小麦冠层含水量反演方程, 本研究考虑采用遗传算法进行全局优化搜索. 遗传算法于1975年被Holland[29]首次提出, 算法核心主要为遗传操作算子, 有选择、交叉、变异算子. 首先对优化的变量进行个体编码, 然后通过种群进化优化相应的目标函数, 得到最优解的个体. 本研究反演方程遗传优化求解具体流程见图5.
图5 遗传优化求解反演方程流程Fig. 5 Technological process for solving the inversion equation by genetic optimization
因此, 需先确定优化的目标函数, 如式 (10)—(12):
式(10)—(12)中:Merror表示遥感影像每个像素的W*和预测值M*差值绝对值, 在此基础上构建目标函数WFitness, 并在后续优化中最大化.
在构建目标函数的基础上, 对求解的小麦冠层含水量和土壤含水量参数进行二值编码, 并合并进一个个体. 对每个个体进行适应度目标函数值的计算, 每个个体对应一组含水量值, 通过式 (12) 可以计算出适应度值, 并保留具有最大适应度值的个体. 其后进行进化操作, 具体包括交叉和变异操作, 可以对种群种的个体进行重组, 生成两倍于上一代个体数目的新种群. 为保持种群数目稳定, 再次进行种群内每个个体的适应度值计算, 并保留前NP个个体形成新种群进入下一步迭代, 直至达到最大迭代次数, 输出当前保留的最大个体, 作为最优的小麦冠层含水量值.
通过双时期法取定像元二分模型的参数值, 如表2所示, 既考虑了小麦生长的连贯性又兼顾了不同的影像特征, 得到了较准确的反演结果. 用实测冬小麦的植被覆盖度值进行验证. 如图6所示, 用于计算植被覆盖度的3种指数 (NDVI、MSAVI、EVI) 能得到较好的反演结果,R2分别为0.898、0.773和0.931, RMSE分别为0.086、0.144和0.070, 均能有效模拟研究区不同生长期冬小麦植被覆盖度的变化规律. 其中, 基于EVI反演小麦覆盖度的结果最佳, 其反演结果将用于后续小麦冠层含水量的反演.
表2 像元二分模型参数取值Tab. 2 Parameter values for the pixel bisection model
图6 基于不同植被指数的冬小麦覆盖度的反演结果Fig. 6 Inversion results of winter wheat coverage with different vegetation indices
反演结果如图7所示, 从图中可以发现在冬小麦生长周期内, 冠层含水量处于50% ~ 90%的范围内, 随着生长周期呈先增长后降低的趋势, 与冬小麦的田间实际生长状况一致. 在小麦生长初期(1—3月), 冠层含水量整体在60% ~ 70%之间, 此时小麦正处于越冬/返青期; 小麦生长最旺盛的时期为4月, 此时小麦正处于拔节期, 冠层含水量整体大于80%; 在灌浆、成熟期 (5月), 小麦水分流失,冠层含水量回落, 整体在70%左右. 对于不同生长期的小麦, 整体冠层含水量反演值与实测值均有较好的一致性, 在小麦的整个生长周期均能有效监测其冠层含水量. 值得注意的是, 在大片种植的田块内部, 小麦冠层含水量反演结果稳定, 而田块边缘存在部分异常值. 考虑到在小麦种植区提取时, 未能精确地将田间小路剔除, 导致部分麦田边缘的反演结果异常. 为评估反演结果的可靠性, 采用田间实测值验证反演精度. 由图8 (b) 可见, 小麦总体反演结果RMSE = 5.6%、R2= 0.567, 与直接利用NDWI指数进行反演的结果RMSE = 31%、R2= 0.432相比, 本文方法误差可降低20%以上, 证明本文反演结果具有一定的可靠性.
图8 小麦冠层含水量估测值与实测值的对比Fig. 8 Comparison between the ground-truth CWC and the predicted CWC for winter wheat
将本文方法与其他文献中的方法进行横向对比, 分别从所提方法中的监测面积、研究数据类型及时期数、小麦生长期、反演精度、植被含水量表示类型和验证点数目共6个方面进行比较, 如表3所示. 可以发现: 首先, 从研究时相来看, 其他方法大多为单时期或部分时期的反演, 且多分布在拔节期以后, 对小麦越冬期和返青期的探究较少. 对于易发生春旱的华北平原而言, 此时期冬小麦需水量大,需对麦田进行实时监测. 本研究几乎涵盖了整个生长期, 能够对小麦进行及时有效的监测和指导灌溉,因此较其他方法更有优越性. 其次, 从反演精度来看, 需结合监测数据的时相数和时间、植被含水量表示类型及验证点数目进行综合评判. 本文与同一反演量纲 (FMC) 的研究进行对比, 虽研究精度略低,但对比方法的研究时相集中在小麦生长旺盛的拔节期、灌浆期和成熟期, 受土壤影响较小, 精度相对较高. 此外, 从图8 (b) 中可以发现反演精度较低的样本点属于小麦生长的前期, 受土壤影响较大. 若只考虑对比方法的研究时相 (拔节期、灌浆期和成熟期), 经统计本研究反演结果为R2= 0.794、RMSE = 3.2%, 较其他方法更优. 因此, 综合来看, 本文为研究冬小麦全生长期冠层含水量反演提供了一种有效的监测方法. 但本研究也存在一些不足之处, 小麦种植区的提取精度及覆盖度的估算结果对冠层含水量反演影响较大. 若能进一步提高两者精度, 可以获得更精准的冠层含水量反演结果. 此外, 受限于影像的质量和数量, 采用landsat-8和sentinel-2数据进行同步反演虽然可以提高时间分辨率, 但两种影像空间分辨率的不同也可能对反演结果造成一定影响.
本文基于多时相的Landsat-8 OLI和Sentinel-2 MSI多光谱数据, 以华北地区的冬小麦田为研究区, 结合不同生长时期的小麦覆盖度与NDWI构建反演方案, 提出基于光学遥感影像反演冬小麦生长期冠层含水量的反演方案. 得到如下结论:
(1) 基于遗传算法全局优化求解具有一定的优势. 本文所构建的冬小麦生长期冠层含水量反演模型, 可以有效进行冬小麦多时期大范围的冠层水分监测.
(2) 本文基于不同植被指数探究了多时期小麦覆盖度反演问题. 研究表明, 采用NDVI、MSAVI、EVI作为像元二分模型的植被信息均能得到较好的反演结果. 其中, 利用EVI反演的结果最为理想,能有效模拟研究区不同生长期冬小麦植被覆盖度的变化规律.
(3) 消除不同生长季土壤比例对小麦冠层含水量的反演精度影响显著. 这对于探索利用中高分辨率光学遥感影像, 开展冬小麦整个生长周期冠层含水量的反演, 进而对小麦作物的生长进行实时监测,都具有重要的实际意义.