张卫正，李旭光，万瀚文，李灿林，张伟伟，金保华，刘 岩

(1.郑州轻工业大学 计算机与通信工程学院，河南 郑州 450002； 2.郑州大学国际学院，河南 郑州 450000)

玉米(ZeamaysL.)是禾本科一年生草本植物，具有较强的耐旱性、耐寒性、耐贫瘠性以及极好的环境适应性，与水稻、小麦并称为世界三大粮食作物[1-2]。玉米的营养价值较高，是畜牧业、养殖业等的重要饲料来源，也是食品、医疗卫生、轻工业等不可或缺的原料之一[3-4]。玉米的表型分析一直是世界各地植物遗传育种研究的热点，表型是研究“基因型-表型-环境”作用机制的重要桥梁，高通量的表型检测有助于辅助作物育种，加速整个环节和进程[5]。

传统的表型提取方法，在获取表型信息时容易损伤植株，并且通常具有劳动力密集、耗时、成本高等缺点。而基于图像的植物表型分析，则可在相对较短的时间内大批量、无损地分析并提取植株性状[6]。可通过将整个植物视为单个对象(整体表型)或作为单个组分(如叶和茎等组分表型)来计算高级表型，以研究玉米等植物的生物物理特征[7]。

植物活力可以通过个体叶片随时间的生长来表征，叶片长度和连接尖端距离是2种重要的表型，有助于确定植物对环境胁迫的反应。叶片曲率是叶子韧性的量度，叶片下垂可作为营养缺乏等反映植株活力的指标。茎叶角是玉米植株主茎与垂直线的夹角，可作为监测植株因风力等自然灾害而倒伏的早期敏感信号。基于计算机视觉的叶片研究可取代人工及机械测量装置，显著降低人工操作的劳动强度[8]。

在玉米植株的茎叶分割方面，苏伟等[9]采用激光雷达获取大田玉米植株的点云，通过Cloud Compare软件裁切工具手动分离出单株玉米，利用玉米叶片和茎秆表面法线的差异，采用DoN(Difference of normals operator)法线差分算法实现茎叶的分割。黄成龙等[7]利用并行细化算法对完整的二值图像提取植株骨架，并去除细小毛刺，采用霍夫变换对玉米骨架图像检测竖直方向上的直线，识别骨架中茎秆部分。根据叶片的联通性和位置信息，对叶片骨架进行识别并用不同的伪彩进行标记，实现玉米植株的茎叶分割。BRICHET等[10]首先从12个侧视图中选择玉米植株展开度最好的视图，然后采用均值平移和HSV阈值算法提取玉米植株，接着采用Scikit-image库的中轴算法提取植株的骨架，使用最短路径算法来检测主茎，从而实现茎叶分割。CABRERA等[11]使用细化算法从二值化图像中提取植株骨架，采用Imagej中的插件Analyze Skeleton进一步获得主茎的中间节点和叶片分支的端点，从而实现茎叶分割。CHOUDHURY等[12]使用快速行进距离变换来对二值图像进行骨架化，之后采用基于阈值的骨架修剪来去除伪分支，然后检测植株骨架中的节点和端点，根据植株的先验知识完成茎叶分割。

将机器视觉和图像处理技术应用于农业研究与生产，具有广泛的应用前景和实用价值。本研究采用已有的骨架提取算法，结合二叉树分析进行玉米植株的茎叶分割，以求处理流程便捷有效。研究将有助于基于计算机视觉的玉米植株以及其他具有相似结构的谷物作物(例如高粱)的表型分析研究。

1 材料和方法

1.1 研究数据

本研究采用内布拉斯加大学林肯分校的植物表型数据集(Nebraska-Lincoln Component Plant Phenotyping Dataset,UNL-CPP)，可从https://plantvision.unl.edu/dataset公开获得。

该数据集包含玉米、高粱、珍珠谷子、普罗米谷子等谷类作物的图像，每棵植株均有0°和90° 2个侧视图。使用内布拉斯加大学林肯分校的Lemnatec Scanalyzer 3D高通量表型分析设备对植株进行拍摄，该设备的可见光摄像头每天采集2个侧视图，本研究选用数据集中的3组玉米植株图像，共计180张图像。该数据集包含原始图像、地面实况图像和带注释的图像，以促进叶片检测和跟踪、叶片分割、整体和组分表型分析领域的研究[13]。

1.2 植株图像裁剪与二值化

UNL-CPPD包含玉米植株从0～27 d生长过程中逐天采集的图像，通常为了便于分析，考虑采用0°侧视图为样本，确保植株完整的展现，在该特定视角下，相机的主光轴垂直植株所在的竖直面。所选取的图像如图1所示。本研究所使用UNL-CPPD中的图像，其图像采集装置经过合理设计，视野清晰，距离适中。观察采集箱的背景板为白色，两侧黑色粗带为采集箱的支撑框架。

图1 UNL-CPPD数据集中的玉米植株图像Fig.1 Image of maize plant from UNL-CPPD data base

借鉴张卫正等[14]的方法，从图像的中轴线分别向两侧逐步检测分析，当检测到出现采集箱的黑色支撑架时，分别将左右两侧的黑色支撑架裁剪掉，如图2所示。由于摄像头与采集箱位置固定，为了便于玉米植株图像的分析和参数提取，统一裁减掉培养盆的上边沿到底部区域(根据先验知识，裁减掉底部的共计260行像素)，如图3所示,从而消除培养盆的干扰。

图2 裁剪掉两侧玉米植株图像黑色支撑架Fig.2 Cut the black support frame of maizeplant image on both sides

采用常用平均值法，将PNG格式的图像转换为灰度图像，即计算PNG图像(与RGB图像都具备3个通道)的3个通道的平均值作为灰度。

本研究采用最大类间方差法(Ostu)算法确定灰度阈值然后进行图像二值化处理。该算法按图像的灰度特性，将图像分成背景和前景两部分。背景和前景之间的类间方差越大，说明构成图像的两部分的差别越大，当部分前景错分为背景或部分背景错分为前景都会导致两部分差别变小。因此，使类间方差最大的分割意味着错分概率最小。Ostu算法不受图像亮度和对比度的影响，计算简单、稳定有效，是实际应用中经常采用的确定阈值的方法[15-16]。

设灰度图像灰度级是L，则灰度范围为[0,L-1]，利用Ostu算法计算图像的最佳阈值为：

t=Max[w0(t)×(u0(t)-u)^2+w1(t)×(u1(t)-u)^2)]

(1)

其中，当分割的阈值为t时，w0为背景比例，u0为背景均值，w1为前景比例，u1为前景均值，u为整幅图像的均值。

使以上表达式值最大的t，即为分割图像的最佳阈值。对原图像f(x,y)使用上面的准则找到阈值t，将图像分割为2个部分，从而实现图像二值化，如图4所示。

图4 玉米植株的二值化图像Fig.4 Binary image of maize plant

1.3 骨架提取

玉米植株的二值化图像中，不可避免会产生噪声及面积较小的背景等干扰，由先验知识发现，整棵玉米植株就是一个连通区域，且该区域的面积最大。采用Matlab中的函数bwlabeln对玉米植株的二值化图像进行连通区域标记，经检测得到3个连通区域，查找并保留面积最大的连通区域，将该连通区域保存为玉米植株，如图5所示。删除剩余连通区域，从而消除噪声及背景干扰。

图5 玉米植株图像中最大的连通区域 Fig.5 The largest connected area in the maize plant image

由于玉米植株的叶片及茎秆颜色变化，植株个别像素点的颜色与背景接近，在图像处理的过程中会产生孔洞。采用Matlab的孔洞填充函数imfill()予以消除。

图像中的目标细化是获得目标骨架的过程，骨架细化算法的目的就是为了方便描述、抽取目标的线性特征[17-18]。一般细化算法要满足以下几点要求：(1)保证原来图像边缘的连通性；(2)原来图像边缘的中心线为骨架图像；(3)尽可能多地将1个像素宽度的线条图像作为细化结果；(4)细化速度要快[19]。本研究采用的骨架细化函数为bwmorph(image,‘thin’,Inf)，其中，image表示被处理的目标图像，‘thin’表示细化操作，Inf表示进行多次的操作，直到图像不能再细化。通过细化获得的骨架模型，如图6所示。

由于叶片较宽或叶片受到外力、病害等因素而出现折弯，经过细化的玉米植株的骨架可能会出现毛刺，本研究设定阈值为30，对骨架的毛刺予以去除。

图6 玉米植株的细化及骨架模型Fig.6 Thinning and skeleton model of maize plant

1.4 二叉树分析

1.4.1 骨架的端点和交叉点查找 分别以构成植株骨架的坐标点为中心，查找该中心的8邻域。当该点的领域有3个值不为0，则该点为交叉点；当该点的领域有1个值不为0，则该点为端点；通过Matlab中的函数bwmorph(image,‘endpoints’)和bwmorph(image,‘branchpoints’)查找并标记骨架线的端点和交叉点，如图7所示。

图7 玉米植株骨架的端点和交叉点Fig.7 End points and branch points of the maize plant skeleton

1.4.2 骨架的二叉树分析 二叉树是每个结点最多有2个子树的树结构。通常子树被称作“左子树”(Left subtree)和“右子树”(Right subtree)[20]。二叉树常被用于实现二叉查找树和二叉堆[21]。本研究通过分析植株的骨架结构，发现可以采用二叉树对骨架进行分析，逐步实现对植株骨架的标记和识别，步骤如下：

步骤1：选取骨架的最低点(即行坐标最大的端点)为起点S，沿骨架线向上查找，检测到的第1个交叉点即为一级节点(Node 1)，从S走到Node 1之间的骨架线即为一级茎段。

步骤2：以Node 1为起点，该点具有3个邻域，其中一级茎段是Node 1的邻域之一，该邻域已经完成标记；分别向Node 1的剩余2个邻域(这2个邻域可称之为2条支路)行走。当沿某条支路查找，遇到端点，则该条支路为叶片1，如果该端点位于Node 1的左侧，则称该条支路为左侧叶片1；如果该端点位于Node 1的右侧，则称该条支路为右侧叶片1。当沿某条支路查找，遇到交叉点，则该交叉点为二级节点(Node 2)，该条支路为二级茎段。

步骤3：由于交叉点有多个，如果该交叉点不是最后一个，则重复步骤2；如果该交叉点是最后一个，则该交叉点分别到最后剩余的2个端点，构成2条支路，分别根据端点相对交叉点的位置，称最后剩余的两条支路为左侧叶片n-1、右侧叶片n。

应用本研究所提出的以上算法对图7中的骨架进行标记和分析，如图8所示。

图8 玉米植株骨架的标记和分析Fig.8 Skeleton marking and analysis of maize plant

针对图7中玉米植株的骨架，其对应的二叉树分析如图9所示。可以看出，玉米植株的结构完全符合二叉树的特征，可以应用二叉树进行分析。图8与图9可实现一一对应，完全满足骨架分析的要求。

图9 玉米植株的二叉树分析Fig.9 Binary tree analysis of maize plant

2 结果与分析

根据以上茎叶标记与分析的结果，将玉米植株骨架和原始图像中的茎叶分割效果展示出来，各个叶片与茎分别用随机颜色显示，如图10所示。

采用本研究所提出的方法对其他2组图像进行骨架提取和端点、交叉点检测，如图11所示，对应的茎叶分割效果如图12所示。

图10 骨架和图像中的茎叶分割效果Fig.10 Stem and leaf segmentation of skeleton and image

图11 样本的骨架提取和端点、交叉点检测Fig.11 Skeleton extraction and end points,branch points detection of samples

图12 样本的茎叶分割效果Fig.12 Stem and leaf segmentation of samples

该图像处理方法可在Thinkpad电脑[产品型号T460P，16G RAM，Intel(R) 双核I7-6700 HQ CPU @2.60 GHz]中进行分析，每幅图像的处理用时小于8 s，可满足科研和生产的要求。

由于数据库中的图像，均是在具备充分的先验知识的情况下所采集，放置玉米植株的采集箱及承载玉米植株的培养盆的大小已知，摄像头与玉米植株样本的距离为常量，可由此计算样本图像中每个像素所代表的尺寸(以mm为单位)，从而进一步计算玉米的株高、叶片长度和茎宽等参数。此外，本研究中的所有图像处理操作均在Matlab软件中开发完成，自动化程度高，具有复现性，可对同一植株的不同生长阶段进行茎叶分割。

3 结论与讨论

骨架提取算法直接影响茎叶分割的效果，本研究中所采用的细化算法快速地将玉米植株的骨架细化到最简状态，并剔除毛刺，保证了植株的结构特征。宗泽等[22]以玉米作物为例，结合最小二乘法和遗传算法提取作物骨架等表型参数，株高误差减小了35%，叶长误差减小了70%，叶倾角误差减小了20%，有效地提高了作物表型参数的测量精度。对比上述研究，本研究中骨架提取算法计算量小，精度满足研究要求。

本研究消除噪声及其他干扰，保留玉米植株的图像，图像错分率可忽略不计，并创新性地提出二叉树分析算法，进行植株骨架的茎叶分割，从下到上依次对植株骨架进行遍历分析，分割的效果具有准确、直观的优点，设置随机色显示不同叶片，便于对分割的精度进行人工评定。玉米植株的茎叶分割所需时间小于8 s，满足实时性的要求。王传宇等[23]提出了一种基于时间序列多图像的玉米植株干旱胁迫表型检测方法，采用统计学习方法对图像序列进行初始背景分割，图像错分率为2.12%，在植株器官分割阶段，融合植株形态先验知识，通过逐行扫描和邻域数量判断拆分植株器官像素，并计算叶片角度比值。对比上述研究，本研究在玉米植株的茎叶分割精度、处理速度及效果展示等方面都具有较大的进步。

本研究从已有的内布拉斯加大学林肯分校的植物表型数据集选取玉米植株图像，该数据集采用标准化种植，图像的采集时间及方式经过精心的设计与规划，相对于现有的玉米植株栽培，具有批量化和集约化的优势。自动进行玉米植株图像的提取和二值化，通过建立植株的骨架模型及去除骨架中的毛刺，实现骨架上的端点和交叉点的检测，采用二叉树分析骨架模型，确立各级节点及各个叶片，从而实现茎叶分割。通过多个样本图像的试验，表明本研究所采用骨架提取和二叉树分析的方法可以实现玉米植株的茎叶分割，具有快捷有效、计算速度快、自动化程度高等优点，满足植物表型提取的科研要求，为玉米等植物表型分析与育种等提供支持，具有一定的实用性。

后续可在本研究的基础上进行茎叶角测量、叶片长度测量、叶片曲率计算、叶片存在交叉情况的叶片分割与识别等研究，为植株活力检测、植物对环境胁迫的响应等整体和组分表型分析提供支持。通过调整相关参数，本研究所提出的方法可扩展到其他具有相似结构的谷物作物(例如高粱)植株的检测与分析，具有一定的推广价值。