基于细粒度校正的育种小区小麦株高无人机测量方法

吴婷婷 刘昕哲 聂睿琪 刘 佳 武 璐 李 涛

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

0 引言

株高作为小麦的重要表型性状[1],是生物量估算模型的重要输入参数[2-4],还可用于小麦产量预测[5]、氮素营养评价[6]、倒伏识别[7]等,对选育高产、稳产且具有抗逆性的小麦品种具有重要意义[8-9]。现有田间小麦株高主要依赖人工标尺测量,效率低下,难以适应高通量表型组学的发展需要[10]。为解决此问题,可以通过效率高、无主观误差的间接法测量小麦田间株高性状。

现有的间接测量方式有利用雷达点云、深度相机、无人机(Unmanned aerial vehicle,UAV)遥感等手段获取田间育种小麦株高[11-13]。苗艳龙等[14]提出一种基于地基雷达生成三维点云的株高测量方法,与人工测量结果对比,均方根误差为21.5 mm。孙翔[12]提出一种基于标尺图像识别的株高测量方法,相对误差小于3.55%。此外,无人机遥感因成本低、灵活性高等优势,在大田作物表型中得到广泛研究。樊意广等[15]基于UAV生成试验区域的数字正射影像(Digital orthophoto map,DOM)和数字表面模型(Digital surface model,DSM)提取冠层光谱特征和株高,实现了植株氮含量的估算。文献[16-17]对基于UAV的多种空间辅助信息组合进行作物株高估算,结果表明最佳株高估算结果中决定系数R2为0.736,RMSE为5.3 cm。OEHME等[18]基于无人机DSM信息对不同生育阶段玉米株高进行估测,测量误差范围为12.07%～19.62%;JAMIL等[9]发现基于UAV的小麦株高估测和地面实测数据相关系数仅为0.49,RMSE高达22.04 cm;HAN等[19]提出一种利用多层地面控制点改进无人机作物株高估算的方法,RMSE降低约20%。

低成本无人机测量作物株高普遍精度不够理想,特别是对小麦这样的低矮作物[20],在进入特定生育期后株高增长缓慢,高度落差不明显,直接采用上述方法会存在较大误差。作者通过全生育期高频次飞行反演育种小区粒度(面积通常是几平方米)的株高性状,以平衡单次测量误差,然而,多次反演出株高数据逆增长现象,显然不符合作物生长机理,且发现这种逆生长差距随时序增长被进一步放大。若随意抛弃逆增长的株高反演结果,则生育期数据空缺易造成整体模型的不稳定,为此,在缺失地面控制点高程信息的条件下,考虑将时序变化中遇到的实际问题转换到空间变化中来分析,本文提出一种在育种小区尺度下基于空间已有株高真实数据的近邻校正法以及一种基于多光谱+RGB数据融合的半经验光谱指数校正法,以实现对田间育种小区小麦株高的无人机遥感精确测量。

1 试验材料与数据采集

1.1 试验设计

试验研究区位于陕西省杨凌农业高新技术产业示范区实验基地(34°17′27.05″N,108°04′31.59″E),海拔435～563 m、年降水量635.1～663.9 mm,年平均气温12.9℃。选取565份(含5个对照品种)具有丰富遗传多样性的小麦种质作为本次研究的试验材料,包括当前国内的主栽品种、高代品系、核心种质、农家种、骨干亲本、不同时期培育的代表性品种以及国外引进的关键材料等。大田田间试验设计采用随机增广试验,共种植640个样本小区,分16个区组,每区组选取西农511、百农207、偃展4110、济麦22和周麦185等5个品种作为对照品种,每个小区播种6行,行长5 m,宽1.2 m,播种密度每平方米270株基本苗,播种期为2021年10月18日,具体研究区规划如图1所示。

图1 研究区规划图及UAV平台Fig.1 Research area map and UAV platform

1.2 无人机数据采集

无人机数据获取平台如图1所示,RGB图像获取平台采用小型多旋翼高精度Phantom 4 RTK型航测无人机(DJI大疆创新公司,深圳),搭载RGB相机对大田农情信息进行采集,相机分辨率2 000万像素,搭配8.8 mm定焦镜头(等效焦距24 mm),多光谱图像获取平台采用大疆Inspire2型无人机搭载MicaSense RedEdge-MX型多光谱相机对农田进行采集,多光谱相机每次拍摄的原始图像尺寸为1 280像素×960像素,相机焦距为5.5 mm,视场角为47.2°,相机配有标准反射板,每次采集数据需要对相机进行校正,多光谱相机各通道中心波长与半波宽如表1所示。

表1 RedEdge-MX型多光谱相机的通道参数Tab.1 Channel of RedEdge-MX multispectral camera nm

2021年10月31日至2022年5月31日,间隔采集试验田的RGB与多光谱影像共31次。为保证数据质量,选取晴天拍摄,拍摄时间为11:00—15:00之间,设定无人机自动巡航的飞行区域为150 m×100 m,飞行高度25 m,飞行速度2 m/s,航向重叠率85%,旁向重叠率80%,相机垂直拍摄,拍照间隔时间1 s。

1.3 地面数据采集

为了获取作物高度,以评估株高估算精度和对应无人机图像质量,在冬小麦的拔节到开花期选取试验田中的107个小区进行了6次完备株高标尺测量,采样日期为:2022年3月7、11、28日和2022年4月5、13、27日;采样时间为10:00—11:00。测量土壤至小麦冠层的高度(不包括麦芒部分)。

2 基于无人机遥感的小麦株高精确测量

2.1 数据预处理

分别将多光谱图像、RGB图像导入到Pix4Dmapper软件中,利用多光谱相机厂家提供的标定系数依次对小麦多光谱单幅影像进行标定,依据对应的位置与姿态数据生成稀疏点云,采用多视点立体匹配算法生成密集点云数据,最后导出为TIFF格式的DOM和DSM数据。生成的多光谱DOM影像空间分辨率约为每像素0.8 cm、可见光DOM影像空间分辨率约为每像素0.56 cm,其投影坐标系设置为WGS 1984坐标系下墨卡托投影分度带的北半区49°带。由于携带多光谱相机的无人机不具有RTK,故需选择标准日期,利用地面控制点(GCP)将其余日期以标准日期为参考进行校准。为避免田间试验破坏性采样和杂草对试验结果精确度造成影响,剔除田间破坏性采样区域保留感兴趣区域,将DOM及DSM(DTM)分割为单独小区,同时利用阈值分割法去除多光谱图像中的杂草信息。

2.2 高频次UAV株高测量的误差机理分析

采用UAV测量小麦株高时,通常将小麦的自然株高(小麦最顶端不包括麦芒部分到土壤的距离)作为度量株高变化的评价指标[16];依据研究侧重点选取日长、积温、生长度日等参数作为时序指标[21-23]。为更好观测细粒度尺度下株高变化,本文选取日历日期作为衡量小麦时序变化的指标。

表2为利用传统方法[24](Traditional method,TM)估测的小麦株高和地面真值,在2022年3月11日、4月27日误差较大,这是由于无人机飞行高度定位误差、小麦冠层高程提取误差而导致。如图2所示,随着时间的增长,小麦冠层发育逐渐封垄使得拍摄的RGB图像特征点变少,同时由于进入多风多雨季节,采样时由于风的影响使得冠层产生扰动,从而导致图像在拼接生成DSM时误差被进一步放大。此外由于无人机处于25 m的较高高度进行群体测量,麦穗间间隙导致小区冠层顶面非水平面,以及部分高杆小麦自然下垂但测量人员会将其扶正后测量,均会导致测量结果比真值偏小,且小麦封垄将进一步导致后期株高监测更为复杂,测量误差也更大。

表2 基于传统方法的无人机测量与地面测量株高结果对比Tab.2 Comparison of plant height based on traditional UAV measurement and ground measurement cm

图2 冠层点云分布示意图Fig.2 Distribution of point cloud in canopy

在空间分布上,试验区小麦种质资源丰富,不同品系之间小麦固有发育速度有较大差异,最慢品种和最快品种的日株高差可达146 cm,这导致无人机每次采样的株高数据之间存在较大梯度。以表2中2022年3月11日为例,虽然相对均方根误差(RMSE100)为59.44%,但相关系数仍然处于较高的显著性水平,进一步说明时序上出现株高逆增长现象,是设备与算法的固有误差被相对放大而导致。

2.3 小麦株高精确测量方法

基于前述分析,本文在缺失地面控制点高程信息的条件下,基于空间分布上已有部分随机采样真值,提出一种在育种小区粒度下的近邻校正法和一种基于多光谱+RGB数据融合的半经验光谱指数校正法。试验区整体平整,南北海拔差异为0.8 m、东西海拔差异为0.3 m,近邻小区南北海拔落差在0.05 m以内、东西海拔落差在0.03 m以内,与各小区群体后期生长发育不齐带来的株高落差相比,可以忽略不计,这为近邻校正法提供现实基础。此外,基于不同波段光谱反射率组合的植被指数(VIs)能够有效反映植物生长状况[25],以及无人机遥感能够提供高分辨率的VIs测量,为光谱指数校正法提供了理论依据,故而可以在育种小区粒度下选取合适的VIs建立回归模型预测小麦株高。图3展示了数据预处理流程与两种模型的建模思路。

图3 总体技术流程图Fig.3 Overall technical roadmap

2.3.1基于高分辨率RGB图像的近邻校正法

本方法在出苗之前采集试验田裸地数据的DTM信息,并采用下述步骤计算小麦株高。

(1)小区群体高程获取

不同采样时间下,由于RTK信号会有轻微偏移,因此将不同日期各小区的遥感影像进行经纬度配准。依据地面控制点进行配准,配准方法采用一阶多项式,其又称为仿射变换,除了平移、缩放和旋转之外,具有纠正地图变形的功能,保证每个日期配准残差小于1个像素。将获取到的每两期图像做差,得到作物株高数据(Wheat plant height, WPH),计算式为

(1)

式中Sday——数字高程数据

T——裸土数字高程数据

(2)结合地埂高程校正

为解决由于设备与算法的固有误差被相对放大而导致株高时序上出现逆增长问题,利用遥感影像中育种小区之间裸露地埂的土壤基底信息进行修正。由于小区分割时感兴趣区域可能包含部分未种植的裸地信息,将每个生长阶段DSM值的95%百分位数视为小麦顶部,将DSM值的5%视为裸土高程DEM,得到WPH的估计值,计算式为

(2)

其中

式中Eday——数字高程数据

α——算法调节因子

(3)依据真值滑动校正

在小麦生长发育后期,叶片会将小区间的地埂覆盖,对无法通过式(1)、(2)修正的小区,采用与其邻近的地面株高真值修正土壤基底信息。具体地,UAV测量的冠层海拔减去地面真值,得到各小区土壤海拔,并将这一海拔共享给周围n个近邻,n个近邻的最终株高则为UAV测量的冠层海拔扣除上述共享海拔,将n+1个小区作为一个滑动窗口依次计算所有小区的株高(图4)。

图4 近邻小区海拔共享Fig.4 Altitude sharing of adjacent cells

(3)

其中

2.3.2基于多光谱和RGB数据融合的光谱指数法

根据采集的五通道多光谱数据,计算植被指数及其与株高的相关系数,优选合适的植被指数与DSM数据参与建模得出株高回归模型。

(1)VIs计算

多光谱遥感影像生成的植被指数地图中像素点灰度变化能够反映株高生长特征的差异,故对各育种小区发育过程中的植被指数求取均值。计算与作物高度相关程度高的15种植被指数[26-28]:归一化植被指数(NDVI)、蓝色归一化差异植被指数(BNDVI)、蓝色归一化差异红边植被指数(BNDREI)、归一化绿蓝差异指数(NGBDI)、标准蓝色指数(NormB)、绿叶指数(GLI)、绿/红比率指数(GR)、过绿指数(ExG)、修正简单比率(MSR)、红边叶绿素指数(CIRE)、改良叶绿素吸收率指数(MCARI)、转化叶绿素吸收反射指数(TCARI)、三角植被指数(TVI)、根号归一化植被指数(RNVI)、非线性植被指数(NLI)等。

V=f1(Blue,Green,Red,Red_edge,NIR)

(4)

式中,V为多种植被指数,Blue、Green、Red、Red_edge、NIR为五通道光谱反射率。

(2)植被指数优选

分别在时序上计算植被指数向量与株高向量的皮尔逊相关系数和空间上各育种小区的植被指数向量与株高向量的皮尔逊相关系数,形成三维相关系数矩阵,并在合适的显著性水平下选出用于下一步建模的VIs[25-27]。

(5)

式中Vday——日期day植被指数向量

Gday——日期day株高向量

(3)VIs-WPH模型构建

将选择出的VIs与DSM数据进行归一化处理,分别采用偏最小二乘g1(V,S,b)、支持向量机g2(V,S,b)、随机森林g3(V,S,b)3种回归算法,计算校正后的株高,再进行模型的选择。

Hday=min{RMSE(g(V,S,b))}

(6)

其中

g(V,S,b)=

{g1(V,S,b),g2(V,S,b),g3(V,S,b)}

式中Hday——日期day小麦株高

S——数字高程数据

b——调节系数

RMSE(g(V,S,b))——计算g(V,S,b)RMSE的函数

2.4 准确性评价指标

采用决定系数R2、均方根误差(RMSE)、相对均方根误差(RMSE100)、相关系数R评价模型的稳健性和准确性。R2越高拟合效果越好,RMSE越小,测量越精确。107个采样点随机分为80%的训练集和20%的验证集,采用五折交叉验证的方法以获得可靠结果。

3 结果与分析

3.1 植被指数与株高的相关性分析

根据相关性系数显著性检验公式可知,自由度为107(即样本量)时,0.02、0.01、0.05、0.001显著水平的临界值分别为0.224 7、0.248 0、0.269 5、0.313 7,即相关性系数绝对值大于0.248 0时达到极显著水平。表3展示了15种植被指数在不同时期与小麦株高的皮尔逊相关系数,大部分植被指数与株高的相关系数从2022年3月28日才开始出现绝对值大于0.313 7,然后随时间推移不断增加,在2022年4月5日达到最高后开始不断减小,值得注意的是过绿指数(Excess green,ExG)在整个时期均表现出较高的显著性水平;在2022年3月28日至2022年4月13日这个时间段大部分植被指数与株高的相关性维持在较高的显著性水平上,说明在挑旗期、抽穗期的植被指数能很好地反映株高变化,而在其它生长期与株高的相关性较低,难以精确反映小区粒度上的株高变化。

表3 不同时期植被指数与株高的皮尔逊相关系数Tab.3 Correlation coefficient between vegetation index and plant height in different periods

3.2 不同模型下的小麦株高时序变化

表4给出了TM、NNCM以及SICM的株高估算结果,可以看出传统测量方法在2022年3月11日和2022年4月27日出现较大偏差,NNCM的测量精度在2022年4月5日达到最高,其R为0.94,RMSE为3.72 cm,RMSE100为5.62%,此时大部分育种小区完成拔节,株高处在40～100 cm范围,动态变化快、空间差异大,有利于数据校正;测量精度在2022年4月13日开始降低,主要原因在于大部分品种的株高已经停止变化,营养用于叶片和生殖器官的发育,导致育种小区间地埂被逐步覆盖,降低了算法的精度。SICM的测量精度在2022年3月28日达到最高,其R2为0.79,R为0.89,RMSE为3.25 cm,RMSE100为6.05%,此时小麦处于快速生长期,很高的R值说明植被指数能够较好地反映小麦株高变化。从时序上看,2022年3月28至4月13日期间NNCM和SICM结果相比于其它时期具有更高的测量精度,说明这2种方法在小麦挑旗期与抽穗期取得了最好的应用效果。为了进一步观测小区粒度株高的反演情况,得到单个小区的模型估测结果与真值的散点图如图5所示,不同列表示不同方法的株高反演情况,不同行表示不同时期的株高反演情况。可以看出,在2022年3月11日、2022年4月27日两个时期无人机传统测量方法估测的小麦株高明显低于地面测量真值;相比于TM,NNCM与SICM观测点拟合直线斜率k在6个时期均接近1,拟合直线截距更接近0,各观测点小区更靠近坐标轴的45°线,说明NNCM与SICM能够很好地测量小麦株高。SICM与NNCM相比在2022年4月27日结果较差,这是由于在小麦株高生长后期植被指数并不能够很好地反映小麦株高的变化,故可以依据不同的实验时期以及实验条件,合理选择SICM模型或NNCM模型。

表4 不同日期不同模型的评价指标计算结果Tab.4 Results of different models on different dates

图5 模型估测结果精度及时序变化分析Fig.5 Accuracy of model estimation results and time series change analysis

3.3 不同生育期的小麦株高空间分布

图6为小区细粒度下采用NNCM绘制的小麦株高时序变化空间分布,不同颜色代表不同株高,2022年3月7日多数小区处于橘黄色,将图中不同小区的颜色与右侧参考值对比,结合图5可以看出在2022年3月7日小麦株高普遍在20～30 cm。此外,原图在QGIS软件中可以查看每个像素点的具体株高。

图6 NNCM测量小麦株高时序变化的空间分布Fig.6 Spatial distribution of wheat plant height

在2022年3月28日试验田北部株高在35～50 cm,南部株高在45～60 cm,表明南部种植的小麦品种相较于试验田北部小麦品种有更快的株高生长速率;在2022年4月13日试验田,南部种植的小麦株高比北部种植的小麦品种高20～40 cm,说明NNCM能够精确反演小区粒度尺度下的株高及不同小区株高在空间中的变化,为后续全基因组关联分析提供精确度高、通量大的小麦田间株高表型数据。

3.4 不同模型的性能对比分析

图7给出了不同时期株高真值GM及使用两种校正方法后的分布结果,对比发现,NNCM与SICM核密度分布曲线与地面真值都很接近,中位数、四分位数、最大值、最小值偏差均不超过0.5%,说明这两种方法解决了细粒度株高精确测量的问题,有效校正了固有误差带来的数据逆增长现象。图7小提琴图的轮廓即核密度曲线,反映了株高性状的分布,可以看出,本试验区株高性状服从正态分布,种质资源丰富,群体间株高性状具有较大的差异性。此外还可以看出,在2022年4月13日株高增长较快的小麦已停止增长,较2022年4月5日采样结果GM中位数增长19 cm、TM模型中位数增长10 cm、NNCM模型中位数增长17 cm、SICM模型中位数增长17 cm,相比于TM模型,NNCM和SICM能够很好地同步跟随GM观测到的株高变化。

图7 实测株高与模型校正后的株高分布Fig.7 Measured plant height and plant height distribution after model correction

4 讨论

本文基于无人机遥感,通过TM、NNCM以及SICM提取小麦冠层株高,并与6个时期的实测株高对比,表明无人机遥感在小麦育种小区粒度下,基于人工测量真值的NNCM和基于多光谱+RGB数据融合的SICM能够精确获取小麦株高。NNCM获取精度略高于SICM,在具有地面随机采样真值的条件下,获取成本低于SICM;在不具有地面随机采样真值的条件下应当选用SICM, NNCM在实际农田的效率可能不如SICM,应当依据实际情况合适选择测量方法。此外,无人机数据采集时,天气变化、风的有无会对数据质量造成较大影响,所以要在天气晴朗、无风或微风的条件下采集数据。下一步将在不同试验田、不同采集设备条件下探究本文所提模型的泛化能力,本文所提思想还可以为小麦其他表型参数在细粒度尺度下精确测量提供借鉴。

5 结论

(1)验证了基于细粒度校正的育种小区小麦株高无人机测量是可行的,本文提出的两种株高精确获取方法能够在育种小区细粒度下精确获取小麦株高,两种方法的RMSE100可达6.00%左右,其中NNCM的RMSE100最低可达5.70%;SICM的RMSE100最低可达6.05%。

(2)NNCM和SICM对真实株高动态变化均具有良好的跟随性,可以精确测量育种小麦小区粒度的全生育期株高,避免了UAV高频次飞行传统测量所出现的株高逆增长现象,在实际生产实践与辅助育种中有较好的应用前景,可为细粒度下小麦株高无人机遥感精确获取提供参考。