向友珍 安嘉琪 赵 笑 金 琳 李志军 张富仓
(1.西北农林科技大学旱区农业水土工程教育部重点实验室, 陕西杨凌 712100;2.西北农林科技大学水利与建筑工程学院, 陕西杨凌 712100; 3.西北农林科技大学综合素质教育学院, 陕西杨凌 712100)
大豆是人类和动物饮食中显著的蛋白质和脂肪酸来源[1],对于众多亚洲国家,大豆对其饮食文化的影响有着十分重要的作用[2]。因此,确保大豆的高质量和高效率生产具有重要意义。
叶面积指数(Leaf area index, LAI)和地上部生物量是作物冠层结构的重要参数,不仅可为作物生长提供动态信息,还决定作物生物物理过程[3-4],同时作为重要的输入参数,在作物生长模型与决策支持系统中发挥着十分重要的作用[5]。通过对LAI和地上部生物量实测值的分析,能够进行作物生长监测和产量预测,也可以进行田间水肥管理[6-8]。因此,快速、准确地测定地上部分的LAI和地上部生物量具有相当重要的意义。当前获取作物LAI和地上部生物量的方法有两种,分别是直接测量法和间接测量法。其中前者的精确度较后者更高,但是直接测量法需要对作物进行破坏性采样,不仅费时费力并且所取样本不一定具有代表性,因此该方法有一定的局限性;与之不同的是,间接测量法是将作物光谱信息与田间实测数据相结合并通过模型对作物生长指标进行估算的方法,且具有快速、高效等特点。无人机多光谱遥感操作简易、灵活性高并且成本较低,这些特点都使其成为近年来进行精准农业研究的一项重要工具以及获取作物冠层参数的重要手段[9-10]。邵国敏等[11]通过配备Red Edge多光谱相机的六旋翼无人机,对不同灌溉条件下的玉米LAI进行了估算,发现在不同的灌溉条件下不同的模型对于估算玉米LAI的效果各不相同,其中在充分灌溉条件下,效果最好的是植被指数的多元线性回归模型,而在水分胁迫条件下,表现更好的则是植被指数的随机森林(Random forest,RF)模型。刘涛等[12]通过得到不同飞行高度的多光谱数据来估计小麦冠层LAI,得到采用偏最小二乘回归技术反演小麦LAI的精度最高。黄林生等[13]通过使用配备MicaSense Red EdgeTM3光谱相机的DJI M600 Pro六旋翼高性能无人机,对水稻LAI进行了估算,提出结合光谱特征、纹理指数和植物覆盖度等多个指标的估算模型具有最佳精度。
作物产量与人民生活水平和国家粮食安全息息相关[14],及时准确地预测粮食产量对农业管理、粮食政策制定、农业保险和农业用水效率评估具有重要的意义。目前,主要通过经验模型、半经验模型和物理模型3种方法来进行遥感估产[15-17]。而无人机遥感可获得比卫星和地面遥感更高的时间分辨率、空间分辨率和光谱分辨率影像[18-20],已成为农业遥感数据的重要途径[21]。无人机遥感的应用对于发展智慧农业、精准农业具有重要的影响。目前,对作物产量的估算已有部分研究。余坤勇等[22]根据光谱数据,利用不同的光谱指数进行玉米产量的预测;STEPANOV等[23]在利用无人机遥感数据进行估产的研究中对不同的光谱指数进行讨论,发现归一化差异植被指数可以更好地用于产量估算。王鹏新等[24]在估算小麦产量时用最小二乘法结合无人机高光谱数据与不同植被指数,发现最优估算指数为EVI2(无蓝带增强植被指数)。
过往的研究大多只通过单一植被指数或光谱反射率分析反演估算的LAI、地上部生物量或者产量,通过不同植被指数的组合作为输入变量进行建模的较少。此外,在利用无人机多光谱数据预测西北地区大豆产量的研究较少。因此本文以西北地区大豆为研究对象,通过无人机获取大豆的田间多光谱数据,分析各生育期植被指数、LAI、地上部生物量和产量之间的关系,采用RF、反向传播神经网络(Back propagation neural network, BPNN)和支持向量机(Support vector machine, SVM)3种机器学习方法建立定量模型。并讨论不同机器学习方法对大豆LAI、地上部生物量和产量估算模型精度的影响以及不同生育期的光谱信息对产量估算模型的精度影响,以期为更加准确、快速地获取大豆的生长指标和产量预测提供理论参考。
试验于2021年6—9月在西北农林科技大学旱区农业水土工程教育部重点实验室节水灌溉试验站(34°20′N,108°24′E,海拔521 m)进行(图1)。该地区气候类型为暖温带季风半湿润气候,是典型的半湿润偏(易)旱区。根据当地气象资料显示降水年际间差异较大(1995—2020年),平均降水量为561 mm,其中全年的降雨主要发生在7—9月,平均占年降水量的一半以上(约55.5%)。试验区土壤为重壤土,0~100 cm土层的田间持水率(质量分数,下同)为23%~25%,凋萎含水率为8.5%,平均土壤干容积密度为1.58 g/cm3。
图1 研究区概况Fig.1 Study area
本试验设置4种覆盖类型:平作秸秆覆盖(SM)、垄覆地膜沟覆秸秆(SFM)、垄覆地膜沟无覆盖(FM)和平作无覆盖方式。设置5个施氮量:0 kg/hm2(N0)、30 kg/hm2(N1)、60 kg/hm2(N2)、90 kg/hm2(N3)和120 kg/hm2(N4),共计20个处理。
设置7个种植密度:1.5×105株/hm2(R0)、2.0×105株/hm2(R1)、2.5×105株/hm2(R2)、3.0×105株/hm2(R3)、3.5×105株/hm2(R4)、4.0×105株/hm2(R5)和4.5×105株/hm2(R6),共计7个处理。
设置各处理小区时采用随机排列原则,并设置2个重复,共54个小区,小区面积为2.4 m×6 m=14.4 m2,试验区周围布设2 m保护带。每个小区磷肥和钾肥保持一致,均为30 kg/hm2。在2021年6月18日按行距50 cm、株距约10 cm人工点播大豆,其他田间生产管理(施药、除草等)均与当地保持一致。于9月30日收获大豆,并回收地膜。大豆全生育期为苗期(2021年6月18—25日)、分枝期(2021年6月26日—7月18日)、开花结荚期(2021年7月19日—8月17日)、鼓粒期(2021年8月18日—9月10日)与成熟期(2021年9月11—30日)。其中分枝期大豆正处于营养生长阶段,是决定整个生育期植株健壮与否、分枝与开花多少的关键时期[25];结荚期是大豆生长最旺盛的时期,营养生长与生殖生长并进,是生物量积累的关键时期[26];鼓粒期是大豆产量和品质形成的重要生育时期,大豆在鼓粒期进行蛋白质和脂肪的储藏[27]。而苗期和成熟期大豆各项生长参数变化不大,因此,本研究选择分枝期、结荚期和鼓粒期3个生育期获取无人机影像数据。
无人机遥感影像获取时间为2021年7月25日(四节期,V4)、2021年8月20日(结荚期,R4)和2021年9月5日(鼓粒期,R6)。采用大疆精灵4多光谱版一体化无人机(Phantom4-M,P4M)获取大豆冠层多光谱影像数据。该设备集成1个可见光传感器通道和5个多光谱传感器通道(蓝光、绿光、红光、红边和近红外),每次拍摄可获得6幅影像。且每幅影像具有2.0×106像素,最高飞行速度14 m/s,最大续航时间27 min,搭配有TimeSync时间同步系统,可获得厘米级定位精度,此外P4M顶部集成光强传感器,可捕捉太阳辐照度数据用于影像的光照补偿,排除环境光对数据的干扰,提高不同时段采集数据的准确度与一致性。试验利用DJITerra软件进行无人机航线规划,航向和旁向重叠率均为80%,影像分辨率分别为1.6 cm/像素。使用Pix4D mapper软件进行影像拼接,可获得不同飞行高度下的RGB和单波段正射影像。将5个波段的多光谱正射影像导入ENVI软件进行波段合成,使用快速辐射校正工具将像素DN值转换为反射率。
在无人机多光谱数据影像采集的当天进行大豆LAI和地上部生物量的数据采集。
LAI测定方法:使用美国LI-COR LAI-2200C型植物冠层分析仪测定。在54个试验小区内各随机进行6次大豆LAI数据采集,以其均值作为各试验小区的实测值。
地上部生物量测定方法:在每个小区选取6株能代表该小区长势的大豆进行破坏性取样,在105℃杀青30 min,随后75℃干燥后使用电子天平称量大豆总质量。最终将测量结果换算成单位面积的干物质量。
产量测定方法:在大豆成熟后,选取每个小区的中间4行测量其产量,并进行考种工作(随机选择10株),籽粒产量按含水率13%计[28]。
根据植被的光谱吸收特性,将遥感影像不同波段反射率之间进行线性或非线性组合形成的无量纲指标参数称为植被光谱指数。该指数是对地表植被状况简单、有效的经验度量,同时能够体现植被在可见光、近红外波段反射与土壤背景之间的差异[10]。植被的LAI、地上部生物量和产量等指标与差异化植被指数(DVI)、比值植被指数(RVI)、归一化植被指数(NDVI)、增强型植被指数(EVI)和土壤调剂植被指数(SAVI)等经验植被指数密切相关。因此,选取15种植被光谱指数构建模型,表达式如表1所示。
表1 经验植被指数及计算公式Tab.1 Empirical vegetation spectral indexes and calculation formulas
基于田间试验共获取162组LAI和地上部生物量样本及54组产量样本和三生育期的光谱数据样本(共162组)。对162组LAI和地上部生物量样本进行从小到大排序,同时也将54组产量样本进行从小到大排序,在这些样本中随机选取2/3的样本作为建模集,剩余1/3的样本作为验证集。表2为建模集和验证集的样本数量及产量的统计特征。
表2 大豆LAI、地上部生物量和产量描述性统计Tab.2 Descriptive statistics for leaf area index, above-ground biomass and yield of soybean
基于计算得到的15种光谱指数,分析全生育期光谱指数与LAI、地上部生物量的相关性,筛选出相关系数(R)最高的5个植被指数,作为模型的输入变量。采用SVM、RF、BPNN进行建模,对大豆的LAI和地上部生物量进行回归预测。同理,将各生育期的光谱指数与产量数据进行相关性分析,在各生育期筛选出与产量相关系数(R)最高的5个植被指数,作为模型的输入变量。采用SVM、RF、BPNN进行建模,对大豆的产量进行回归预测。
支持向量机模型的核函数类型设置为“poly”,根据最小交叉验证误差原则分别确定参数惩罚系数C和γ为20和0.02;在RF模型的构建中,经过参数优化和多次训练后,LAI模型和地上生物量模型中的决策树数均设置为600。反向神经网络模型隐蔽层传递函数设置为“TANSIG”,基于数值优化理论的Levenbeger-Marquardt(Train-LM)算法被用作网络训练函数,经过多次训练,确定中间层的神经元数量为15个[39]。
采用均方根误差(RMSE)、决定系数(R2)与平均相对误差(MRE)对模型拟合结果进行评价[37],模型的预测精度与决定系数(R2)成正相关;RMSE与MRE越小,说明模型的性能更稳定,预测结果也更为集中。
利用数据集中实测的LAI值分别与15个植被指数进行相关性分析,计算结果如表3所示。其中,RVI的相关系数最高,为0.706;由大到小排序分别筛选出5个植被指数作为LAI估算模型的输入,分别为RVI、NDVI、GNDVI、NDRE和CIRE。地上部生物量分别与15个植被指数的相关性分析结果如表3所示,OSAVI与地上部生物量的相关系数最高,为0.676,同理,选择OSAVI、SAVI、GOSAVI、EVI和MTCI作为地上部生物量估算模型的输入。
表3 植被指数与LAI和地上部生物量相关系数计算结果Tab.3 Calculation of correlation coefficients between vegetation index and leaf area index and above-ground biomass
各生育期15组植被指数与产量的相关性分析结果如表4所示,其中VIopt在V4期与产量相关性最佳,相关系数高达0.813;而在R4和R6生育期与产量相关性最佳的植被指数为CIRE,相关系数分别为0.849和0.827。与LAI和地上部生物量筛选植被指数的方法同理,在V4选取植被指数VIopt、RVI、GCI、TVI和MTVI作为该生育期的产量估算模型的输入;选取植被指数CIRE、NDRE、GCI、MTCI和GNDVI作为R4的产量估算模型的输入;选取植被指数NDRE、CIRE、GCI、MTCI和GNDVI作为R6的产量估算模型的输入。
表4 植被指数与产量相关系数Tab.4 Correlation coefficient between vegetation index and yield
以2.1节筛选出的植被指数组合作为自变量,以大豆LAI作为响应变量,分别采用SVM、RF、BPNN构建大豆全生育期LAI估算模型,从R2、RMSE、MRE 3方面综合评定模型精度,不同建模方法对于大豆叶面积的预测结果如图2所示。结果表明:基于SVM构建的大豆LAI估算模型建模集和验证集的R2分别为0.506和0.593,RMSE分别为1.054 m2/m2和1.006 m2/m2,MRE分别为26.454%和23.563%; 基于RF构建的大豆LAI估算模型建模集和验证集的R2分别为0.783和0.801,RMSE分别为0.713 m2/m2和0.675 m2/m2,MRE分别为17.693%和18.684%; 基于BPNN构建的大豆LAI估算模型建模集和验证集的R2分别为0.588和0.708,RMSE分别为0.961 m2/m2和0.756 m2/m2,MRE分别为23.564%和19.741%。
图2 基于SVM、 RF和BPNN的大豆LAI估算模型Fig.2 Estimation model of soybean LAI based on reverse SVM, RF and BPNN
同理,不同建模方法对于大豆地上部生物量的预测结果如图3所示。结果表明:基于SVM构建的大豆地上部生物量估算模型建模集和验证集的R2分别为0.481和0.451,RMSE分别为2 296.584 kg/hm2和 2 213.534 kg/hm2,MRE分别为28.880%和26.649%;基于RF构建的大豆地上部生物量估算模型建模集和验证集的R2分别为0.765和0.744,RMSE分别为1 550.173 kg/hm2和1 548.140 kg/hm2,MRE分别为21.673%和18.769%; 基于BPNN构建的大豆地上部生物量估算模型建模集和验证集的R2分别为0.557和0.536,RMSE分别为2 096.711 kg/hm2和1 973.511 kg/hm2,MRE分别为27.127%和25.642%。
图3 基于SVM、RF和BPNN的大豆地上部生物量估算模型Fig.3 Estimation models of soybean above-ground biomass based on reverse SVM, RF and BPNN
以2.2节筛选出的各生育期植被指数组合作为输入变量, 以大豆产量作为响应变量,分别采用SVM、RF、BPNN构建大豆各生育期产量估算模型,从R2、RMSE、MRE 3方面综合评定模型精度,不同建模方法对于大豆产量的预测结果如图4所示。
图4 基于SVM、RF和BPNN的大豆产量估算模型Fig.4 Estimation models of soybean grain yield based on reverse SVM, RF and BPNN
结果表明:在V4生育期,基于SVM构建的大豆地上部生物量估算模型建模集和验证集的R2分别为0.634和0.615,RMSE分别为289.474、388.397 kg/hm2,MRE分别为6.894%和9.735%;基于RF构建的大豆地上部生物量估算模型建模集和验证集R2分别为0.781和0.741,RMSE分别为225.671、332.876 kg/hm2,MRE分别为5.433%和7.731%;基于BPNN构建的大豆地上部生物量估算模型建模集和验证集的R2分别为0.752和0.709,RMSE分别为238.264 kg/hm2和338.495 kg/hm2,MRE分别为5.529%和7.990%。
在R4生育期,基于SVM构建的大豆地上部生物量估算模型建模集和验证集的R2分别为0.704和0.679,RMSE分别为268.616、341.065 kg/hm2,MRE分别为4.224%和7.246%; 基于RF构建的大豆地上部生物量估算模型建模集和验证集的R2分别为0.847和0.818,RMSE分别为188.617、287.539 kg/hm2,MRE分别为3.707%和7.128%; 基于BPNN构建的大豆地上部生物量估算模型建模集和验证集的R2分别为0.805和0.792,RMSE分别为271.815、296.111 kg/hm2,MRE分别为4.937%和7.399%。
在R6生育期,基于SVM构建的大豆地上部生物量估算模型建模集和验证集的R2分别为0.687和0.636,RMSE分别为281.901、372.816 kg/hm2,MRE分别为6.191%和8.432%;基于RF构建的大豆地上部生物量估算模型建模集和验证集R2分别为0.802和0.776,RMSE分别为214.021、313.524 kg/hm2,MRE分别为4.808%和8.983%;基于BPNN构建的大豆地上部生物量估算模型建模集和验证集的R2分别为0.773和0.741,RMSE分别为227.772、315.486 kg/hm2,MRE分别为5.617%和8.916%。
目前无人机多光谱遥感技术已经较为成熟,且在作物生理生态及生化参数反演方面的应用也越来越广泛。当前使用无人机多光谱遥感技术对作物的生长生理指标的估算研究主要是通过计算植被指数并以此构建统计模型或机器学习模型,不同研究地域、植被种类和估算方法得到的研究结果基本相同[40-43],但精度存在一定的差异。
本文通过对植被指数与目标参数的相关性分析,发现植被指数GCI与LAI、地上部生物量和产量的相关性均较好,这可能是因为GCI对植物中叶绿素含量变化的响应比较敏感[12]。叶绿素含量可以反映植被的生理状态,叶片的生长及生物量的积累需要叶绿素的驱动进而进行光合作用[44],这直接导致叶绿素与LAI、地上部生物量和产量的相关性均较高,所以通过GCI的变化能很好地实时监测大豆的长势。本研究还发现,随着生育期的推进,大部分植被指数与产量的相关性呈先增大后减小的趋势,同时,在输入变量和机器学习模型相同的情况下,R4的产量预测模型精度更高,这可能是因为R4生长时期是作物生长最旺盛的时期,其LAI、叶绿素含量和光合速率都达到峰值,而在V4,由于大豆还处于营养生长阶段,作物体内叶绿素含量水平较低,光谱镜头对该信息的捕捉效果受限;同理,在R6随着叶片衰老和脱落同样导致作物体内叶绿素含量变低,这与陶惠林等[45]的研究结果相似。
不同的建模方法显著影响监测模型的预测精度[46]。本文通过筛选15种植被指数,同时使用3种机器学习方法构建大豆LAI、地上部生物量和产量估测模型。经对比分析发现,基于RF模型的估测精度均高于SVM模型和BPNN模型,说明RF在估算大豆LAI、地上部生物量和产量方面更有优势,这与陈晓凯等[47]的研究结果相类似,主要原因是RF模型对自变量共线性问题敏感度较差,具有强抗干扰性和抗过拟合能力,对噪声异常值容忍性高,且无需检验变量的正态性和独立性,可以有效地提高模型的精度[48];而SVM模型的核心问题是确定核函数与惩罚因子,由于其抗干扰性较差,核函数与惩罚因子等参数选择受限,导致预测效果受限[49];BPNN模型的估计精度较RF模型低,可能是因为本研究中样本数较少导致泛化能力低造成的[50]。因此,RF模型可以作为大豆LAI、地上生物量和产量建模的首选方法。
此外,本研究结果还发现在进行大豆LAI与地上部生物量估算时,在运用同一种建模方法时,前者的估算模型精度均高于后者,这一结果与文献[51-52]的研究结果一致。造成这一现象的原因可能是LAI与生物量相比更为直观,是作物重要的形态参数。目前,基于经验植被指数建立作物生长估算模型和产量估算研究中仍有一些问题有待解决,如不同区域、相同区域不同作物、甚至相同区域相同作物但不同时期的不同植被指数的相关系数变化及模型反演精度可能有所不同。除此之外,本研究后期也可以尝试使用波段更多的高光谱仪进行光谱信息的采集,以避免多光谱受部分光谱信息缺失的不利影响。
(1)光谱指数GCI与大豆LAI、地上部生物量和产量的相关性均较好,通过GCI的变化能实时监测大豆的长势。
(2)基于RF模型构建的大豆LAI和地上部生物量预测模型的精度显著高于SVM与BPNN模型。
(3)在R4生育期,基于RF模型构建的大豆产量预测模型的精度最高,故在R4生育期使用RF模型能最准确地估算产量。