基于梯度提升树算法的夏玉米叶面积指数反演

2019-06-04张宏鸣韩文霆刘全中宋荣杰侯贵河

农业机械学报 2019年5期

张宏鸣刘雯韩文霆刘全中宋荣杰侯贵河

(1.西北农林科技大学信息工程学院，陕西杨凌 712100； 2.西北农林科技大学机械与电子工程学院，陕西杨凌 712100)

0 引言

叶面积指数(Leaf area index，LAI)是生态系统能量流动中极其重要的植被特征[1]。在LAI定义的不断发展中，为了既可以适用针叶植被又可以适用阔叶植被，最终被定义为单位地表面积上绿叶表面积总和的一半[2]。传统获取LAI的方法主要依靠人工实地测量[3]，不仅浪费人力，而且对作物破坏性极大。随着遥感技术的不断发展，获取LAI的方法开始转向于遥感反演[4-6]。目前，遥感数据主要来源于卫星和无人机，但卫星遥感属于高空遥感技术，在获取遥感图像过程中易受到大气因素的干扰[7]。且卫星遥感的研究区域面积适合在公顷以上，如果进行小范围的区域研究，其分辨率受到一定的限制。因此，近年来产生了利用无人机搭载多光谱或高光谱相机的方式采集遥感图像，该技术弥补了卫星遥感的不足，提高了图像的地面分辨率，可快速准确获得作物叶面积指数，为农情信息监测提供了有效手段[8-9]。

中国作为农业大国，玉米是不可或缺的农作物。长期以来，我国在作物的叶面积指数反演方面进行了相关研究[10]。早期在进行LAI反演时都选择用单一植被指数作为输入变量[11-12]，但单一的植被指数存在不同程度饱和性[13]，对于LAI的反演受到了一定的限制。之后的研究中发现，红波段和近红外波段的反射率与作物叶片的特征不仅最为密切相关且有别于其他地物[14]，可组合出多种的植被指数，因而是当前用于LAI反演的常用波段[15]。为了可以融合多种与LAI相关性较强的植被指数及其他相关数据，衍生出了机器学习的方法。机器学习就单因子训练模型进行了改进，利用多种相关因子的非线性拟合来构建较为精确的模型，使得验证集与实际值之间的预测误差最小，泛化能力最强[16]。

谢巧云等[17]利用高光谱数据进行作物LAI遥感反演，结果表明，非线性支持向量机模型最适宜用于研究区域冬小麦LAI反演；王丽爱等[18]将支持向量回归和反向传播神经网络算法作为比较模型，证明随机森林模型具有强学习能力和预测能力；张春兰等[19]基于无人高光谱模型研究了4个生育期小麦的LAI，发现随机森林模型LAI的反演精度较高，且适用性较强；HOUBORG等[20]发现，利用机器学习的方法，可有效分析和利用数量大、维度高的观测数据。以上研究成果表明，机器学习中的回归模型广泛应用于作物遥感反演，且取得了较好的研究成果[21-22]。但回归模型都存在小样本数据出现过拟合的问题，且无法判断输入因子中的主要贡献因子。梯度提升树算法(GBDT)可以解决以上问题，该算法通过每一棵树学习之前所有树的残差和结果，一步步迭代构建弱学习器，纠正原模型误差，有效提高了预测精度[23]，可作为反演夏玉米LAI的有效方法。

本文以2018年内蒙古自治区鄂尔多斯市达拉特旗实验基地为研究区域，采集夏玉米无人机多光谱影像，提取关键生育期V3～R5(出苗期至初入蜡熟期)的8种植被指数，实测LAI及玉米株高。将株高和8种植被指数作为输入变量，LAI作为输出变量输入SVM算法、RF算法、GBDT算法3个模型中训练学习，以验证GBDT算法反演模型在夏玉米LAI反演中的适用性。

1 材料与方法

1.1 研究区概况

实验基地(40°25′46.99″ N， 109°36′35.68″ E)位于内蒙古自治区鄂尔多斯市达拉特旗昭君镇，属于典型的温带大陆性气候，干旱多大风，主要粮食作物为小麦和玉米。实验区所种作物为夏玉米，属于春播中晚熟型，一年一熟。播种时间为2018年5月中旬，收割时间为9月。选取株高为60 cm，叶片宽度约两指宽(3 cm左右)时进行田间玉米LAI实验。由于当地的降水量无法满足玉米生长期内的需水量，主要的供水方式为喷灌机灌水。据此，将研究区中的④区利用喷灌机作水分胁迫，使之除大气降水外，额外供水量少于其余区域。研究区无人机影像如图1所示，实验区域面积约为1.13 hm2,以喷灌机为中心，分为5个扇形区域，每个扇形区域内划分3个4 m×4 m的实验样方，即研究区总计15个小样方。

图1 研究区无人机影像和分区示意图Fig.1 UAV image and partition map of study area

1.2 数据的获取与处理

田间夏玉米LAI的测量使用LAI-2200C型植物冠层仪进行。数据于2018年6月26日开始采集，2018年8月25日结束，共9次实验，期间覆盖夏玉米的V3(出苗期)、V6(拔节期)、VT(抽雄期)、R1(吐丝期)、R3(乳熟期)、R5(初入蜡熟期)6个关键生育期。前8次实验数据用于建立模型，最后1次实验数据用于模型预测及验证。LAI采集方式采用ABBBB，即测量时，取4次冠下B值，一次冠上A值。由于玉米属于行栽作物，测量方式通常采用对角线法：B值取样点位于两垄之间均匀分布，A值取样点避免阳光直射，位于冠层上方。每个小样方重复采集4组数据，取均值。考虑正午时分阳光强烈，测量时需配备遮光帽,以减少视野遮盖帽周围的反射光对数值的影响。夏玉米的株高在测量LAI时一同测量同一样方对角线的点，在对角线3点的区域内每一点测量5组数据，取均值。

遥感图像数据采集为六翼无人机RedEdge五波段搭载多光谱相机，相机焦距为5.5 mm，视场角为47.2°，图像分辨率为1 280像素×960像素。相机配备了光强传感器和两个3 m×3 m的灰板。光强传感器可对无人机航拍过程中外界光线的变化对光谱影像造成的影响进行校正，而灰板由于具有固定的反射率，可对航拍影像进行反射率的校正，从而生成反射率影像图，进行植被指数的提取。相机的波段信息以及灰板的反射率如表1所示。实验时晴朗无云，平均气温28.6℃，平均相对湿度61.96%，平均风速1.12 m/s，微风。时间在11：30—14：30。多光谱无人机飞行高度为70 m，飞行方向为南北方向，航向、旁向重叠度分别为80%和70%。将每次实验按照固定航线拍摄的多张图像，以日期为索引导入到瑞士Pix4D公司的Pix4D mapper软件中，以实时动态(Real time kinematic,RTK)测量的方法获取地面像控点，导入小图对应的POS数据，在软件中进行初始化处理，几何校正，构建三维模型，提取纹理以及构造地物特征，最终生成高清正射多光谱影像。由于拼接预处理后的原始图像包含除研究区域以外很大的区域，为了更加突显遥感影像的作物特征，需在ENVI软件中裁剪处理。根据可见光影像中的样方对应裁剪出多光谱影像的15块实测区域，每一块实测区域一一对应实验样地的每一小样方。取裁剪后每一小样方的对角线3点，生成这3点对应的各项植被指数，最终以均值来确定每个小样方的植被指数。

表1 RedEdge多光谱相机参数及灰板对其中心波长的反射率Tab.1 Multispectral camera parameters and reflectivity of gray plate to its center wavelength

由于正午时分和黄昏时分叶片的蜷缩程度受温度影响会有一定的区别，且在不同的时间段接收的光谱信息也有很大的差异，因此多光谱影像数据和地面数据采集需在同一天的同一时间段。

1.3 植被指数的提取

植被对于不同波段入射光子的吸收作用和散射作用不同，形成了特殊的光谱响应特征。由于植被和农作物在红光波段强吸收和在近红外波段强反射的特性，大量的研究证明这两个波段与作物覆盖度和LAI具有很好的相关关系[24]。因此本文借鉴前人研究，选择红光、绿光和近红外波段组合出与LAI相关性较强的8种植被指数进行LAI反演，各计算公式[25-26]见表2。

表2 植被指数计算公式Tab.2 Formulas of vegetation index

注：RNir、RRed、RGreen分别为灰板对RedEdge相机近红外、红波段、绿波段的平均反射率。

1.4 模型构建

1.4.1支持向量机

支持向量机(Support vector machine，SVM)中的支持向量回归(Support vector regression，SVR)是早期机器学习中常用的回归模型。SVR就是寻找一个最优的回归平面，让集合中所有的数据到这个回归平面的距离最近[13]。本文选择了SVM中优化的LIBSVM库进行夏玉米的LAI预测。LIBSVM[27]具有操作简单、快速有效和可处理高维空间数据等特点，常用来做分类和回归。支持向量回归与传统的回归模型相比，优点是此算法可容忍基于模型的输出f(x)与真实的输出y之间有ε的偏差，也就是当|f(x)-y|>ε时才计算损失。但SVM算法对非线性问题没有通用的解决方案，难以找到合适的核函数。

1.4.2随机森林

随机森林(Random forest，RF)是决策树的集成算法，通过对大量分类树的汇总以提高模型的预测精度，是取代神经网络等传统机器学习方法的新算法。随机森林算法中包含多个决策树来降低过拟合风险[28]，思想是Bagging算法和随机特征选取[20]。随机森林通过对训练样本重新采样的方法得到不同的训练样本集；在新的训练样本集上分别进行训练学习，由于每个学习器相互独立，所以此类方法更容易并行；最后合并每一个学习器的结果，从而得到最终的学习结果。随机森林的一个缺点是在噪声较大的分类或者回归问题中容易过拟合。

1.4.3梯度提升树

梯度提升树(Gradient boosting decision tree，GBDT)是集成学习中重要的一种算法，是基于Booting算法的一种改进[29-30]。Booting算法的工作原理是在初始训练集样本上给每个训练样本赋予相同的权值，在每一次训练之后对于出错的样本进行增加错分点的权值，在经过多次迭代后，生成相应的多个基学习器，之后对于这些基学习器进行组合，最终通过加权或投票得到模型。而梯度提升树回归与分类算法的区别是输入的训练数据是残差，即将上一次的预测结果带入残差中求出本轮的训练数据，而不是损失函数的梯度[22]。

梯度提升树具有可灵活处理各种数据、预测准确率高、使用健壮的损失函数和对异常值具有很强的鲁棒性等优点，可有效进行回归预测。GBDT回归算法如下：

(1)输入训练样本

D={(x1,y1),(x2,y2),…,(xm,ym)}

(2)初始化弱学习器

式中L——损失函数

c——样本y的均值

(3)计算负梯度

式中T——最大迭代次数

(4)利用(xi,rti)拟合一棵分类回归树(Classification and regression tree，CART)中的回归树，从而得到第t棵回归树，其对应的叶子节点区域Rtj(j=1,2,…,J)，J为回归树t的叶子节点的个数。

(5)计算最佳拟合值

(6)更新强学习器

(7)得到强学习器表达式

(8)输出为强学习器。

梯度提升树算法中，所产生的树是回归树而不是分类树，GBDT的树会累加之前所有树的结果，这种累加的实现只能用CART回归树实现。

1.5 模型评价

以决定系数(R2)和均方根误差(RMSE)来进行模型精度的评价。

为保证模型的有效性和稳定性，本文基于实验分区来划分不重复的训练集和验证集，将样本3次分组，重复放入模型训练学习。样本组1的训练集为区域①、②、④，验证集为区域③和⑤；样本组2的训练集为区域②、③、⑤，验证集为区域①和④；样本组3的训练集为区域①、③、⑤，验证集为区域②和④。

2 结果与分析

2.1 植被指数相关性分析

基于前人对于作物LAI反演进展的研究，分析光谱特征信息发现作物LAI对于红光参数与近红外参数较为敏感，为后期的植被指数的选择提供了方向。本文对3组样本中8种植被指数和株高与LAI进行相关性分析，结果如表3所示， LAI与8种植被指数和株高在P<0.01水平呈极显著相关，训练集相关系数均不小于0.660，验证集相关系数均不小于0.668。NDVI、OSAVI、RDVI、RVI、SAVI、EVI2、MASVI、TVI与LAI在总样本中相关系数平均为0.777、0.765、0.751、0.769、0.747、0.743、0.744、0.715(P<0.01)，因此可选择此8种植被指数作为构建LAI反演模型的变量。而本实验中新加入的实测株高在P<0.01水平下与LAI的相关系数均值也达到了0.769，说明株高与夏玉米的叶面积指数有着较强的相关性，可以选择与8种植被指数一同作为SVM、RF和GBDT模型的输入变量，进行夏玉米LAI的预测研究。

2.2 夏玉米LAI反演模型的构建

由于实地采集的样本相对较少，因此选取3个区夏玉米生长关键生育期的9次实验数据进行反演模型的训练，8种植被指数和同期的株高共9个因子一同作为训练模型的输入变量，夏玉米LAI作为输出变量，分别使用SVM算法、RF算法和GBDT算法来构建夏玉米LAI反演模型。SVM算法模型用LIBSVM库来实现，核函数选择径向基，其余参数根据网格搜索法来确定最优参数；RF算法模型根据多次实验，确定决策树的数量为100，节点分割变量为3；而本文构建GBDT算法模型，用Python语言编写回归程序，根据输入样本组的不同，防止出现过拟合现象，需多次实验确定每组训练集的子模型数目(n_estimators)、损失函数(loss)、树的最大深度(max_depth)等参数。

表3 植被指数与LAI相关性分析结果Tab.3 Correlation analysis result of LAI and vegetation index

注：** 表示在P<0.01水平上极显著相关。

以样本组2为例说明参数选择过程：调参前，样本组2由于选择的树深和子模型数目等过大引起过拟合(图2a)。调参时，选择最小绝对偏差(lad)为损失函数，n_estimators为500，max_depth为4，步长(learning_rate)为0.01，叶节点最小样本(min_samples_leaf)为6，作为样本组2的模型参数(图2b)。由图2可以看出，在提升迭代次数(Boosting iterations)为500时，训练集偏差和验证集偏差(Deviance)分别最小，500为最优的提升迭代次数，过大或过小都会造成预测精度的降低。

图2 样本组2调参对比Fig.2 Contradistinction of parameter adjustment for sample group 2

GBDT算法模型训练结束后，会出现所使用输入变量的相对重要度(Relative importance)，便于理解哪些因素对于预测结果有关键影响力。同样也可以判别出几种由红、近红外波段组合出的与LAI强相关的相似植被指数中，对于夏玉米LAI反演结果的影响力强弱。由图3可知，3组样本中株高在相对重要度中比例占据第一，是影响模型精度的主要因素；而MASVI在3组样本中相对重要度比例相对比较低，即MASVI对LAI反演模型的影响力较弱；综合3组的变量重要性分析，株高对于LAI反演模型的结果贡献较大，MASVI对于LAI反演模型的结果贡献较小。

图4 夏玉米LAI实测值与预测值关系Fig.4 Relationship graphs of measured and predicted summer maize LAI

2.3 夏玉米LAI反演模型比较分析

先将3组训练集分别以3种算法进行训练，不断进行调参和多次迭代训练，得到优化的模型；进而分别将3组样本组中独立于训练集的验证集，作为训练模型的最后验证；最后将训练集和验证集模型预测得到的LAI与实地测量的LAI分别进行散点图分析，拟合成线性的回归线(图4)。

由图4可知，GBDT算法对于连续值的预测效果较SVM算法和RF算法好。GBDT算法在每一样本组中都体现出强大的学习能力，训练集的决定系数R2分别为0.815 4、0.746 5、0.847 5，对应的RMSE为0.000 7、0.000 8、0.000 6；验证集的决定系数R2分别是0.571 0、0.755 8、0.644 1，对应的RMSE为0.002 7、0.001 5、0.001 6。

2.4 夏玉米LAI空间分布

依据GBDT算法对3个样本组反演结果的分析，选取效果最佳的样本组2模型作为反演模型，对8月23日的数据进行研究区内夏玉米LAI反演，该研究区的夏玉米LAI空间分布如图5所示。

图5 研究区夏玉米LAI空间分布图Fig.5 Spatial distribution map of summer maize in study area

由图5可知，8月23日的夏玉米LAI在①、②、③、⑤区主要集中在1.8左右，④区的LAI在1.4左右。与8月23日的实测LAI相比较，总研究区域实测均值为1.851 4，实测最大值为2.190 0，实测最小值为1.315 8。5个区域实测均值分别为2.092 6、1.795 8、2.087 0、1.406 9、1.888 1。由以上数据可知LAI反演数据与实测数据基本相符。整体上①、②、③、⑤区LAI较高，这是因为①、②、③、⑤区喷灌正常，玉米长势均匀；研究区④夏玉米LAI较低，主要是因为④区进行了一定的水分胁迫，导致夏玉米生长较其余4个区缓慢，区域内作物稀疏；综上，基于GBDT算法模型反演的LAI与研究区的LAI空间分布相对一致，进一步证明了该模型的合理性和可靠性。

3 讨论

LAI反映的是叶片的疏密程度，随着作物的生长，高度不断增长，叶片的疏密程度按照稀疏-稠密-稀疏变化，高度表现出对于LAI的显著影响。本文经过输入变量相关性分析和输出结果相对重要度分析两方面佐证了株高对LAI具有极显著关系，将遥感数据(植被指数)与实测数据(株高)相结合一同进行LAI的反演，可得到较好的反演效果。因此，LAI反演中不应只局限于遥感图像提取的植被指数，还可考虑其他对于LAI有影响的因子。

机器学习算法与作物遥感反演密切相关，一个回归能力很强的机器学习算法模型，可以融合众多因子(多种植被指数和株高)共同反演LAI，大幅度提高结果的精度。本文将GBDT算法应用到农作物邻域，反演夏玉米LAI。基于模型稳定性的考虑，分为3组样本重复训练，以8种植被指数和株高构建GBDT算法反演模型，并与SVM算法和RF算法R2、RMSE相比较，结果表明GBDT算法构建的模型两项指标均高于SVM算法的模型，同样也较同出一派的RF算法构建的模型性能有了进一步的提升。在后期的研究中可以将GBDT算法应用到玉米叶绿素、玉米生物量等作物相关参数的反演，以扩大精准农业技术支持的范围。

GBDT算法的优势在于将若干个弱学习器组合成强学习器，结果是多棵回归树的累加之和。由此算法构建的模型可以灵活处理各种数据，不论是本文中的LAI连续值，还是后续研究中作物冠层温度的离散值；在小样本的回归问题中，GBDT算法可通过设置不同的损失函数以及在相对较少的调参时间下，提高反演精度。GBDT算法不仅减少了SVR模型因选择核函数和其他参数造成的时间复杂度的浪费，而且也解决了RF算法对待各输入因子是同一权值，无法判断其中每一因子的贡献率的问题。因此，GBDT算法在回归问题中有很强的应用价值。但还需要注意两点：①利用GBDT算法进一步判断了在红波段和近红外波段组成的几种相似植被指数的不同影响度，突出了算法的优势，但在精度方面的后续研究中还有进一步提升空间。②GBDT算法中的基学习器之间存在依赖关系，一般难以进行并行计算。本文尚未考虑各个基学习器之间的并行操作，在今后的研究中应着重考虑如何实现部分的并行操作，进一步提高反演模型的效率。