基于无人机多时相植被指数的冬小麦产量估测

程 千 徐洪刚 曹引波 段福义 陈 震

(中国农业科学院农田灌溉研究所，新乡 453002)

0 引言

收获前，及时准确地估测农作物产量有助于相关粮食政策的调整及制定。遥感技术是监测作物生长状况的有效手段，遥感估产则是将遥感信息作为输入变量，通过建立遥感估产模型来估测作物产量。随着科技的发展，无人机能够以低成本和更实用的方式获取高时空分辨率的遥感数据，并将其应用于作物生长监测[1-6]和产量估测中[7-9]，弥补了卫星遥感运行周期长、空间分辨率低和受气象条件影响大等不足。

多光谱传感器包含的近红外波段能够较好地捕捉作物整个生育期内叶片结构变化产生的光谱反射特征，其与其他波段构建的植被指数包含更多的作物生长信息，有助于产量估测。通过建立植被指数与产量间的统计关系进行遥感估产的方法，简单易行，目前已有很多学者利用无人机平台搭载多光谱传感器对不同作物进行产量估测[6,10-13]。

在冬小麦生长阶段，冠层反射率在不断发生变化，冬小麦生殖生长阶段的光谱指数与产量具有更高的相关性，可更好地用来估测产量，如开花期[14]或灌浆初期[15]的光谱植被指数等。随着生长过程的延续，产量与单一时期光谱指数之间的相关性达到非常显著的水平，但是无法在单个时期内确定产量形成的过程。相反，由多个生育期生长信息构成的数据集能更好地反映作物生长变化过程，可以提供更多有用的信息，有利于提高算法的估测精度。文献[16-17]利用多个生育期的植被指数的累积值进行作物产量估测，提高了估测精度。本文尝试将新生育期的植被指数作为样本数据新的属性值添加到上一个生育期的数据集中，形成新的多时相植被指数数据集，从而对冬小麦产量进行估测。

本文采用无人机搭载多光谱传感器，对冬小麦多个生育期进行遥感监测，提取多种光谱植被指数，并利用多时相的植被指数构建数据集，训练机器学习算法，进行冬小麦产量的估测分析。

1 材料与方法

1.1 试验区概况

试验区位于河南省新乡市七里营镇中国农业科学院新乡综合试验基地(113°45′38″N，35°8′10″E)，该试验区地形为黄河北部冲积扇平原，海拔73 m，土壤类型为粘壤土。气候类型属于温带季风气候，四季分明，雨热同期，多年年平均气温14.5℃，年平均降水量560 mm左右，主要集中在6—9月，平均无霜期199 d，试验研究区位置见图1。本次研究共有3种灌溉水量处理(Irrigation treatment，IT)：240 mm(IT1)、190 mm(IT2)和145 mm(IT3)，每种灌水处理包含60个小区；试验选取了30个小麦品种，每个小区1个品种，随机分布种植，则每个灌水处理对应30个品种，重复两次，试验设计如图1所示。

1.2 无人机平台与遥感数据获取

采用大疆经纬M210型四旋翼无人机平台，无人机净质量4.69 kg，最大承载质量1.45 kg,续航时间20 min左右。无人机搭载MicaSense RedEdge-MX型多光谱相机拍摄多光谱影像。多光谱相机质量232 g,图像尺寸1 280像素×960像素,具有蓝(Blue，B，(475±20)nm)、绿(Green，G，(560±20)nm)、红(Red，R，(668±10)nm)、红边(Red-edge，RE，(717±10)nm)和近红外(Near infrared，NIR，(840±40)nm)共5个通道。

冬小麦抽穗之后，进入生殖生长阶段，各器官营养逐步转化为籽粒的形成，小麦叶片开始衰落，此时虽然冠层光谱反射率不断降低，但与产量的相关性逐步升高[12]。因此，从冬小麦开花期到灌浆期，选取多次采集的多光谱影像，进行产量预测。具体采集日期为：开花期(S1，2020年4月14日)、灌浆初期(S2，2020年4月23日)、灌浆中期(S3，2020年4月30日)和灌浆末期(S4，2020年5月10日)。

影像采集时使用大疆GSPro地面站对试验区覆盖范围规划航线自主飞行。为保证影像质量，减少大气影响，选择天气晴好、无云低风速且11:00—14:00时间段飞行。飞行高度30 m，影像重叠率为航向85%和旁向80%,地面空间分辨率为2 cm。

1.3 数据处理

1.3.1影像预处理

采用Pix4Dmapper软件对无人机的多光谱影像进行拼接处理，并导入每次飞行前后拍摄的校准板的影像进行辐射校正，生成该试验区地物反射率的正射影像图。将正射影像图导入开源软件QGIS，进行几何校正、试验小区矢量文件制作、多光谱植被指数计算、各小区多光谱植被指数统计等，得到各小区多光谱植被指数均值作为产量估测的输入变量。

1.3.2多光谱植被指数选取

目前已有多种多光谱植被指数用于作物产量估测研究，如表1所示，大致分为以下几类[18]：①基于红波段和近红外波段运算的多光谱指数：差异植被指数(Difference vegetation index，DVI)、增强植被指数(Enhanced vegetation index，EVI)、比率植被指数(Ratio vegetation index，RVI)、改善简单比率指数(Modified simple ratio，MSR)、归一化植被指数(Normalized difference vegetation index，NDVI)、植被色素比率(Plant pigment ratio，PPR)、结构不敏感色素指数(Structure insensitive pigment index，SIPI)和改进的非线性植被指数(Modified nonlinear vegetation index，MNVI)。②优化土壤背景的红波段和近红外波段运算的多光谱指数：土壤调节植被指数(Soil-adjusted vegetation index，SAVI)和优化土壤调节植被指数(Optimization of soil-adjusted vegetation index，OSAVI)。③利用红边波段特征的多光谱指数：红边叶绿素指数(Red-edge chlorophyll index，CIRE)、红边归一化植被指数(Red-edge NDVI，NDVIRE)、改良叶绿素吸收率指数(Modified chlorophyll absorption ratio index，MCARI)和转化叶绿素吸收反射指数(Transformed chlorophyll absorption reflectance index，TCARI)。④基于从绿波段到近红外波段光谱反射率曲线形状的指数：三角植被指数(Triangular vegetation index，TVI)和改善三角植被指数2(Modified triangular vegetation index 2，MTVI2)等。

表1 植被指数计算公式Tab.1 Vegetation index calculation formula

1.4 数据分析方法

水分亏缺条件下冬小麦生长特征出现差异，遥感影像的像元灰度变化能够反映生长特征的变化。通过多光谱植被指数评价水分亏缺条件下，冬小麦生育期的生长差异，并采用皮尔逊相关系数和决定系数评价不同生育期的多光谱植被指数与冬小麦产量的线性相关程度,选择合适的生育期和植被指数进行产量估测。

分别采用偏最小二乘回归(Partial least squares regression, PLSR)、支持向量机回归(Support vector regression, SVR)、随机森林回归(Random forest regression, RFR)等统计学习算法进行冬小麦产量估测；并采用决定系数R2、均方根误差(Root mean squared error, RMSE)和平均绝对误差(Mean absolute error, MAE)等指标评价不同统计学习算法的估测精度。以上数据分析均使用Python语言中的Scikit-learn[34]机器学习包完成。

2 结果与分析

2.1 植被指数与产量相关性分析

不同生育期多光谱植被指数与冬小麦产量的相关系数如图2所示，可以看出，进入开花期后，光谱指数与产量出现相关性，随着冬小麦生长，相关性逐渐增强，灌浆末期多个植被指数与产量相关系数达到0.7。灌浆期后大部分植被指数与产量呈正相关关系，且植被指数之间相关性较强，而TCARI和SIPI与产量及其他植被指数呈较强的负相关关系。

不同生育期多光谱植被指数与冬小麦产量的线性回归决定系数如表2所示，可以看出，随着冬小麦生殖阶段不断发展，决定系数不断增大，在灌浆末期达到最大。除PPR外，其他植被指数均能在一定程度上反映产量的变化特征。可见本文选取的多光谱植被指数与冬小麦产量有较强的相关性，能够用来估测冬小麦产量。

表2 不同生育期植被指数与产量决定系数Tab.2 Coefficient of determination between vegetation index and yield in different growth stages

2.2 多时相植被指数估测产量

选取不同生育期多光谱植被指数，采用PLSR、SVR和RFR等统计学习算法估测冬小麦产量，结果如图3和表3所示。单时相时(如S1或S2等)，将文中选取的16种植被指数分别作为自变量的1个维度，形成自变量向量，即每个样本数据有16个属性值。当构建多时相植被指数数据集时(如S1+S2等)，将新生育期(S2)的16个植被指数合并到前1个生育期(S1)的自变量向量中，形成新的多时相植被指数数据集，每增加1个时相，自变量向量增加16个维度，即每个样本数据增加16个属性值。

在训练统计学习算法时，为综合提高算法精度和泛化能力，常采用交叉验证的方法，由于试验有3种灌水处理，为减少灌水处理对算法的影响，保证数据分布的均匀性，本文采用10次5折交叉验证，取均值作为每种算法的评测指标值，图3和表3中各评测指标值针对测试集的数据，更能代表模型的估测能力。

由图3和表3可知，当使用单个生育期的植被指数估测冬小麦产量时，由开花期到灌浆初期、灌浆中期和灌浆末期，不同算法的预测精度不断提高，这与植被指数与产量相关性变化规律一致。

表3 不同回归模型估测产量统计分析Tab.3 Validation statistics of different regression models for grain yield prediction

使用生育期(S1+S2)的多时相植被指数比单个时相(S2)的植被指数估测产量的精度高，PLSR模型估测精度R2提高约0.021；SVR模型R2提高约0.015；RFR模型R2提高约0.051。同理，生育期(S1+S2+S3)的多时相植被指数比单个时相(S1或S2或S3)的植被指数估测产量的精度高，说明在灌浆末期(S4)到来前，这种数据集构建方式能够提高产量估测精度，对产量作出较好的提前估测。灌浆末期(S4)的植被指数与产量相关性更好，估测精度更高，此时构建的多时相植被指数，不同模型均有较高的估测精度。

对比3种统计学习算法，针对本文的数据集，当数据属性值增多时，PLSR估测能力低于SVR和RFR。PLSR模型针对生育期(S1+S2+S4)数据测试集的估测精度最高，平均R2为0.459，RMSE为1 822.746 kg/hm2；SVR模型针对生育期S4数据测试集的估测精度最高，平均R2为0.540，RMSE为1 676.520 kg/hm2；RFR模型针对生育期(S3+S4)数据测试集的估测精度最高，平均R2为0.560，RMSE为1 633.896 kg/hm2。而且，RFR模型针对生育期S4后不同时相植被指数组合的数据测试集，均表现出较高的估测精度。

图4为不同灌水处理时，不同品种冬小麦产量估测值与实测值对比，图中估测值针对某一生育期的全部数据样本，且估测模型分别选取估测精度最高的训练模型，PLSR模型针对生育期(S1+S2+S4)，SVR模型针对生育期S4，RFR模型针对生育期(S3+S4)。可以看出，3种模型的估测值均反映不同灌溉处理条件下产量变化规律，并且RFR模型的估测值与实测值更为接近。

正常灌溉处理(IT1)时，冬小麦产量平均达到8 700 kg/hm2，水分亏缺处理(IT2和IT3)时，产量逐步降低，且水分亏缺严重时，产量下降较多，相比正常灌溉，冬小麦减产将近1/2。对比冬小麦产量估测值与实测值，IT1灌水处理时，估测值偏小；IT2和IT3灌水处理时，估测值与实测值较为接近；但每个灌水处理条件下估测值的分布范围均比实测值的小。由于试验区冬小麦品种较多，同一灌水处理条件下，特别是正常灌水处理时，产量相差较大，使得数据分布较为分散。

3 讨论

本文将新生育期的植被指数作为新的属性值，添加到上一个生育期的数据集中，形成新的多时相植被指数数据集，提高了冬小麦产量估测精度，这与文献[11,13,17]的结论一致。而且文献[11]研究中同时使用了株高参数，并发现灌浆期株高比植被指数对冬小麦产量有更高的估测精度，说明不同类型传感器获取的冠层参数，经过数据融合能进一步提高作物产量估测精度。可见，无论是多时相数据或多源数据，均提高了样本数据属性值的维度，有助于产量估测模型的训练。

采用统计学习模型进行回归预测时，样本数据的选择对模型的训练效果有较大影响，特别是样本数据量较小时，训练数据集的划分更为重要。交叉验证法适用于样本数量不大的情况，而且对模型偏差和方差能够兼顾。本文采用P次K折交叉验证法，减小训练数据随机选取对模型预测精度带来的影响。P次为重复K折交叉验证的次数，通过P次重复减小随机误差，本文P选取10。K折即把数据集划分为K等份，其中K-1份作为训练集，剩余1份作为测试集，由于样本总量为180个，单生育期时数据属性为16个，多时相时数据属性最多有64个，为保证不同生育期训练模型的偏差和方差均较低，K值不宜过大或过小，本文K选取5作为数据集划分标准。

冬小麦产量构成因素较多，与多光谱植被指数之间存在较为复杂的非线性关系，PLSR模型虽然包含了主成分分析过程，可以有效解决自变量之间的多个共线性问题，但当自变量与因变量之间存在复杂的关系时，其仍有较大的局限性。SVR模型和RFR模型对复杂的非线性关系均有一定的回归能力，并且RFR模型具有更好的抗噪能力，本文冬小麦试验中有30个小麦品种，数据分布较为离散，针对灌浆末期的不同时相植被指数，RFR模型均表现出较高的估测精度。但是当样本数据量过大，或者样本数据特征过多时，传统统计学习算法的预测能力将下降，需要训练更优秀的模型，如深度神经网络模型等。

4 结论

(1)随着冬小麦的不断生长，植被指数与产量的相关性不断增强，灌浆末期多个植被指数与产量相关系数达到0.7，基于此生育期植被指数估测冬小麦产量的精度也较高。

(2)采用多时相植被指数能进一步提高冬小麦产量的估测精度，采用冬小麦开花期和灌浆初期的多时相植被指数进行估产比采用单个生育期的植被指数估测产量的精度高，PLSR模型估测决定系数R2提高约0.021，SVR模型R2提高约0.015，RFR模型R2提高约0.051。灌浆末期的多时相植被指数，3种模型均有较高的估测精度。PLSR模型估测精度最高时R2、RMSE分别为0.459、1 822.746 kg/hm2，SVR模型估测精度最高时R2、RMSE分别为0.540、1 676.520 kg/hm2，RFR模型估测精度最高时R2、RMSE分别为0.560、1 633.896 kg/hm2，相比之下，RFR模型具有更好的稳定性。

(3)长时间的水分亏缺导致冬小麦产量下降，与正常灌溉水量相比，长时间的水分亏缺将造成冬小麦减产约1/2，此时需及时制定灌溉计划。