基于超体素聚类和局部特征的玉米植株点云雄穗分割
1.
2.
3.
Tassel Segmentation of Maize Point Cloud Based on Super Voxels Clustering and Local Features
1.
2.
3.
通讯作者: 苗 腾(1985-),男,副教授,研究方向为数字植物技术和表型检测技术。电话:
收稿日期: 2021-02-01 修回日期: 2021-02-23 网络出版日期: 2021-06-01
基金资助: |
|
Received: 2021-02-01 Revised: 2021-02-23 Online: 2021-06-01
作者简介 About authors
朱超(1990-),男,博士,研究方向为作物表型检测技术。E-mail:
针对当前三维点云处理方法在玉米植株点云中识别雄穗相对困难的问题,提出一种基于超体素聚类和局部特征的玉米植株点云雄穗分割方法。首先通过边连接操作建立玉米植株点云无向图,利用法向量差异计算边权值,并采用谱聚类方法将植株点云分解为多个超体素子区域;随后结合主成分分析方法和点云直线特征提取植株顶部的子区域;最后利用玉米植株点云的平面局部特征在顶部子区域中识别雄穗点云。对3种点云密度的15株成熟期玉米植株点云进行测试,采用F1分数作为分割精度判别指标,试验结果与手动分割真值相比,当点云密度为0.8、1.3和1.9个点/cm时,雄穗点云分割的平均F1分数分别为0.763、0.875和0.889,分割精度随点云密度增加而增高。结果表明,本研究提出的基于超体素聚类和局部特征的玉米植株点云雄穗分割方法具备在玉米植株点云中提取雄穗的能力,可为玉米高通量表型检测、玉米三维重建等研究和应用提供技术支持。
关键词:
Accurate and high-throughput maize plant phenotyping is vital for crop breeding and cultivation research. Tassel-related phenotypic parameters are important agronomic traits. However, fully automatic and fine tassel organ segmentation of maize shoots from three-dimensional (3D) point clouds is still challenging. To address this issue, a tassel point cloud segmentation method based on point cloud super voxels clustering and local geometric features was proposed in this study. Firstly, the undirected graph of the maize plant point cloud was established, the edge weights were calculated by using the difference of normal vectors, and the spectral clustering method was used to cluster the point cloud to form multiple super voxel sub-regions. Then, the principal component analysis method was used to find the two end regions of the plant and based on the observation of the straight direction of the bottom stem regions, the top and bottom regions were distinguished by the point cloud linear features. Finally, the tassel points were identified based on the plane local features of the point cloud. The sub-regions of the top region of the plant were classified into leaf regions, tassel regions, and mixed regions by plane local features of the point cloud, the tassel points in the tassel sub-region, and the mixed region were the finally segmented tassel point clouds. In this study, 15 mature maize plants with 3 point cloud densities were tested. Compared with the ground truth segmented manually, the average F1 scores of the tassel segmentation were 0.763, 0.875 and 0.889 when the point cloud density was 0.8/cm, 1.3/cm, and 1.9/cm, respectively. The segmentation accuracy of this method increased with the increase of plant point cloud density. The increase of point cloud density and the number of point clouds mainly affected the calculation results of point cloud plane features in tassel segmentation. When the number of point clouds was small, the top leaf point cloud was relatively sparse. Therefore, the difference between the plane feature of the leaf point and the plane feature of the tassel point was not obvious, which led to the increase of the misclassification of the point cloud. However, the time complexity of the algorithm was O(n3), so the increase in the density and number of point clouds would lead to a significant increase in the running time. Considering the segmentation accuracy and running time, the research obtained the best effect on the mature maize plants with a point cloud density of 1.3/cm and an average number of 15,000. The segmentation F1 score reached 0.875 and the running time was 6.85 s. The results showed that this method could extract tassels from maize plant point cloud, and provided technical support for the research and application of high-throughput phenotyping and three-dimensional reconstruction of maize.
Keywords:
本文引用格式
朱超, 吴凡, 刘长斌, 赵健翔, 林丽丽, 田雪莹, 苗腾.
ZHU Chao, WU Fan, LIU Changbin, ZHAO Jianxiang, LIN Lili, TIAN Xueying, MIAO Teng.
1 引 言
基于二维图像的视觉技术率先被应用于雄穗表型提取中。Gage等[7]开发了一种基于图像的玉米雄穗表型分析系统,可对离体的雄穗进行标准图像获取,进而提取穗长、穗宽等表型参数。一些研究工作聚焦于在玉米植株或冠层图像中识别雄穗。Ye等[8]结合方向梯度直方图特征和支持向量机分类器对雄穗区域进行识别,再采用显著性检测方法实现雄穗像素分割。Kurtulmus和Kavdr[9]结合颜色特征和支持向量机分类器进行雄穗定位,并采用层次聚类方法实现雄穗分割。Lu等[10]采用多个角度的图像提升雄穗分割的精度,并实现了对雄穗颜色、分支数等参数的粗略测量。Lu等[11]构建了雄穗识别深度学习网络TasselNet,实现了在玉米冠层图像中对雄穗的准确计数。二维图像只是三维空间的一个投影,利用二维图像可以较好地进行雄穗计数工作,但对穗长、穗宽等参数的提取结果并不准确。
目前,在玉米植株或群体点云中对雄穗表型参数进行提取仍具挑战,其中,对雄穗点云的准确分割是技术难点。针对该问题,本研究提出一种基于超体素聚类和局部特征的玉米植株点云雄穗分割方法,可直接在玉米植株点云上对雄穗点云进行分割和识别,为玉米高通量表型检测、三维几何重建等研究和应用提供技术基础。
2 材料与方法
2.1 试验材料与点云数据获取
田间试验于2019年5月至10月在沈阳农业大学玉米试验田(116.16°E,39.56°N)进行。选取种植行距60 cm、株距25 cm的成熟期玉米植株为研究对象,使用FreeScan X3手持激光扫描仪进行玉米植株三维点云获取。将田间成熟期玉米植株移植到盆中,并转移到室内进行点云获取。由于FreeScan X3激光扫描仪需要借助激光感光片进行三维定位,因此,在采集数据时,先将感光片粘贴在辅助支架上,并将支架围在被测植株周围,保证激光线能同时扫描到植株和感光片,数据采集过程如图1(a)所示。利用CloudCompare软件对获取的原始点云数据进行预处理,手动删除盆、移动支架的点云数据,仅保留玉米植株点云数据。
为验证本研究分割算法对不同密度点云植株的分割效果,选取15株不同形态的成熟期玉米植株点云作为测试样本,每个玉米植株点云均被采样到3个点云密度(如图1(b)),分别为0.8、1.3和1.9个/cm,点云密度表示每厘米存在多少个点;3个密度的植株点云平均数量约为6000、15,000和30,000个点,雄穗点云平均数量约为140、360和890个点。利用CloudCompare软件对被测玉米植株点云进行手动分割获得雄穗点云的分割真值。
图1
图1
点云数据获取和预处理
Fig.1
Point cloud data acquisition and preprocessing
2.2 雄穗点云分割方法
本研究方法以玉米植株散乱点云为输入,通过点云超体素聚类、植株顶部子区域提取和雄穗点云识别3个主要步骤,实现对雄穗点云的识别和分割。
2.2.1 点云超体素聚类
本研究首先通过点云超体素聚类,将植株点云分解为多个点云子区域,这些子区域由具有几何相似性的局部区域点云整合形成,并且这些子区域较好地保留了植株点云器官间的边界信息,便于后续进行雄穗点云的分割。本研究采用谱聚类算法实现点云超体素聚类,假设输入的三维玉米点云为P ={pi|i=1,2,…,n},点数量为n,对P进行超体素聚类的具体步骤如下。
(1)边连接操作。对P中的所有点均进行边连接操作形成无向图。对于任意一点pi,首先通过K近邻空间搜索方法[26]获得P中与点pi欧式距离最小的k个近邻点,形成点集Ai ={aj| j=1,2,…,k,aj ∈ P}。之后遍历Ai中的所有点,判断每个点是否能与pi连接形成边。如果Ai中的点aj在P中的k个近邻点不包含pi,则pi与aj不能连接成边;反之则可连接形成1条边,且该边的权重值e通过
(2)构建归一化拉普拉斯矩阵。构建植株点云P的邻接矩阵W,该矩阵为n×n大小的方阵,其第i行、k列的矩阵元素wik表示点pi与pk之间的连接关系,如果pi与pk之间无边连接,则wik=0;如果pi与pk之间有边连接,则wik等于该边的权重值。构建n×n大小的度矩阵D,该矩阵第i行、k列的矩阵元素值dik通过
从
(3)点云聚类。采用奇异值分解方法计算归一化拉普拉斯矩阵L的特征值和特征向量,并选出最小的K个特征值对应的特征向量ui(i =1, 2, …, K),每个特征向量ui均为n维列向量,K为用户输入的类别数量参数,本研究设K=80。将上述K个特征向量组合成n×K大小的特征矩阵F =[u1, u2, …, uK],并将F按行标准化。如果将F中的每个K维行向量看作1个数据样本,则F矩阵中共有n个数据样本。采用K-means聚类方法将这n个数据样本分为K个类别,其中第i个数据的类别即为植株点云P中第i个点pi的类别。最后将类别相同的点划分为同一子区域,则植株点云最终被划分为K个子区域点集。用符号Ø表示子区域点集,第m个子区域用符号Øm(m=1,2,…,K)表示。本研究以点云法向量作为特征,将法向量相似的邻近点云聚类成同一子区域。
不同株型植株的点云超体素聚类结果如图2所示,从图中可以看出,植株的每个器官均被分解成多个边界光滑的子区域,且每个子区域的整体形态相对简单,便于用简单的几何特征进行表示和识别。例如,叶片上的子区域近似平面,茎上的子区域近似圆柱。但在个别器官交界处的位置,会存在分属不同器官的点云被误划分为同一子区域的现象,如部分叶领点云与茎形成同一子区域、部分雄穗点云与顶端叶片形成同一子区域。
图2
图2
不同株型植株的点云超体素聚类可视化图
Fig.2
Visualization of point cloud super voxels clustering for different plant types
2.2.2 植株顶部子区域提取
(1)拟合植株点云的主成分向量。采用主成分分析法提取植株点云P的主成分方向。首先求植株的中心点pc。
之后利用
其中,矩阵C为3×3对称实方阵。
采用奇异值分解算法对C进行特征值分解,得到3个非负特征值λ1、λ2、λ3 (λ1≥λ2≥λ3)以及对应的3个特征向量I1、I2、I3。特征向量I1、I2、I3分别为植株点云的第1、2、3主成分向量。
(2)提取候选的顶部子区域集合。以植株点云第1主成分向量I1为方向向量,过植株中心点pc的直线参数方程L(t)为:
利用
找到所有点投影坐标的最大值t1和最小值t2,之后利用
定义子区域集合V1和V2,并遍历所有点云子区域。如果子区域中存在投影坐标>t3的点云,则将该子区域加入集合V1。如果子区域中存在投影坐标<t4的点云,则将该子区域加入集合V2。V1和V2均为候选的植株顶部子区域集合,表示玉米植株的两端,其中一个为植株底部区域,另一个为顶部区域。提取候选顶部区域的过程如图3(a)所示。
图3
图3
植株顶部子区域提取
Fig. 3
Sub-region extraction at the top of the plant
(3)识别顶部子区域。成熟期玉米植株底部区域的茎为圆柱形,具有较明显的直线线性特征。与之相反,植株顶部区域包含了雄穗,呈现一定的扫帚状,直线线性特征较弱。基于上述特点,本研究通过计算子区域中点云的直线线性特征从V1和V2中选出顶部子区域集合。
定义1个子区域点集的直线线性特征为集合内每个点的直线线性特征的平均值。为计算子区域点集Ø中一点p的直线线性特征,首先在Ø中查找与p点距离小于r1的所有点,组成点集B1。其中r1为用户可调节参数,本研究设r1=10 cm;之后采用主成分分析法,计算B1点集的3个特征值λ1、λ2、λ3 (λ1≥λ2≥λ3),则p点的直线线性特征f1通过
利用上述方法分别计算V1、V2中每个子区域的直线线性特征,并求V1中所有子区域直线线性特征的平均值fv1以及V2中所有子区域直线线性特征的平均值fv2,比较fv1和fv2的大小,较小值对应的子区域集合为植株顶部区域,较大值对应的区域集合为植株底部区域。V1、V2集合的子区域直线线性特征热力图如图3(b)所示,从图中可以看出,顶部区域的雄穗子区域具有非常小的直线线性特征值,导致该候选区域的直线线性特征的平均值大幅度降低。而茎则具有较大的直线线性特征值,从而使底部区域的特征值保持在较高的数值。
2.2.3 雄穗点云识别
上节提取的顶部子区域集合(图4(a))包含了3类子区域,分别为只包含雄穗点的雄穗子区域、只包含叶片点的叶片子区域以及同时包含雄穗点和叶片点的混合子区域。为实现雄穗点云的精确分割,需要在顶部区域中识别出雄穗子区域,并在混合子区域中识别出雄穗点云。由于叶片比雄穗具有更明显的平面特征,因此本研究通过子区域内点云的平面特征来识别雄穗点云和叶片点云,具体步骤如下。
图4
图4
雄穗点云识别示意图
Fig. 4
Schematic diagram of the tassel point cloud recognition
(1)提取子区域内点云的平面特征。为计算子区域点集Ø中一点p的平面特征,在Ø中找到与p点距离小于r2的所有点,组成点集B2,其中r2为用户可调节参数;之后计算B点集的3个特征值λ1、λ2、λ3 (λ1≥λ2≥λ3),则p点的平面特征f2(p)通过
子区域点集Ø的平面特征f2(Ø)为集合内每个点的平面特征的平均值。植株顶部点云平面特征的热力图如图4(b)所示,特征值越大,表示平面特征越弱。
(2)点云分类。计算所有子区域平面特征的平均值
(3)去噪。在混合子区域中,叶片边缘的点云的平面特征也较小,因此在分类时容易被分类成雄穗点。但通常情况下,这些误分的点与整个雄穗点云的距离较远,因此本研究采用基于统计分析的点云去噪方法[30]从雄穗点云中排除这些误分点。
2.3 分割精度评估
本研究使用精确度(Precision)、召回率(Recall)、F1分数(F1 score)对分割精度进行评估。对于雄穗点云,用TP表示属于雄穗的点云被正确分割为雄穗的点数量,用FP表示不属于雄穗的点云被误分割成雄穗的点数量,用FN表示属于雄穗的点云被误分割为其他器官点云集合的点数量,则精确度、召回率和F1分数的计算如公式(
其中,Pr表示精确度;Re表示召回率;F1表示F1分数。一般来说,分割精确度和召回率之间是矛盾的,因此常采用F1分数来对分割算法效果进行综合考量,越高的F1分数说明分割效果越好。
3 结果与分析
3.1 分割效率和精度
在配置为2.2 GHZ CPU、DDR32 G内存的笔记本工作站上进行了算法测试。本研究在进行点云超体素聚类时,利用奇异值分解方法对n×n大小的归一化拉普拉斯矩阵进行特征向量分解,导致整个算法的时间复杂度达到O(n3)。因此随着点云密度和点云数量的增长,算法的运行时间会显著增加。当点云密度为0.8个/cm时,成熟期玉米植株点云个数平均在6000个左右,此时算法的平均总运行时间为2.26 s;点云密度为1.3个/cm时,植株点云数量平均约为15,000,平均总运行时间为6.85 s;当点云密度增大到1.9个/cm时,植株平均点云数量增加到30,000个左右,此时算法的总运行时间增长到27.37 s。本研究方法的具体运行时间如表1所示。
表1 算法运行时间
Table 1
点云密度/(个·cm-1) | 超体素聚类平均运行时间/s | 植株顶部子区域提取平均运行时间/s | 雄穗点云识别平均运行时间/s | 平均总运行时间/s |
---|---|---|---|---|
0.8 | 1.58 | 0.53 | 0.15 | 2.26 |
1.3 | 4.61 | 1.77 | 0.47 | 6.85 |
1.9 | 17.78 | 8.11 | 1.48 | 27.37 |
本研究算法的分割精度如表2所示。点云密度为0.8个/cm时,雄穗的平均F1分数、平均精确率和平均召回率分别为0.763、0.746和0.832;点云密度为1.3时,上述3个指标分别为0.875、0.898和0.867;点云密度为1.9个/cm时,3个分割指标分别为0.889、0.950和0.849。从F1分数指标可以看出,随着点云密度增加,本研究的分割精度逐渐增加,但点云密度1.3个/cm和1.9个/cm之间的分割效果差距不大。
表2 雄穗分割精度
Table 2
点云密度/(个·cm-1) | F1分数 | 精确率 | 召回率 | ||||||
---|---|---|---|---|---|---|---|---|---|
最小值 | 平均值 | 最大值 | 最小值 | 平均值 | 最大值 | 最小值 | 平均值 | 最大值 | |
0.8 | 0.531 | 0.763 | 0.911 | 0.370 | 0.746 | 0.993 | 0.607 | 0.832 | 1.000 |
1.3 | 0.714 | 0.875 | 0.976 | 0.568 | 0.898 | 1.000 | 0.613 | 0.867 | 1.000 |
1.9 | 0.737 | 0.889 | 0.998 | 0.802 | 0.950 | 1.000 | 0.584 | 0.849 | 1.000 |
点云密度和点云个数的增加,主要影响雄穗分割步骤中的点云平面特征计算结果。当点云个数较少时,顶部叶片点云也相对稀疏,叶片点的平面特征和雄穗点的平面特征差距不明显,造成点云类别误判的情况增大,分割出的雄穗点云中包含较多叶片点云,分割精确度不高;当点云个数增加时,叶片点云的平面特征更加明显,因此,雄穗点和叶片点的区分度加大,分割出的雄穗点云中包含的叶片点云减少,分割精确度显著提高。从表1的精确度指标可以看出,随着点云密度的增大,精确度显著提高,而召回率则变化不大,F1分数的增长主要来自于精确度的增加。
图5
图5
不同F1分数的雄穗分割可视化结果
Fig. 5
Visualization results of tassel segmentation with different F1 score
3.2 算法参数
影响本研究算法的参数主要为点云超体素聚类步骤中的K参数、植株顶部子区域提取中的r1参数和雄穗点云识别中的r2参数。K参数会影响点云超体素聚类得到的子区域个数,K值越大,子区域个数越大。本研究在后续植株顶部子区域提取和雄穗点云识别操作中,均以子区域作为点云集合,计算点的直线线性特征和平面特征。K值越小,每个器官被分成的子区域个数越少,每个子区包含多个器官点云的概率会增大,造成直线线性特征和平面特征对茎、叶片、雄穗点的区别度下降,降低分割精度。但K值越大,超体素聚类的计算时间会显著增加。本研究中,通过试验确定参数K=80,该值对于不同株型、不同密度的玉米成熟期点云均适用,可将每个器官分解成多个边界光滑且形态简单的子区域,便于利用直线线性特征和平面特征对点云类别进行判断。
r1参数会影响点的直线线性特征计算结果。如果r1值较小,底部子区域茎点云在计算一个点的直线特征时,搜索得到的邻域点仅包含茎圆柱面上的一个侧面,使计算得到的直线特征值减少。相反,r1值越大,茎点云的直线特征值越大,底部子区域和顶部子区域的直线特征对比越明显,但过大的r1值会增多邻域点个数,降低特征值的计算效率。本研究设r1=10.0 cm,由于玉米茎的直径通常不会大于10.0 cm,因此该值可以使每个茎点的邻域点集合均为一个较完整的圆柱状点云集合,使茎点云的直线特征值较大,同时,也不会过多地降低计算效率。r2参数会影响点的平面特征计算结果,进而改变雄穗点云的识别结果。玉米叶片的平面特征具有局部性,即整个叶片表面呈现较明显的弯曲特征,但局部小区域是较平展的。随着r2值的增大,每个点搜索得到的邻域点范围会增大,则叶片点的平面性特征值会增大,降低叶片点和雄穗点之间的平面特征对比度,因此r2值不宜过大。由于顶部子区域包含的叶片均为上部叶的叶尖区域,该区域叶片宽度通常不会大于6.0 cm,同时顶部叶片经常沿着叶脉卷曲,因此r2值要小于3.0 cm,这样才能尽可能地保证叶片点云的平面特征较明显。但过小的r2值会使邻域点个数过少,很难描述局部区域的几何特性,降低平面特征计算的准确性,同时,雄穗的每个分支在极小区域内也可能会呈现一定的平面特征,因此要保证r2值不能小于雄穗分支的宽度值。玉米雄穗每个分支的宽度通常会小于1.0 cm,因此r2值要大于1.0 cm。
本研究在设定r2具体参数值时,针对每一个点云密度的植株,均从r2 = 1.0 cm开始进行参数选择测试,按照间隔0.1 cm逐渐增大r2,直到r2 = 3.0 cm。选择分割效果较好的r2参数进行计算,当多个参数值分割效果差距不大时,选择其中的最小值作为r2参数,从而尽可能地减少邻域点数量,提高平面特征计算效率。本研究通过该方式确定r2值,当点云密度= 0.8个/cm时,r2设为2.4 cm;点云密度= 1.3个/cm时,r2设为1.8 cm;点云密度=1.9个/cm时,r2设为1.0 cm。随着点云密度的减少,r2参数值明显增大,这是由于点云密度较少时,为了找到足够多的邻域点描述局部平面特征,需要增大邻域搜索半径。当给定一个新密度的植株点云时,可以利用本研究给出的r2值作为参考提高参数选择效率。例如,如果植株点云密度为1.6个/cm,则r2值可在点云密度为1.9和1.3个/cm的参数值之间进行测试,即r2值从1.0 cm开始进行参数选择,按照间隔0.1 cm逐渐增大r2到1.8 cm。
3.3 展望
本研究算法目前不能较好地处理雄穗被叶片包裹住或距离较近的情况,此时叶片和穗之间的边界点云的法向量和平面局部特征非常相似,造成分割精度下降。Jin等[31]开发了基于中值归一化向量的区域生长茎叶分割算法,利用欧式距离作为特征,采用启发式区域增长策略从器官底部的种子点开始,逐步将器官点云进行完整分割,该方法证明了距离特征在玉米器官分割时的有效性。在后续的工作中,以本研究分割出的穗点云作为种子点,结合基于距离条件限定的区域增长方法,有望优化最终的穗点云分割结果。
本研究算法只能对单株玉米植株点云进行雄穗分割。为获得单株点云,需将玉米植株移植入盆并采用激光扫描仪手动进行数据获取,数据获取效率低。MVS Pheno平台可快速、自动地获取不用角度的玉米植株图像,并利用运动恢复结构方法获得质量较好地玉米单株点云[32]。下一步可将本研究方法与MVS Pheno平台进行结合,加快雄穗表型检测速度。
4 结 论
本研究提出一种基于超体素聚类和点云局部特征的玉米植株三维点云雄穗分割方法,验证了在玉米成熟期植株点云上识别、分割雄穗的可能性。通过分析本研究试验结果,得出以下结论。
(1)本研究方法的分割精度随植株点云密度的增大而增加。当点云密度为0.8、1.3和1.9个/cm时,雄穗点云分割的平均F1分数分别为0.763、0.875和0.889。
(2)本研究算法的时间复杂度为O(n3),随着点云密度和点云数量的增长,算法的运行时间会显著增加。
(3)综合考虑分割精度和运行时间,本研究对点云密度为1.3个/cm、点云的平均数量约为15,000的成熟期玉米植株点云处理效果最好,分割F1分数可达到0.875,运行时间为6.85 s。
结果表明,本研究方法具有在玉米植株散乱点云中提取雄穗点云的能力,可为玉米高通量表型检测、玉米三维重建等研究和应用提供技术手段。