杜蒙蒙,Ali Roshanianfard,刘颖超
1. 河南科技大学农业装备工程学院,河南 洛阳 471003 2. Faculty of Agriculture and Natural Resources, University of Mohaghegh Ardabili, 56199-11367, Ardabil, Iran 3. 上海交通大学机械与动力工程学院,上海 200240
小麦是我国第二大粮食作物。2020年我国小麦产量达1.34亿吨,占全球小麦总产量的17.68%。随着选育优良品种、土地整理、优化水肥运筹等技术的推广实施,近年来我国的小麦产量得到了大幅提升。国家统计局数据显示,2010年—2020年我国小麦产量在波动中保持增长趋势; 数据同时表明,我国小麦产量的逐年增幅比例已由高速增长转入低速增长阶段,持续提高小麦产量的任务十分艰巨[1-2]。另一方面,随着工业化和城镇化深入发展,近年来我国小麦种植面积逐年减少,2020年我国小麦种植面积约2338万公顷,较上一年度减少34.8万公顷,小麦进口依存度呈逐年增长态势,亟须提高小麦产量以保障粮食供应能力。小麦产量受品种遗传特性、环境条件及栽培措施等重因素影响[3],由有效穗数、穗粒数和粒重等共同决定。有效穗数是影响小麦产量的主导因素[4],适当增加小麦分蘖数是提高小麦有效穗数、改善小麦群体茎蘖结构的有力措施,与小麦高产、稳产密切相关。小麦茎蘖数及有效分蘖成穗数的多少与栽培技术措施紧密相关。在小麦分蘖期内(Growth Stage 21~29)适时、适量追施氮肥可显著促进分蘖的发生与发育[5]。然而,冬小麦分蘖过程时间跨度长,田块尺度上小麦茎蘖密度(单位面积内的小麦茎蘖数)的空间差异显著。忽视小麦茎蘖密度的空间差异状况,经验性地实施均一化追施氮肥作业往往导致氮肥利用率差与肥效偏低: 茎蘖密度较低的区域因肥力不足造成麦株长势瘦弱、分蘖少; 而茎蘖密度偏高的区域则因氮肥过度施用造成麦株过剩发育、田间郁蔽、有效分蘖成穗率降低,易倒伏[6-7]。因此,基于小麦分蘖的实际发生与发育状况,适时适量进行变量施氮管理,对小麦分蘖促、控同抓,可有效调控小麦茎蘖群体,提高小麦产量,实现化肥减量增效的目标。在小麦分蘖期内通过变量施氮作业对小麦茎蘖群体进行精准调控的技术关键,在于准确获取田块尺度的小麦茎蘖数空间分布数据。传统的通过人工田间调查统计获取小麦茎蘖密度的方法,在时效性与精准度上存在不足,工作量大、效率低,而且稀疏的点源统计数据往往无法精准反映田块内部的小麦茎蘖密度的空间差异信息。然而,图像处理与遥感技术为高效获取大田作物长势信息提供了技术途径。
国外,Jin等[8]在小麦出苗期(Growth Stage 10~18),将SONY ILCE α5100数码相机搭载在无人机上进行超低空遥感作业(航高3, 5与7 m),获取小麦苗期的高分辨率RGB遥感图像,通过二值化处理与图像分割,在毫米级分辨率的遥感图像中识别出小麦特征,建立并验证了小麦出苗密度的估算模型。Sripada等[9]在冬小麦起身期(Growth Stage 30),利用无人机遥感平台获取彩色-近红外遥感图像,针对冬小麦冠层生物量与多种植被指数分别进行了相关性分析,最终得出最佳追施氮肥剂量与代表冬小麦冠层生物量的相对绿色植被指数之间的二次回归模型,为冬小麦追施氮肥提供理论指导。Basso等[10]在冬小麦起身期(Growth Stage 30), 进行变量追施氮肥作业,并在施肥前后分别获取RapidEye卫星遥感图像,基于卫星图像提取出24种不同的植被指数,结合利用地面设备测量的小麦冠层叶绿素含量,对氮肥的有效利用率进行了分析。Walsh等[11]在冬小麦起身期(Growth Stage 30),利用手持式光谱仪GreenSeenker获取麦田光谱数据,计算出归一化差分植被指数值,并根据土壤水分含量对植被指数值进行修正,进一步提高了小麦产量预测模型的精度。Rezzouk等[12]在小麦抽穗期与乳熟期(Growth Stage 55,75), 分别获取航空遥感图像与近地拍摄图像,探讨了不同遥感方式在小麦表型研究中的适用性,结果表明在小麦乳熟期获取的基于航空遥感图像的植被指数与小麦产量具有较高的相关性。
国内,张猛等[13]在冬小麦返青期与起身期(Growth Stage 30),利用手持式作物光谱测量仪在高于小麦冠层50 cm处获取麦株光谱数据,结合4波段光谱数据计算并分析了9种植被指数与小麦茎蘖数之间的回归关系,为便携式小麦茎蘖数在线检测仪器设计提供了依据。吴焕丽等[14]在冬小麦返青期至抽穗期(Growth Stage 30~55), 手持Cannon EOS 7D数码相机在距离地面1.3 m处获取小麦冠层图像,基于改进K-means图像分割算法提取出麦田的植被覆盖度。黄芬等[15]在冬小麦拔节期(Growth Stage 31),手持Olympus E620数码相机获取小麦冠层图像,并将颜色空间下的基础颜色分量与小麦叶面含氮量进行了回归分析,为构建小麦叶面含氮量估算模型提供了方法借鉴。刘明星等[16]在冬小麦拔节期至乳熟期(Growth Stage 31~75),利用手持高光谱设备FieldSpec获取冠层高光谱数据,结合简单产量预测模型对冬小麦叶面积指数进行遥感检测与动态模拟,提高了冬小麦地上生物量的估算精度。Jiang等[17]在冬小麦拔节期、抽穗期与乳熟期(Growth Stage 31, 55, 75),利用无人机平台获取多光谱遥感图像,通过对遥感图像进行去噪、渐晕与辐射校正等预处理,建立了冬小麦的叶面积指数反演模型。杜颖等[18]在冬小麦乳熟期(Growth Stage 75),手持SONY DSC-H9数码相机在距离小麦顶部1米处获取小麦冠层图像,经超像素分割、样本标记、植被指数分析等处理,建立麦穗自动计数模型,为小麦大田估产提供参考。
综上所述,国内外关于冬小麦长势的遥感研究在生育期上主要集中在起身期(Growth Stage 30)及之后,在研究内容上主要集中在预测产量与反演叶面氮含量,在研究方法上主要集中在手持光学设备进行近地检测,小麦茎蘖密度预测精度普遍低于100株·m-2,难以满足变量追施氮肥作业的精度要求。在小麦分蘖期内,对小麦茎蘖密度进行田块尺度上的遥感检测与定量反演的研究较为少见。因此,为满足变量追施氮肥作业对田块尺度的小麦茎蘖密度专题图的需求,本研究使用大疆Mini 2航拍无人机,在小麦分蘖期获取试验田的可视光波段遥感图像,基于无人机定量遥感、图像处理、植被指数分析和人工神经网络预测等技术,在田块尺度上完成小麦茎蘖密度的定量反演,为变量追施氮肥作业的高效实施提供数据支持。
如表1所述,本研究采用大疆Mini 2无人机,搭载CMOS (complementary metal oxide semiconductors)传感器,在小麦分蘖期获取可见光波段的遥感图像。小麦试验田位于河南省开封市(北纬34°44′—35°01′,东经114°41′—115°15′),年均降水量670 mm,年有效积温4 592 ℃,无霜期219 d,属半干旱大陆性季风气候。试验田为多年耕种熟地,长约50 m,宽约10 m[如图1(b)红色方框标志部分所示],土壤成土母质为黄河冲积物,以潮土为主,肥力中等,前茬作物为花生。供试小麦品种为“矮抗58”,2020年10月9日使用小麦八行(行距为0.2 m)单箱播种机完成小麦播种,播种量为165 kg·hm-2。基肥为磷酸二铵、氯化钾(施用量均为225 kg·hm-2)与尿素(150 kg·hm-2)混合施用,于小麦播种前撒施后翻埋。
表1 试验设备Table 1 Experiment equipment
无人机遥感试验于2020年12月15日进行,无人机飞行速度约为5 m·s-1,飞行高度约40 m,航拍相机的快门速度77.78 ms,ISO100,获取影像空间分辨率为8 mm,具有蓝、绿、红3个波段。采用“1米双行”定点测量的方法[如图1(c)所示],调查试验田内部距离麦垄1 m左右处、位于田畦中心部位的[如图1(b)红色短线所示]、具有上不少于3个叶片的茎蘖的数量,根据麦田20 cm行距以及式(1)计算单位面积内的小麦茎蘖数(小麦茎蘖密度)地面真值,共计20处,如表2所示。
图1 试验区域及小麦茎蘖密度地面真值调查点的位置(a): 试验地块位置; (b): 小麦试验田及小麦茎蘖密度地面真值调查点; (c): 调查“1米双行”区域内的小麦数Fig.1 Study area and positions of ground truths of wheat tiller density(a): Experimental field; (b): Experimental field and positions of ground truths of wheat tiller density;(c) Counting wheat tiller numbers within the zone of “1 meter and 2 rows”
表2 小麦茎蘖密度地面真值与相应的VI值Table 2 Ground truths of wheat tiller density and corresponding VI values
Td=5×Tm
(1)
式(1)中,Td与Tm分别为小麦茎蘖密度与 “1米双行”人工调查区域内的小麦茎蘖数量。
商业彩色相机的CMOS成像传感器普遍具有350~1 100 nm的光谱响应能力,相机基于拜耳模式滤光阵列,将摄入镜头的地物光谱信息划分至蓝、绿、红三个通道,并通过光电感应装置将光信号转换为电信号,量化电脉冲信号(0~255),输出具有蓝、绿、红3个波段光谱信息的彩色图像。图像校正是利用遥感图像进行定量分析的首要步骤,根据成像传感器响应特性进行图像校正,可以修正光学镜头的非均匀性引起的径向畸变与切向畸变。使用Matlab相机标定工具箱,拍摄棋盘格标定板获取相机的内部参数矩阵(像素尺寸、图像原点相对于光心成像点的纵向偏移量及横向偏移量)、外部参数矩阵(世界坐标系转换至相机坐标系的旋转变换矩阵与平移变换矩阵)、以及相机的径向畸变系数与切向畸变系数,最终完成无人机遥感图像校正,并提取蓝、绿、红3个波段的遥感图像,分别如图2(a),(b)和(c)所示。
图2 校正后的无人机遥感图像(蓝、绿、红3波段)(a): 蓝色波段; (b): 绿色波段; (c): 红色波段Fig.2 Calibrated drone remote sensing images of blue, green and red bands(a): Blue band; (b): Green band; (c): Red band
植被指数(vegetation index, VI)将两个或两个以上波段的地物反射率或辐射值进行组合运算,以增强植被的某一特性或细节。在定量遥感领域,VI被广泛应用在作物长势监测与植被覆盖度(fractional vegetation coverage, FVC)估算等方面。1960年代末Jordan基于多光谱图像的近红外波段与红光波段的光谱反射信息,构建比值植被指数(ratio vegetation index, RVI),用于评估森林生物量。针对RVI在FVC低于50%时受大气影响较大的问题,Rouse等在1970年初基于近红外波段与红光波段的光谱反射信息,提出归一化差分植被指数(normalized difference vegetation index, NDVI)。截止目前,有超过150余种VI被应用在遥感分析的各个领域,而NDVI的应用最为广泛,通常被用于进行叶面积指数定量反演、植被光合作用强度预测等方面[19]。然而,NDVI对土壤亮度较为敏感,在FVC较大时易出现溢出效应。另外,获取近红外波段光谱信息往往依赖于多光谱传感器,其摄取的各波段影像需经严格配准与标定之后才能用于构建VI,而且多光谱传感器价格昂贵,严重阻碍其在农业定量遥感领域的推广应用。
另一方面,商业彩色相机在可见光波段范围内(400~800 nm)针对健康的绿色植被与土壤背景具有较强的光谱辨识能力(如图3所示)。植被在绿色波段范围内具有较高的反射率,而在红色波段范围内具有明显较低的反射率; 与之相反,土壤在绿色波段范围内具有较低的光谱反射率,而在红光波段范围内具有较高的反射率。利用植被与土壤在不同波段范围内的反射率差异,构建VI并提取无人机遥感图像中的大田植被特征在理论上具有较高可行性。基于植被与土壤在可见光波段的光谱响应特性,选取可以较好地突出植被特征、减轻光照强度不均对遥感数据造成影响的四种比值类型VI,具体为可见光波段差分植被指数(visible-band difference vegetation index,VDVI)、归一化绿红差分指数(normalized green-red difference index,NGRDI)、归一化绿蓝差分指数(normalized green-blue difference index,NGBDI)、绿红比值指数(ratio green-red index,RGRI),用以定量反演小麦茎蘖密度指标,分别如式(2)—式(5)所示。
图3 可见光波段范围内植被与土壤的光谱响应特性Fig.3 Spectral response characteristics of vegetation and soil in visible light domain
VDVI=(2G-R-B)/(2G+R+B)
(2)
NGRDI=(G-R)/(G+R)
(3)
NGBDI=(G-B)/(G+B)
(4)
RGRI=G/R
(5)
式中B,G,R分别为无人机遥感图像中蓝色波段、绿色波段与红色波段的地物反射率。
FVC指植被在生长区域的地面垂直投影面积占研究区域面积的百分比,是反映小麦茎蘖个体/群体长势的重要指标。遥感检测是大范围获取FVC的主要技术手段,选取与FVC具有良好相关关系的VI, 并建立VI与FVC的转换关系。因此,本文选用VDVI,根据式(6),计算并获取小麦分蘖期的FVC专题图,用于小麦茎蘖密度测算。
FVC=(VDVI-VDVIb)/(VDVIv-VDVIb)
(6)
式(6)中,VDVI,VDVIb及VDVIv分别为任意像素点的VDVI值、表征土壤背景的VDVI值(取在小麦试验田内部人工随机选取的10处土壤特征的VDVI平均值)及表征植被的VDVI值(本文取在小麦试验田内部人工随机选取的10处植被特征的VDVI平均值)。
由于小麦分蘖期茎蘖数增多,FVC过饱和趋势明显,单独使用FVC难以实现对小麦茎蘖数的高精度定量反演。因此,为提高小麦茎蘖数定量反演模型的精度与鲁棒性,基于小麦分蘖期的无人机遥感图像,利用VI法计算FVC,同时综合VDVI,NGRDI,NGBDI与RGRI四种植被指数,结合人工田间统计的小麦茎蘖密度地面真值数据,建立BP(backward propagation)神经网络预测模型,在田块尺度上完成小麦茎蘖密度的定量反演。
BP神经网络由Werbos于1974年提出,是一种按误差反向传播进行模型训练的多层前馈神经网络,具有非常强大的非线性映射能力,在满足合适数目的神经元的前提下,可无限逼近任意复杂度的非线性函数。BP神经网络通常由输入层、隐含层与输出层构成,分别具有n,q和m个神经元。
输入层向量为X=(x1,x2, …,xn)T,隐含层向量为Z=(z1,z2, …,zq)T,输出层向量为O=(o1,o2, …,om)T,期望输出向量为Y=(y1,y2, …,ym)T。输入层与隐含层、隐含层与输出层的连接权值、阈值分别为V=(v1,v2, …,vq)T与W=(w1,w2, …,wm)T。BP神经网络将训练样本从输入层经由隐含层传输至输出层,比较实际输出值与期望输出之间的误差,将误差进行反向传播,并根据误差信息调整连接各神经元的权值与阈值,若干次迭代后获取目标函数值小于预定期望值时停止训练,最终建立预测模型。
以小麦茎蘖密度地面真值调查点处1 m×1 m观测区域的FVC,VDVI,NGRDI,NGBDI及RGRI平均值作为5节点输入层,小麦茎蘖密度地面真值Td为单节点输出层,根据式(7)拟定隐含层神经元的节点数为10,建立一个单隐含层、多输入、单输出的3层BP神经网络预测模型,如图4所示。
图4 BP神经网络训练模型Fig.4 Training model of BP neural network
(7)
式(7)中,l,m,n,a分别为隐含层神经元的节点数、输入层神经元的节点数、输出层神经元的节点数为1~10之间的经验常数。
VI专题图指每个像素只有一个亮度值的栅格图像,像素亮度值通常根据遥感图像特定波段的反射率或亮度值计算获取。利用图2所示的无人机遥感图像蓝、绿、红3波段的图像分量,在ENVI遥感图像处理软件环境下,基于“波段运算”,根据式(2)—式(5)获取相应的VDVI,NGRDI,NGBDI与RGRI专题图。进一步地,基于VDVI专题图,在小麦试验田地块内部随机选择并提取10处植被特征与10处土壤特征的VDVI值,分别计算植被特征与土壤背景的VDVI平均值为0.389与-0.069; 根据式(6),通过“波段运算”,获取小麦分蘖期的FVC专题图,如图5所示。
图5 植被指数专题图Fig.5 Vegetation index maps
在小麦试验田内部随机选择并提取10处植被特征与10处土壤特征的FVC,VDVI,NGRDI,NGBDI及RGRI的像素值,获取植被与土壤特征在各VI专题图中的统计特征值,如表3所示。10处植被特征的FVC平均值与标准差分别为0.77与0.052; 而10处土壤特征的FVC平均值与标准差分别为0.158与0.016,两类图像特征在FVC专题图中的像素值具有较大差异。类似地,两类图像特征在VDVI,NGRDI,NGBDI及RGRI专题图中的像素值的平均值与标准差分别为0.389与0.071、-0.069与0.089; 0.204与0.026、-0.047与0.011; 0.375与0.043、0.059与0.008; 1.516与0.083、以及0.911与0.020。由表3及图 5可知,植被特征与土壤特征在各VI专题图中区分明显: 植被特征在各VI专题图中均具有较高的像素值,其色调在图中表现明亮; 相反地,土壤特征在上述VI专题图中均具有较低的像素值,其色调在图中表现灰暗。植被特征与土壤特征在FCV专题图中的像素值的平均值具有较大差异,两类图像特征的FVC平均值相差0.612。另一方面,VDVI对植被特征及土壤特征具有较强的区分能力,植被特征在VDVI专题图中为远大于零的正值,而土壤特征则为略小于零的负值。虽然NGRDI与VDVI具有类似的性质,但植被特征的NGRDI平均值小于VDVI平均值,表明其提取植被特征的能力相对较弱。另外,NGBDI提取植被特征的能力介于NGRDI与VDVI之间,但其土壤特征的平均值为略大于零的正值,表明其提取土壤特征的能力较弱。RGRI与FVC的性质相近,两类图像特征的RGRI平均值相差0.605。
表3 植被与土壤在FVC,VDVI,NGRDI,NGBDI与RGRI专题图中的统计特征值Table 3 Statistic values of vegetation and soil in FVC,VDVI,NGRDI,NGBDI and RGRI maps
因此,由植被特征与土壤特征在各VI专题图中的统计特性可知,特定区域内植被越多,具有较大VI值的像素点则越多,其区域VI值的总和(或平均值)则越大。针对图1(b)所示的小麦茎蘖密度地面真值调查点,分别统计各感兴趣区域(region of interest, RoI)内所有像素点的FVC,VDVI,NGRDI,NGBDI及RGRI值,计算VI平均值,进而定量反演小麦茎蘖密度值。
图6 BP神经网络模型训练结果Fig.6 Training result of BP neural network model
采用表2中的5组试验数据对基于BP神经网络的小麦茎蘖密度预测模型进行精度验证,调用已完成模型训练的神经元权值W1,W2与阈值V1,V2,输入5组FVC,VDVI,NGRDI,NGBDI及RGRI试验数据,获取小麦茎蘖密度的预测值,如表4所示。根据式(8)及式(9)可得小麦茎蘖密度定量反演的BP神经网络模型的均方根误差(root mean square error, RMSE)及平均绝对百分比误差(mean absolute percentage error, MAPE)分别为19及3.62%,因此,该BP神经网络模型具有较高的小麦茎蘖密度预测精度。
表4 小麦茎蘖密度地面真值与预测值Table 4 Ground truths and prediced values of wheat tiller density
(8)
(9)
式(9)中,Td与Tdp分别表示小麦茎蘖密度的地面真值与使用上述的BP神经网络预测模型获取的预测值。
进一步地,在Arcmap地理数据处理软件环境下,将小麦茎蘖密度反演试验田的各VI专题图进行1 m×1 m网格化处理,通过“非重叠块统计”,分别获取各网格的FVC,VDVI,NGRDI,NGBDI及RGRI像素平均值,共计150组。同样地,调用上述小麦茎蘖密度定量反演的BP神经网络预测模型,代入FVC,VDVI,NGRDI,NGBDI及RGRI,获取相应的150组小麦茎蘖密度预测值,在Arcmap地理数据处理软件环境下,生成田块尺度的小麦茎蘖密度定量反演专题图,如图7所示。
图7 田块尺度的小麦茎蘖密度定量反演专题图Fig.7 Map of quantitatively inversed wheat tiller density at field level
由图7可知,该试验田的小麦茎蘖密度最小值、最大值及平均值分别为455,865及675株·m-2。小麦茎蘖密度低于500株·m-2的地块面积占比6.67%,低于越冬前500万株·hm-2的要求,适当提高该区域的追施氮肥施用量可有效提高小麦茎蘖数量。小麦茎蘖密度介于501~800株·m-2之间地块面积占比74.67%,具有较为合理的小麦茎蘖群体结构,按照150 kg·hm-2的常规追施氮肥施用量即可保障后期具备足够数量的有效分蘖成穗数。而地块面积占比18.66%的小麦茎蘖密度高于800株·m-2,需适当降低追施氮肥的施用量,避免小麦个体发育过剩导致的群体结构恶化、田间郁闭、光合作用与氮代谢能力减弱,减低有效分蘖成穗率及整体产量。
田块尺度的小麦茎蘖密度定量反演建立在高精度的遥感图像获取与分析的基础上。利用航拍无人机获取大田图像,为在田块尺度上进行作物长势监测提供了有效途径。植被指数可以突出遥感图像中的植被特征,被广泛应用在叶面积指数及植被覆盖度的定量反演中。选用VDVI,NGRDI,NGBDI及RGRI四种比值类型的植被指数,以突出遥感图像中的植被特征,并减轻光照强度对遥感图像各波段像素值的影响,在田块尺度上对小麦分蘖密度进行定量反演,取得较高的预测精度。
然而,航拍相机可见光波段的遥感图像未经严格辐射定标,各波段的图像数据无法直接用于定量遥感的反演分析。此外,光照强度直接决定无人机图像各波段像素值的大小,并间接影响遥感图像的对比度与清晰度。但是,使用两个及两个以上波段的像素值或辐射亮度值之差获取的差值类型植被指数(过绿植被指数等),以及使用两个及两个以上波段的像素值或辐射亮度值的加权和类型植被指数(绿度植被指数),均不能解决光照强度对无人机遥感图像的像素值造成的影响,在小麦茎蘖密度定量反演的模型中选用此类型的植被指数将对模型的预测精度与模型的适应性造成不可避免的影响。
(1) 在小麦分蘖期内获取可见光波段的无人机遥感图像,建立植被指数与小麦茎蘖密度地面真值之间的预测模型,在田块尺度上定量反演小麦茎蘖茎蘖密度,取代精度差、效率低的传统人工田间调查方法,满足现代农业对大田信息获取的高效率、高精度的要求。验证试验数据表明: 小麦茎蘖密度定量反演BP神经网络模型的均方根误差及平均绝对百分比误差分别为19及3.62%,与小麦茎蘖密度地面真值的平均值635相比较,该模型具有较高的小麦茎蘖密度预测精度。
(2) 田块尺度的小麦茎蘖密度反演专题图的统计数据显示: 小麦茎蘖密度低于500株·m-2、介于501~800株·m-2之间、以及高于800株·m-2的地块面积分别占比6.67%,74.67%和18.66%,为变量追施氮肥提供数据支持。
利用商业航拍无人机获取小麦分蘖期内的可视光波段遥感图像,在田块尺度上实现对小麦茎蘖密度的定量反演,是实施变量追施氮肥作业、推进精准农业大田信息获取及农业数字化建设的现实需求。