基于无人机多光谱数据的玉米苗株估算

赵晓伟， 黄杨， 汪永强， 储鼎

(1.黑龙江省测绘科学研究所，哈尔滨 150081； 2.中国科学院东北地理与农业生态研究所，长春 130102)

0 引言

玉米产量取决于多种因素，包括遗传因素、生长环境、播种方式等。玉米的出苗情况同样是影响玉米产量的关键因素，对出苗情况的研究能够提早预防因土地缺苗而造成的产量损失[1]。出苗情况的计算主要依赖于苗株数的估算，出苗情况的高低是作物高产、稳产的前提条件，缺苗直接影响作物的品质产量，所以对作物的长势评估可以影响农田施肥和管理[2]。因为无人机(unmanned aerial vehicles，UAV)遥感技术的快速发展，UAV技术在农业科学中的应用越来越广泛，提取农田信息的方式愈加成熟、多样、灵活[3]。UAV技术不仅可以节约人工实地检查出苗情况的人力物力成本，而且减轻了人工对作物的损害，也可以为大范围内的出苗信息获取提供优质服务[4]。

目前，国内外关于UAV影像在植被识别和分类方面开展大量研究。Mitch等[5]使用UAV可见光影像，基于植被颜色指数及纹理合理地分类了植被类型，证明了颜色指数对提取植被信息的可行性； 丁雷龙等[6]使用归一化绿-红差值指数(normalized green-red difference index，NGRDI)、过绿指数(excess green index，ExG)、过绿减过红指数(excess green minus excess red index，ExG-ExR)和绿叶指数(green leaf index，GLI)等4种RGB颜色指数，识别了植被区域和非植被区，并对几种颜色指数做出了适用性分析和因素干扰性评价，得出4 种颜色指数识别植被的精度均超过90%的结论； 汪小钦等[7]基于可见光波段UAV影像，提出了可见光波段差异植被指数对植被信息提取的精度可达 90%以上。但是对于哪种颜色指数更适用于分割类似于玉米农作物和土壤背景的相关研究仍很薄弱，依然没有统一的定论。

随着深度学习的发展，遥感技术与机器学习方法相结合的方式，能够实现农作物长势的快速高效评估，逐渐替代了原始人工评估出苗情况的笨重方法，已经在作物识别和长势监测中展现优势[8]。Zheng等[9]综合了深度学习算法和UAV多光谱遥感手段，分别构建了支持向量机(support vector machine,SVM)和支持向量回归(support vector vegression,SVR)的棉苗株数估算模型，对比了2种模型的精度，为估测棉苗株数提供了有价值的参考; 戴建国等[10]同样利用UAV可见光遥感影像和形态学特征参数，构建了SVR模型，提取出棉花苗情信息。但是二者研究的作物都是棉花，棉花和玉米无论从形态还是耕种地点、方式、季节均差异较大，SVR模型是否仍然适用于估算玉米株数，值得商榷。刘帅兵等[11]虽然基于UAV数据通过颜色变换将玉米与土壤分离，通过Harris,Moravec和Fast角点检测算法有效地识别出了玉米株数，但是未结合机器学习的算法，且仅仅适用于实验小区内，是否适用于东北大面积的玉米耕地仍然未知。

针对上述不足，本文基于UAV多光谱数据，计算不同颜色指数，利用OTSU阈值分割算法，找到最佳颜色指数提取苗对象。组合不同数量的形态学特征参数，找到最佳的特征参数组合。结合SVR模型，预测出玉米苗数，统计出苗情况。在东北大面积种植玉米的条件下，展示了UAV多光谱数据对玉米苗期评估长势的适用性，为提高田间尺度的植被精细化管理提供了保障。

1 研究区概况及数据源

1.1 研究区概况

研究区位于吉林省长春市开安镇，地理坐标为E125°5′～125°8′，N44°30′～44°50′，地处松辽平原腹地，属于大陆性季风气候区，主要粮食作物为玉米，属于农业发达地区。该地区的农业规模化和机械化程度较高，是开展智慧农业、UAV遥感技术的良好基地，能够满足农业及遥感领域科学研究的需求。该区域采用机械播种方式，播种间隔为0.2 m，垄距约0.6 m，图1为在研究区裁剪后的多光谱影像，红色矩形框为人工测量苗株的样方。

图1 研究区及监测样方

1.2 UAV数据获取与预处理

UAV多光谱影像于2020年5月28日拍摄，拍摄时天气情况良好，无风少云。UAV平台为大疆精灵4，水平飞行速度为3 m/s，飞行高度为30 m，航向重叠为60%，旁向重叠为60%，本次试验共采集260幅玉米苗多光谱影像，数据为 24位真彩色jpeg格式。使用PhotoScan软件对影像进行拼接，然后将影像进行裁剪，去除边缘空洞部分和异常值，得到如图1所示的影像数据，拼接后的UAV影像长约150 m，宽约120 m，空间分辨率为0.5 cm。

在研究区均匀选取28个样方，如图1所示的红色方形框，用来调查样方中玉米苗的株数。监测样方为3 m×3 m的方形，方向分别为顺垄和垂直垄方向。研究区的播种垄宽为0.6 m，所以尽量保证垂直垄方向能够包含5根垄，顺垄方向的株距为0.2 m，所以理论上每个样方的苗数为75株。因整幅UAV影像覆盖面积较大，难于处理，本研究按照不同尺度的监测单元对整幅影像做分割处理，分割成渔网状，分别对每个监测网格内的苗株数做预测。监测单元边长为0.5 m,1 m,2 m，3 m，4 m，5 m，6 m，7 m的正方形。3 m×3 m监测单元内的株树为实地量测样方，因UAV影像的分辨率可达0.5 cm，所以在监测样方尺度下，玉米苗清晰可见，通过目视判读的方法进行估算，同时用实测样方内株数对目视判读的株数做验证，其他尺度内的真实株数采用目视解译的方法进行估算。

2 研究方法

2.1 玉米苗对象提取

如何将玉米苗从与土壤背景中提取出来是获取玉米苗对象的关键。本研究的UAV影像数据主要由玉米苗和土壤构成，褐色的为土壤，绿色的为玉米苗。因此可以通过线性组合红(R)、绿(G)和蓝(B)3 种颜色成分使绿色植被与褐色土壤的差异最大化，有效分离玉米苗与土壤背景[12-13]。此外，合理的选择颜色指数对UAV影像的精准分割十分重要[14]。本研究初步选取的颜色指数为ExG、绿-蓝差值指数(green-blue difference index，GBDI)、ExG-ExR，NGRDI和GLI[15-16]等5个指数作为候选颜色指数。通过特征分析和SVR模型精度的结果对比得出最优颜色指数。首先将UAV影像数据转成双精度(double)类型颜色参数，分别将R，G，B像素值除以255做归一化处理，得到取值范围为0～1之间的r，g，b值。然后计算5个颜色指数，计算公式分别为：

ExG=2g-r-b，

(1)

GBNI=g-b，

(2)

ExG-ExR=3g-2.4r-b，

(3)

NGRDI=(g-b)/(g+b)，

(4)

GLI=(2g-r-b)/(2g+r+b)。

(5)

通过5种颜色指数实现阈值分割玉米苗与土壤背景。阈值的精准确定是从土壤背景中分离出玉米苗对象的关键因素，最大类间方差法(OTSU)作为阈值分割方法，比其他图像分割方法更加直观、简单[17-18]，本研究选取OTSU算法可以到达自动确定阈值的目的。OTSU[19]算法是一种自适应的阈值提取方法，利用的图像灰度特性能够找到目标对象和背景类间方差S的最大值，最佳阈值即为S最大时对应的值。具体算法如下： 首先，分别得到玉米苗和土壤背景的像元比例(阈值为T)w0和w1、平均灰度g0和g1; 其次，计算图像总平均灰度gm，计算二者的方差S，公式分别为：

gm=w0g0+w1g1，

(6)

S=w0(g0-gm)+w1(g1-gm)

；

(7)

最后，求得使S最大的阈值T，即为最佳阈值。

2.2 模型构建

以图2为例，该样方点内共有37个苗对象。红色矩形框为各个苗对象的外接矩形，可计算出各个苗对象外接矩形的长、宽、周长、面积等特征参数。

2.3 去除杂草干扰

观察UAV影像数据发现，研究区内的垄间杂草容易被划分为连通区域。为去除杂草，本研究首先将获取的玉米影像旋转为田垄竖直的方向，结合特征参数提取苗对象外接矩形的中心点坐标，导入ArcMap软件中，给定投影，生成矢量点数据。对点数据做缓冲区分析，因播种株距为0.2 m，所以缓冲距离设置为0.2 m。生成缓冲区后，顺垄方向连接相邻点，生成线矢量文件，将其当作垄的中心线。判断各点距离最近中心线的距离，如果距离大于0.2 m，被认为为杂草，该方法能有效去除垄间面积较大的杂草，如图3所示，红色矩形框内为杂草。

2.4 基于 SVM 建立苗数统计模型

在计算出特征参数后，下一步是将特征参数和人工实际测量苗数做拟合，但是如何选择特征参数以及特征参数数量是提高拟合精度的关键。本研究采用深度学习算法中的SVM试图解决该问题。

SVM[20-21]被广泛应用于机器学习、人工智能、大数据等领域，其本身是解决二分类问题，而SVR是SVM的重要应用。SVR追求的最优超平面能够让所有样本点和超平面(函数)的“距离”最小。而回归的本质是找到数据的内在关系，对于SVR而言，不管数据有多少种类别，SVR都能求出超平面，拟合数据，建立模型，给定输入参数，就能够求得一个新值[22]。因此，本研究计算出各个样方点内的所有苗对象特征参数，结合SVR来估算类似于图2(右)所标记的所有苗对象包含的植株数。将二值图像中所有苗对象(单株和连通区域)的特征参数作为SVR模型的分类变量。并且根据原始实地人工测量苗数，统计出相应苗对象的植株数量作为模型识别结果，统计共计1 003条样本数据。70%样本点用于模型训练，30%用于模型测试。

训练过程中采用对数据进行归一化处理，采用默认的RBF核函数。首先，利用交叉验证方法寻找惩罚因子c，c的大小决定模型是否过拟合或欠拟合； 其次，寻找RBF核函数中的方差s，s的大小决定支持向量的多少; 再用最佳的参数训练模型; 最后，对数据做反归一化处理。

3 结果分析

3.1 玉米苗对象提取结果分析

在本研究中，研究区选取28幅样方地点图像，每个样方为9 m2，水平方向为平行垄，长度为3 m，垂直方向为垂直垄，长度为3 m。在ArcMap软件中打开带投影的UAV数据，找到对应的样方坐标位置，在影像数据上量测对应的样方面积，截取图像中玉米苗植株、土壤的样本对象。然后进行颜色指数计算，采用OTSU阈值分割方法，分别计算5种颜色指数下，每个指数分割玉米苗与土壤背景的最佳阈值，得到5种颜色指数下玉米苗与土壤背景的二值化影像数据。以监测样方23为例，图4为该样方不同指数下的数据结果对比。

(a) 原始影像(b) ExG二值化影像(c) NGRDI二值化影像

(d) GBDI二值化影像(e) GLI二值化影像(f) ExG-ExR二值化影像

各个指数的数据影像通过腐蚀、膨胀及形态学的开运算、闭运算等操作步骤，并剔除杂草噪声后效果如图5所示。

(a) 原始二值化(b) 形态学操作(c) 缓冲区方法去杂草操作

如图5(a)所示，红色矩形框内为玉米苗植株周边的噪声，面积较小，基于形态学操作可以剔除。但是对于面积较大的杂草，如图5(b)红色圆形框内的杂草，该操作去除杂草的效果不佳，将杂草作为玉米苗的错误识别会降低出苗情况计算的精度。而本文去除杂草噪声的效果如图5(c)所示，可以看出该方法能够有效剔除垄沟周边较大面积的杂草，提高识别玉米苗对象的精度。

3.2 基于SVR的最优模型

本研究基于5种指数，随机选取所有样方点70%的数据作为训练集，30%数据作为测试集。分别对各个指数，所有特征参数进行SVR模型测试。如图6所示为基于ExG指数估算株数的训练集和测试集预测结果，综合对比各个指数SVR模型的预测精度(表1)。

(a) 训练集(b) 测试集

从图6可以看出，基于ExG指数估算株数无论在训练集和测试集，样本点主要集中在1株，侧面证明了本实验田的玉米苗以单株植被为主的特点，且整体偏差较小。样本点随着株数的增大而减少，说明本实验田玉米苗的连通区域多以2株为主。综合表1来看，SVR模型估算值对真实值的模型拟合精度都较高，R2均在0.81以上。总体相对误差较小，RMSE均在0.02以下。表1可以看出，基于ExG指数的SVR模型精度最高，本研究采用该指数进行下一步研究。

表1 各个指数训练和测试精度

3.3 最优特征参数组合

由3.2节的结果可知，在所有特征参数参与计算的前提下，基于ExG指数的SVR模型精度最高，但R2都没有超过0.9，精度仍然没有达到预期的效果。而正确选择特征参数会影响模型精度[21]，因此评估各特征参数的重要性，优化选择特征集，能够达到简化模型、降低SVR的过拟合风险，提高模型精度的目的。本研究将选用相比其他指数精度较高的ExG指数，对所有样方内苗对象的特征参数和玉米苗实测株数做Pearson系数[4]的相关性分析，将Pearson系数大于0.7的特征参数(图7)，作为输入模型的候选参数。逐渐增加特征参数的个数，采用数学理论上的组合方式，选择最优的特征参数以及特征参数的数量，最终作为输入变量进行建模。

图7 特征参数实测株数的Pearson系数

如图7所示，本研究最优参数共有7个，分别是面积A、周长B、矩形长D、矩形面积F、矩形周长G、椭圆长轴长度H、形状因子Q，共有127个组合方式，将不同的组合方式分别作为基于ExG指数下SVR模型的输入参数。计算玉米苗训练集和测试集实测和预测苗数的拟合精度如表2所示。

表2 特征参数组合

从表2中可以看出，经过选择后的特征参数组合精度明显有所改善，其中，最优组合基本都包括多边形面积A参数。这也侧面印证了该参数对模型的重要性很大，对精度评价的所占比重较大。而且在选择特征参数时并不是参数选择得越多，模型精度越高。针对本研究最佳的特征参数组合为6种特征参数： 面积A、周长B、矩形长D、矩形周长G、椭圆长轴长度H、形状因子Q。过多或过少的参数都会使SVR模型过拟合或欠拟合，只有适当的组合特征参数才能使SVR模型精度最高。

3.4 模型评价与出苗情况计算

本研究基于ExG指数，选取面积A、周长B、矩形长D、矩形周长G、椭圆长轴长度H、形状因子Q等参数的组合方式构建SVR模型进行测试。将估算株数和实测株数的测试和训练集做散点图分析，结果如图8所示，准确率为达到96.54%。整体上来看，SVR模型估算值对真实值的模型拟合精度较高，总体相对误差更小。模型在训练集和测试集中的R2相差不大，体现了基于UAV数据的SVM算法在玉米苗计算株数方面存在优势。

图8 实测株数和预测株数的比较分析

从图8中可以看出，当玉米苗为单株植被时，更容易发生过高估计。该误差可能是由玉米叶的不均匀性且玉米叶较大而导致冠层重叠，进行引起形态学特征误差，容易将单株植被当作连通区域。当玉米苗为2株以上时，总体上呈过低估计，该误差可能是因为连通区域各个玉米苗距离较近，所以各玉米苗的冠层重叠度高，容易被错误分割成数量较少的苗对象。

进一步分析模型识别的预测株数和实测株数，以确定SVR模型的适用性。按照1～7的顺序，逐个统计出株数中被高估和低估的数量，统计结果如表3所示。从表3中可以看出，玉米苗数量为单株时，误差最小，模型精度最高。随着株数逐渐增大，误差变大，精度变差。所以，该模型更适用于玉米苗长势初期，植被冠层较小，连通区域较小的情况下，更能体现本模型的优越性。在1～3株的情况下，总体而言，误差较小，估测苗数的准确性较高。所以只要在合理的UAV数据获取苗对象数据的时间内，比如在玉米苗生长初期，幼苗多以单株为主，冠层覆盖度不大，连通区域面积不大的情况下，模型具备一定的适用性，所以科学合理的UAV拍摄时间是十分重要的，同时也说明了模型在实际应用中有植被生长窗口期的局限性。

表3 实测株数和预测株数的过高过低估计

3.5 监测单元尺度与玉米株数制图

本研究对整幅UAV影像做分割处理操作，但是分割尺度的大小是否对精度有影响，研究做了进一步分析，图9为不同尺度下监测单元的玉米苗对象识别效果。如图9所示，0.5～4 m的分割尺度下，可以有效识别出玉米苗对象，且识别效果差距不明显，如图9中红色椭圆形所示，识别的玉米苗对象随着尺度变大，形状并未改变，只是在相同的屏幕内形态由大变小。因此，在0.5～4 m的尺度下，基于SVR模型对苗株数估算的精度未受到影响，但是随着监测单元尺度增大，识别效果逐渐变差。其中，如图9中5～9 m尺度下的蓝色位置所示，当监测单元尺度大于5 m后，一些长势较差的幼苗无法识别，这势必影响株数估算精度。综上，考虑到实测样方尺度为3 m，以及运算效率，本研究监测像元的尺度设置为3 m更为合理。同时，在研究区内随机选取一个顺垄方向20 m，垂直垄方向24 m的矩形地块，将该块地按照1 m×1.2 m的分辨率分割成多个小格子。基于SVR模型识别结果来计算该地块的玉米苗株数，该区域如图10(a)所示，基于SVR模型识别计算的玉米苗株数如图10(b)所示。图10(b)可直观地显示该区域玉米苗的出苗状况，且小格子颜色代表玉米苗株数的高低，颜色越红，代表出玉米苗株数越多。从图10(b)可知，多数小格子填满了该地块的分辨率，玉米总体的玉米苗株数多集中在7株数以上、长势较好。图10(a)可以看出，有一块无苗断垄的空缺地，对比图10(b)，同区域内的小格子为深蓝色，玉米苗株数为0，这可能是人为漏播所导致的结果。

图9 不同尺度监测单元的苗对象识别效果

(a) 随机矩形地块(b) 矩形地块玉米株数

4 结论

本文基于UAV多光谱数据，通过对比不同颜色指数，借助OTSU阈值分割算法获取玉米苗对象，确定最佳颜色指数的二值化苗对象。优化组合出形态学特征参数，确定了最佳特征参数组合。借助支持向量机回归(SVR)模型，预测出玉米株数，统计了出苗情况。研究成果能够在短时间内迅速、快捷、准确地得知玉米出苗情况和长势趋势。

ExG，GBDI，ExG-ExR，NGRDI和GLI等颜色指数均能用于UAV影像识别玉米，但是针对本研究区，基于ExG颜色指数的UAV影像识别玉米株数的精度最高。不同的特征参数组合会直接影响SVR模型的拟合精度，使SVR模型过拟合或欠拟合。针对本研究的最佳特征参数组合为6种特征参数： 面积A、周长B、矩形长D、矩形周长G、椭圆长轴长度H、形状因子Q。该组合特征参数的SVR模型精度最高。

利用形态学特征参数构建的SVR模型能够有效解决玉米植株粘连在一起时统计苗株困难的问题。模型分类精度达到96.54%，统计误差为0.6%。当玉米苗为单株植被时，更容易发生过高估计。而当玉米苗为2株以上时，总体上呈过低估计。在1～3株的情况下，误差较小，估测苗数的准确性较高。随着苗株数逐渐增大，模型精度逐渐变差。当幼苗多以单株为主，冠层覆盖度不大的情况下，模型的准确性最高。

本研究对杂草的处理仍然不够完善，对于垄间杂草处理效果较好，但对于一些垄上小面积杂草未能进行有效的处理，这对模型精度的干扰性未知。同样，实测样方需要人力和物力，受制于一些客观条件的限制，本研究的样本点为1 003个，样本点较少，如果增加样本点的数量，势必会提高SVR模型的预测精度。