基于扇环形区域图像分割的小麦秸秆截面参数测量方法

徐胜勇　彭程里　陈　可　王令强　任喜峰　段宏兵

(1.华中农业大学工学院， 武汉 430070； 2.农业部长江中下游农业装备重点实验室， 武汉 430070；.华中农业大学植物科学技术学院， 武汉 430070)

0　引言

作物秸秆是一种重要的可供开发利用的生物质资源，其综合利用对稳定农业生态平衡、促进农民增产增收、缓解能源与环境压力具有重要作用[1]。秸秆理化结构特性不仅与倒伏性能及其综合利用有关，而且其维管束系统担负输送水分、无机盐和有机养料的重要作用，是影响产量、品质和养分利用效率等重要的农艺性状[2]。作物秸秆的力学特性与其微观结构密切相关。王庭杰等[3]认为影响玉米茎秆抗压强度的主要因素为厚壁/半径、机械组织比例、维管束个数等。刘唐兴等[1]分析的结果表明茎秆的木质部、皮层和维管束对茎秆抗倒性有重要影响。陈桂华等[4]研究结果表明，水稻单茎抗推力和茎粗、大维管束数目及小维管束数目呈显著正相关。赵春花等[5]分析牧草茎秆力学性能与显微结构的相互关系，发现牧草茎秆的强度和刚度主要取决于各自机械组织的厚度、维管束的数量以及各组织及其细胞之间的连接形式和连接强度。ZHANG等[6]对番茄植物茎的微观结构进行观察和分析，发现其微结构主要由表皮、皮质、维管束和髓腔等组成，并测试分析了茎部不同部位的水分、拉伸和剪切机械性能。姚金保等[7]发现小麦秸秆的力学特性和抗倒性与其机械组织细胞层数和厚度、维管束数目和面积以及髓腔大小密切相关。对秸秆微观结构的精确测量，可为深入分析作物秸秆力学性能与显微结构的相互关系，进而指导强抗倒性品种的选育以及设计开发新机具、新工艺提供理论依据。

当前，在农业领域使用显微成像与分析技术研究微观结构已成为趋势，常用的成像技术包括：扫描电子显微镜(SEM)图像[8-9]；核磁共振成像(MRI)[10]；X射线微计算机断层扫描(CT)[11-13]；光学显微镜成像[14-16]。

光学显微镜和数码成像装置构成的显微镜成像系统，操作方便、直观、检定效率高，同时适用范围非常广泛，配以测量软件可以构建高通量的显微图像分析系统[17-22]。光学显微图像分析系统的研究为农业领域的微观结构观测开启了新途径。

人工测量作物秸秆微观结构的方式费时费力，严重制约了相关研究的深入开展，急需高通量、高精度的自动化分析系统。徐胜勇等[2]研究了一种小麦秸秆截面显微结构测量系统，可以精确地测量高质量小麦茎秆切面显微图像中的厚壁厚度、髓腔大小和整体尺寸，并可以识别维管束的数量。该测量系统证明了运用显微图像处理技术对小麦秸秆的微观性状参数进行分析和测量是高效可行的方法。秸秆标本显微图像的采集主要包括切片、染色和拍照3个环节，任何一个环节的失误都将导致显微图像质量下降。针对文献[2]的显微图像分析系统测量精度不高，尤其是维管束的识别较难。本文基于厚壁和维管束在小区域的灰度特征，研究一种高鲁棒性的图像分析方法，用于精确分割厚壁和维管束并测量一系列关键参数。

1　材料和方法

1.1　小麦秸秆截面的显微图像采集

2017年5月于华中农业大学实验田采集15个品种的成熟期小麦秸秆，每个品种取10根秸秆，截取第1节、第2节，每节手工制作2个切片(厚度约0.02 mm)，共计600个切片样本。先使用5%浓度的对钾苯酚-乙醇溶液浸泡1 min、再使用95%的浓盐酸浸泡1 min，对切片进行染色。清洗染色切片并晾干，置于载玻片上，放置在光学显微镜下观察和拍照，图片保存为JPG格式。显微镜型号为Nikon SMZ800，相机型号为Nikon DS-Ri 1，图像分辨率为1 024像素×1 280像素。根据显微图像中的比例尺可以将像素距离转换为物理距离，本文中2个相邻像素点的物理距离为1.25 μm。本文算法开发和测试的平台为通用计算机(i5 4200M 2.5 GHz、8 GB DDR3内存、250 GB三星固态硬盘、NVIDIA GeForce GT 755M显卡)、Windows 7专业版SP1和Matlab 2016a。

1.2　扇环形分区域图像分割方法

成熟期小麦秸秆质地坚韧、含水率各异，制作标本切片时容易破坏截面的组织结构，染色中药剂的浓度、浸泡时间、对焦和光源都可能导致标本图像质量的下降。图1a为高质量的显微图像，维管束的内部结构清晰可见；图1b为低质量的显微图像，刀片造成了一些划痕，截面的左右两部分染色不均匀，维管束内部结构不可见且曝光过度。文献[2]的方法不能测量图1b之类的低质量显微图像，需要研究鲁棒性更强的图像处理算法。

图1　小麦秸秆截面显微图像Fig.1　Wheat stalk micrographs

染色是秸秆中的木质素与化学物质发生化学反应的过程，木质素含量越高，该区域被染色得越深。在小麦秸秆中，厚壁和维管束的木质素含量较高，薄壁组织也含有少量的木质素。在染色的秸秆截面显微图像中，厚壁和维管束的整体亮度显著低于相邻的其他组织。因此，依据显微图像的灰度特性进行阈值分割分离出厚壁和维管束区域是可行的。但如果使用一个全局阈值，低质量的显微图像中某些薄壁组织可能被错误地分割为厚壁或维管束。因此，有必要按照小区间的灰度特征进行阈值分割，以得到精确的分割结果。

扇环形分区域图像分割方法的整体思路是：将茎秆截面按照一定的规则划分为很多小区域，再对各小区域单独进行阈值分割，最后对分割结果进一步处理，分离出完整的厚壁和维管束。如图2所示，小麦秸秆截面的内外侧轮廓构成了一个近似的环形区域A。先粗略定位厚壁，依据定位信息从A中划分出包含有厚壁的环形区域A1；再粗略定位维管束，划分出紧邻A1内侧的包含有维管束的环形区域A2。以圆心为端点，对A的外侧圆周画射线，射线将A1和A2分割为很多个扇环形区域。其中，厚壁所在的A1区域被均匀划分为n份；对于维管束所在的A2区域，先对维管束进行粗略定位，再依据定位信息将A2区域划分为m份，每个子区域包含1个或2个完整的维管束。秸秆截面被划分为大量的扇环形区域之后，对每个区域独立进行阈值分割，再进行拼接等后续处理，得到完整而精确的厚壁和维管束，继而以此为基础精确测量截面的各种参数。

图2　小麦秸秆截面的扇环形区域划分示意图Fig.2　Sector ring region-based division of wheat stalk micrograph

1.3　厚壁的检测

厚壁定位与分割的算法流程和结果如图3所示。对小麦秸秆截面的彩色显微图像进行灰度化，得到Ig。再使用Otsu算法进行阈值分割，得到一幅包含全部染色组织的二值图像Ib。将该二值图像反色，再使用形态学填充算法填充空洞，得到小麦秸秆截面的整体图像Iwh。统计Iwh中白色像素点个数Swh，即为截面的像素面积。小麦秸秆截面的外轮廓近似于圆形，可以估计小麦秸秆截面的整体像素宽度H，计算公式为

(1)

对多个品种大量小麦秸秆样本的厚壁进行测量，统计结果显示厚壁的厚度均小于0.045H，即厚壁分布于秸秆截面外侧宽度0.045H的圆环区间。根据这个特性，Iwh减去尺寸为0.045H的圆形结构算子腐蚀一次的Iwh，余下的部分即为厚壁分布的圆环区域A1。将A1划分为n个扇环(n为可以被360整除的整数)，在灰度图像Ig中对每个扇环形小区域独立进行阈值分割，合并结果得到一幅包含厚壁和少量组织粘连的二值图像。用尺寸为0.045H的圆形结构算子对该二值图像进行形态学“开”运算以断开粘连，然后检测和保留最大连通分量，得到精确的厚壁区域。此时，厚壁中还存在有少量的空洞和断裂，用尺寸为0.01H的圆形结构算子进行形态学“闭”运算，最终得到分割效果良好的厚壁二值图像Ith。统计该二值图像中灰度为255的像素点的数量，即为厚壁面积Sth。厚壁为近似的圆环，厚度Hth的计算公式为

(2)

如果不使用扇环形分区域阈值分割方法，则由于厚壁灰度分布的不均匀，将导致较多的分割错误，如图3中小图所示。

1.4　维管束的检测

1.4.1维管束分布的圆环区域定位

维管束的数量和面积是秸秆截面参数测量中最重要而复杂的问题。本文提出一种迭代腐蚀法定位维管束分布的圆环区域A2，再进行后续处理得到维管束区域。如图4所示，秸秆截面整体图像Iwh减去厚壁图像Ith，得到近似圆形区域B，该区域中不包含厚壁。将包含全部染色组织的二值图像Ib反色，再使用尺寸为0.005H的圆形结构元素先后进行形态学“闭”运算和“开”运算，以填补圆环中的孔洞并消除小面积的噪声，得到二值图像Ic，其中的圆环形区域C包含有厚壁、维管束和薄壁组织。统计该图像中灰度255的像素点数量，即为小麦秸秆截面的面积Sb。维管束就分布于B的外轮廓到C的内轮廓之间的环形区域。定义一个尺寸为0.005H的圆形结构元素对B进行t次腐蚀，直到腐蚀结果与C的交集为空集，此时对B腐蚀的结果为Bt。从B中减去Bt得到圆环区域，即为维管束分布的环形区域A2。

图3　厚壁的检测流程与结果Fig.3　Procedures and results of sclerenchyma detection

图4　迭代腐蚀法定位维管束分布的圆环形区域Fig.4　Ring region of all vascular bundles distribution located via iterative-erosion method

1.4.2扇环形小区域的划分

图5　维管束所属区域的扇环形区域划分Fig.5　Ring regions segmentation of vascular bundle

扇环形小区域划分的流程如图5所示。首先对原灰度图像进行尺寸0.005H×0.005H的中值滤波以消除孤立的小区域噪声，再进行尺寸0.02H×0.02H的中值滤波以消除薄壁组织形成的网状条纹并消除维管束内部的小孔洞，得到灰度图像Inw。将均分为k个扇环形区域的A2作为掩膜作用于Inw，计算每个小区域的灰度均值。以逆时针顺序统计的区域灰度均值为纵坐标，以区域的序号k为横坐标，可以构建关于扇环形分割的扇环区域序号-灰度均值坐标系(图5b)。依次将2个相邻的极大灰度均值区域之间的区域合并成一个新的区域，最终A2被划分为m个新区域(m为灰度均值极大值的个数)，每个新区域将包含1个或2个维管束的全部或绝大部分。

1.4.3扇环形区域阈值分割和结果处理

图6　维管束区域图像分割过程与结果Fig.6　Procedures and results of image segmentation of vascular bundle region

扇环形区域A2中包含有维管束、薄壁组织和高亮度的背景。如图6所示，将划分为m个扇环形区域的A2作为掩膜作用于Inw，对每个小区域都使用Otsu算法进行图像分割，将秸秆截面组织区域从高亮度的背景中分离出来；再次使用Otsu算法对该秸秆截面组织区域进行图像分割，分离出维管束。合并m个扇环形区域的阈值分割结果，得到以维管束为主体、包含噪声的二值图像。使用半径为0.01H的圆形结构算子进行“开”运算，消除图像中的小块噪声，并断开维管束上的粘连。再对每个连通分量独立使用半径为0.01H的圆形结构算子进行“闭”运算，以消除维管束中的孔洞和断裂；再用一个半径为0.02H的圆形结构算子做开运算消除毛刺。统计全部连通域的面积平均值，剔除小于该均值20%的连通域，保留的结果即为完整的维管束。最后统计连通区域个数，即为维管束数量K；统计全部连通域的像素点数量，即为维管束的总面积Sk。相比均匀分区域阈值分割的方法，本文的分割方法具有更高的精度和稳定性，如图6中小图所示，分割数m的取值对于分割结果的影响非常显著，很容易造成误分割。

2　实验结果与分析

2.1　算法性能测试实验

为评价本文算法的有效性和正确性，测量了50幅小麦秸秆显微图像(高质量30幅、中等质量15幅和低质量5幅)。以人工测量结果为真实值，将本文算法测量结果、文献[2]的测量结果与其进行比较。比较了截面直径H、截面面积Sb、厚壁厚度Hth、厚壁面积Sth、维管束面积Sk、维管束数量K。本文算法的典型参数为：n=180、k=90，文献[2]典型参数为：模板式样为1、缩放级数为7、旋转角度为2°。人工测量方法为：将一幅显微图像复制多份，使用Windows画图工具软件打开。根据测量目标，人工用画笔分别在复制的图像中标注出秸秆截面的外轮廓、内轮廓、厚壁区域和维管束区域等，删除目标区域之外的背景，最后统计目标区域的参数作为真实值。使用均方根误差(RMSE)和平均相对误差(MRE)来描述测量精度，R2是决定系数。显微图像处理中的度量单位是像素，表1中，截面直径、截面面积、厚壁厚度、厚壁面积、维管束面积的单位是像素，维管束数量的单位是个。实际应用中需要转换为物理尺寸，单位为微米(μm)或平方微米(μm2)，本文实验得到的尺寸转换因子分别为1.25 μm/pixel和1.56 μm2/pixel。

本文测量方法结果：截面直径的RMSE和MRE分别为6.9像素和0.7%；截面面积的RMSE和MRE分别为11 220.6像素和3.6%；厚壁厚度为0.40像素和1.6%；厚壁面积为1 580.2像素和2.5%；维管束面积为2 873.7像素和5.4%；维管束数量的均方根误差为0.39个，平均相对误差为1.1%。实验结果表明，本文算法具有很高的测量精度和稳定性。相比文献[2]，每一种参数的测量精度都提高了。尤其是厚壁厚度的检测误差从9.9%降低到了1.6%，这归因于分区域阈值分割方法减少了秸秆切片阴影引起的误差。维管束数量的识别率从93.2%提高到了98.9%。分析实验结果发现，文献[2]的方法对于高质量的显微图像中维管束的检测率达到100%，但图像质量越差，其检测到的维管束数量越少。这是因为低质量图像中，维管束的“三眼”结构已完全被破坏，模板匹配的方法已不能正常运行。本文算法对于维管束数量的检测也存在1.1%的误检率，当相邻的维管束之间相连的薄壁组织灰度非常近似于维管束时，2个维管束会被检测为1个，造成检测数量比实际数量少；当染色非常不均匀导致部分薄壁组织的灰度值近似于维管束的时，薄壁组织被错误判断为维管束，造成检测数量比实际数量多。

表1　测量结果精度Tab.1　Measurement result of two methods

2.2　参数设置对算法性能的影响

秸秆截面的扇环形分区域图像分割方法中，扇环形区域的划分方式对算法的准确度和速度有直接的影响。随机选取了100幅显微图像，以最具代表性的厚壁面积、维管束面积和数量这3个待检参数为例，测试了参数设置的变化对算法性能的影响。定义图像分辨率为1 024像素×1 280像素、n=180、k=90时测试图像的算法执行时间平均值为基准时间。当其他参数改变时，算法执行时间均值与基准值的百分比称为耗时比。以下实验结果均是在仅改变某一个参数、其他参数不变的条件下获得的。

2.2.1n对厚壁面积测量的影响

当n取值为1～360之间可以整除15的整数时，统计本文算法的正确率和耗时比(基准时间为61 s)，结果如图7a所示。随着n的增大，正确率持续增大，耗时比与n取值呈现线性增加关系。当n取值大于90时，正确率稳定在98%以上。由此可见，环形区域A1的划分数量越大，灰度不均匀的影响越小，同时耗时会线性增加。

图7　参数设置对算法性能的影响Fig.7　Effect of parameters setting on performance of algorithm

2.2.2α对维管束数量和面积测量的影响

根据实验统计，一个维管束所占角度为6°～8°。k的取值从36增加到360，相应地扇环形区域的角度α从10°递减到1°(7°因没有合适的k值，被排除)，统计α在该区间每隔1°的正确率和耗时比(基准时间为49 s)，结果如图7b所示。对于维管束数量，当α取值在2°～8°之间时，检测正确率均超过98%；当α大于9°时，正确率迅速下降。对于维管束面积，正确率基本随着α取值的增大而减小；当α大于6°时，正确率迅速下降。耗时比与α呈近似的线性关系，拟合方程为y=-17.73x+180.52，R2=0.891 6。

3　结束语

针对小麦秸秆中各组织的结构和颜色特征，设计了针对性的图像处理算法，实现了对截面直径、截面面积、厚壁厚度、厚壁面积、维管束面积和维管束数量等关键参数的精确测量。使用该算法测量一幅图像平均耗时约110 s，而技术人员人工测量一幅图像平均耗时超过600 s，且工作强度非常大。本文提出的扇环形分区域分割算法对于切片、染色和拍摄等环节造成的低质量显微图像具有很好的处理效果。在厚壁的分割中，如需更高的精度，只需提高划分区域的数量，同时耗时也线性增加。在维管束的分割中，不能单纯地依靠提高划分区域的数量实现检测精度的提升，当划分的扇环区域所占角度接近并小于维管束所占角度时，可以保持高正确率；过多或过少的划分区域数量，都会导致维管束的粗定位信息不准确，以致于导致后续处理的错误。

1刘唐兴，官春云，黎移新. 甘蓝型油菜主茎显微结构与抗倒性关系的初步研究[J]. 中国农学通报,2011,27(5): 139-143.

LIU Tangxing,GUAN Chunyun,LI Yixin. Prelimilary study on the relationship between stem microstructure and lodging resistance in rapeseed (BrassicanapusL.)[J].Chinese Agricultural Science Bulletin,2011,27(5):139-143.(in Chinese)

2徐胜勇, 段宏兵, 李东臣,等. 小麦茎秆截面参数显微图像测量系统[J/OL]. 农业机械学报, 2017, 48(7): 46-52.http:∥www.j-csam.org/jcsam/ch/reader/view_abstract.aspx?flag=1&file_no=20170706&journal_id=jcsam.DOI:10.6041/j.issn.1000-1298.2017.07.006.

XU Shengyong,DUAN Hongbing, LI Dongchen, et al．Measurement system for parameters of wheat stem section based on microimage processing[J/OL].Transactions of the Chinese Society for Agricultural Machinery, 2017,48(7): 46-52.(in Chinese)

3王庭杰，张亮，韩琼，等.玉米茎秆细胞壁和组织构建对抗压强度的影响[J]. 植物科学学报,2015,33(1): 109-115.

WANG Tingjie，ZHANG Liang，HAN Qiong，et al. Effects of stalk cell wall and tissue on the compressive strength of maize[J].Plant Science Journal,2015,33(1):109-115.(in Chinese)

4陈桂华，邓化冰，张桂莲，等. 水稻茎秆性状与抗倒性的关系及配合力分析[J]. 中国农业科学, 2016, 49(3): 407-417.

CHEN Guihua, DENG Huabing, ZHANG Guilian, et al. The correlation of stem characters and lodging resistance and combining ability analysis in rice[J]. Scientia Agricultura Sinica, 2016, 49(3): 407-417.(in Chinese)

5赵春花，韩正晟，师尚礼，等. 新育牧草茎秆收获期力学特性与显微结构[J]. 农业工程学报, 2011,27(8): 179-183.

ZHAO Chunhua, HAN Zhengsheng, SHI Shangli, et al. Mechanical properties and microstructure of new species forage stems in harvesting period[J]. Transactions of the CSAE, 2011, 27(8): 179-183. (in Chinese)

6ZHANG Xiliang, GUO Qiang, XU Yunfeng, et al. Mechanical testing of tomato plant stem in relation to structural composition[J]. Agricultural Research, 2016, 5(3): 236-245.

7姚金保，马鸿翔，姚国才，等. 小麦抗倒性研究进展[J]. 植物遗传资源学报, 2013, 14(2):208-213.

YAO Jinbao，MA Hongxiang，YAO Guocai，et al.Research progress on lodging resistance in wheat (TriticumaestivumL.) [J].Journal of Plant Genetic Resources, 2013, 14(2): 208-213.(in Chinese)

8KALIYAN N, MOREY R V. Natural binders and solid bridge type binding mechanisms in briquettes and pellets made from corn stover and switchgrass[J]. Bioresource Technology, 2010,101(3): 1082-1090.

9何春霞, 傅雷鸣, 熊静, 等. 不同表面处理对麦秸秆结构和性能的影响[J]. 南京农业大学学报, 2016,39(2): 325-331.

HE Chunxia, FU Leiming, XIONG Jing, et al. Study on the microstructure and properties of wheat straw with different treatments[J].Journal of Nanjing Agricultural University, 2016, 39(2): 325-331．(in Chinese)

10MUSSEA M, QUELLEC S, CAMBERT M, et al． Monitoring the postharvest ripening of tomato fruit using quantitative MRI and NMR relaxometry[J].Postharvest Biology and Technology, 2009,53(1-2): 22-35．

11GULATI T,DATTA A K. Coupled multiphase transport, large deformation and phase transition during rice puffing[J]. Chemical Engineering Science, 2016,139: 75-98.

12MAIRHOFER S, JOHNSON J, STURROCK C J, et al. Visual tracking for the recovery of multiple interacting plant root systems from X-ray μ CT images[J].Machine Vision and Applications, 2016, 27(5): 721-734.

13TING V J L, SILCOCK P, BREMER P J, et al. X-Ray micro-computer tomographic method to visualize the microstructure of different apple cultivars[J]. Journal of Food Science, 2013,78(11): 1735-1742.

14张哲, 赵静, 田津津,等. 冷冻-复温过程中葡萄细胞结晶变化研究[J/OL]. 农业机械学报, 2016,47(5): 211-217.http:∥www.j-csam.org/jcsam/ch/reader/view_abstract.aspx?flag=1&file_no=20160529&journal_id=jcsam.DOI:10.6041/j.issn.1000-1298.2016.05.029.

ZHANG Zhe, ZHAO Jing, TIAN Jinjin, et al． Research on crystallization change of grape cells during freezing-thawing process[J/OL].Transactions of the Chinese Society for Agricultural Machinery,2016,47(5): 211-217.(in Chinese)

15SIRICHAROEN P, SCOTNEY B, MORROW P, et al. Automated wheat disease classification under controlled and uncontrolled image acquisition[C]∥12th International Conference on Image Analysis and Recognition, ICIAR 2015,2015:456-464.

16BERNKLAU I, LUCAS L,BECKER T, et al. Protein network analysis—a new approach for quantifying wheat dough microstructure[J]. Food Research International, 2016,89: 812-819.

17PIECZYWEK P M, ZDUNEK A, UMEDA M. Study on parameterisation of plant tissue microstructure by confocal microscopy for finite elements modelling[J]. Computers and Electronics in Agriculture, 2011,78(1): 98-105.

18刘杏娥，胡飞，汪佑宏，等. 撑绿杂交竹主要解剖特征的研究[J]. 安徽农业大学学报, 2012,40(1): 1-4.

LIU Xing’e, HU Fei, WANG Youhong, et al. Main anatomical characteristics ofBamdusapervariadilis×Dendrocalamopsisdaii[J].Journal of Anhui Agricultural University, 2012,40(1): 1-4.(in Chinese)

19DEVAUX M F, BOUCHET B, LEGLAND D, et al. Macro-vision and grey level granulometry for quantification of tomato pericarp structure[J]. Postharvest Biology and Technology, 2008,47(2): 199-209.

20DHONDT S, HAERENBORGH D V, CAUWENBERGH C V, et al. Quantitative analysis of venation patterns of Arabidopsis leaves by supervised image analysis[J]. The Plant Journal, 2012,69(3): 553-563.

21PIECZYWEK P M, ZDUNEK A. Automatic classification of cells and intercellular spaces of apple tissue[J]. Computers and Electronics in Agriculture, 2012,81: 72-78.

22SIRAULT X R R, CONDON A G, WOOD J T, et al. “Rolled-upness”: phenotyping leaf rolling in cereals using computer vision and functional data analysis approaches[J]. Plant Methods, 2015,11(1): 1-11.