基于水平集和视觉显著性的植物病害叶片图像分割

周强强，王志成，赵卫东，陈宇飞

（1.同济大学 电子与信息工程学院，上海201804；2.江西农业大学 计算机与信息工程学院，南昌330027）

植物病害传染会引起严重经济损失并影响农业生产和食品安全，因此检测、防治植物病害成了农业生产中一项非常重要的工作.但传统植物病害的识别诊断方法主要是靠专业技术人员的田间观察或实验室化验分析.这些方法不仅耗时费力，而且也会由于缺乏足够的经验导致判断的不准确.近些年来基于计算机视觉的植物病害自动检测和识别研究，已经证明在植物的大面积栽种监控方面是行之有效的，能在植物生长的早期阶段提供识别和治理植物病害的重要信息.实际上，植物病害的主要来源是叶片，比如叶斑病，稻瘟病、褐斑病、叶枯病等，大约80%～90%的植物病害在叶部［1］.因此国内外出现了许多基于叶片症状的植物病害图像识别方法研究.其中植物病害图像分割是对植物病害图像进行处理的基础，是实现植物病害情况自动检测和田间管理自动化的前提.而植物病害图像成分复杂，病叶上的病斑排列无规则，颜色深浅不一，且存在一定的随机噪声，因此分割问题也是植物病害图像研究分析的一个重点.文献［2］中，通过先将RGB红绿蓝色彩模式的病害叶片图像进行特定颜色空间的转换，然后根据直方图中的亮度局部最大化通过阈值分割找出病理位置.文献［3］中开发了水稻病害诊断系统，该系统使用局部熵阈值法和最大类间方差法（OTSU）等图像分割方法进行相关目标的分割.文献［4］借助水平集方法对叶子病斑图像的RGB分量图像颜色信息取加权值，以差分图像能量作为能量函数最终值，以适应不同的病害种类.此外还有许多基于人工神经网络，支持向量机和线性判别分析的研究应用到植物叶片病害的分类和识别中［5-6］.

已有的大多数分割方法应用于植物的叶片病害图像大多存在边缘模糊或边缘不连续等缺点，给后续的特征提取带来不便，并且由于一些叶片图像背景的复杂性和非均匀性，初始参数设置过多等原因，也会影响上述方法的效果.本文所提出的基于水平集和视觉显著性的植物病害叶片的彩色图像分割方法可以较好地解决这种问题，并且由于引入了视觉注意机制，还提高了基于水平集分割算法的自适应性和运行效率.

1 基于视觉注意力机制的初始轮廓

活动轮廓模型在演化中的初始轮廓选取往往是从整个图像的边缘开始迭代进行的，这样易使轮廓演化的结果受到背景信息的干扰，导致无法得到满意的目标结果.而且大多数活动轮廓模型仅仅只考虑了图像灰度信息，当图像中存在噪声或目标被部分遮挡等情形时，仅靠这些模型是很难使得轮廓收敛在目标物体上，这时如果还能够加入目标的颜色、纹理、形状等方面的先验知识，可使轮廓在这些先验知识的约束下准确地演化到目标物体的周围.

而视觉显著性检测是一种模拟生物视觉注意机制的选择性注意模型，其利用图像的颜色、灰度、对比度和方向等信息模拟人类视觉体系，一般用于自然场景图像或视频中感兴趣区域的定位、预测和转移.正好可以借此来构建轮廓先验.于是考虑引入图像的显著性分析，可以预判出目标物体所在的位置即物体的大致轮廓，这样既避免了背景信息对目标的干扰，也使得轮廓演化的次数得以大幅减少，从而提高了传统活动轮廓模型的效果和效率.

本文采用基于小波变换使用底层特征的显著性检测算法，相比经典的Itti模型和SR（谱剩余）模型，同样具有计算量较小的优点，还具有小波变换的独有优势是它能够同时提供多尺度上的空域和频域分析，可以将图像分解为更低分辨率水平上的低频轮廓信息和原始信号在水平、垂直和对角线方向的高频细节信息，且可以对图像作多次分解多级子带信号［7］.

基本算法描述如下：

步骤1 将RGB的颜色图像转化成为更符合人眼感知特性的CIE Lab颜色空间图像，同时应用一个m×m的2D高斯低通滤波器消除高频噪声，即：式中：g表示滤波器；Ic为输入的彩色图像；Igc为输出图像；“＊”是卷积运算；m为模板大小值，一般取3.

步骤2 对图像通过小波变换（wavelet transform，WT）（这里采用Daubechies小波）形成多层级的子带；

式中：N是小波分解过程的最大尺度数，尺度指标s∈｛1，…，N｝，c是对应图像Igc的所属颜色通道，c∈｛L，a，b｝.KcN是每个通道上最大尺度的逼近输出，Hcs，Vcs，Dcs分别是给定c和s图像上的水平，垂直和对角细节上的小波系数.

步骤3 逆小波变换创建特征图.

步骤4 计算局部特征的全局分布以得到全局显著图.

由式（1）可知，图像（x，y）位置处可以依据所有的特征图表示成一个3×N大小的特征向量v（x，y）.一个给定位置处的特征概率可以用一个符合正态分布的概率密度函数来定义：

式中：Σ=E［（f（x，y）-μ）（f（x，y）-μ）T］是一个n×n的协方差矩阵；n=3×N；|Σ|表示行列式；μ是均值向量，μ=E［v］.

结合式（2）就可计算得全局显著图sG：

这里gk×k是一个k=5的高斯低通滤波，以得到一个平滑显著图.

步骤5 局部显著图sL（x，y）的计算.通过线性融合（x，y）处每个尺度上的特征图得到局部显著图如下：

步骤6 联合局部和全局的显著图.

式中：s′（x，y）是最终显著图，s′L（x，y）和s′G（x，y）分别是归一化到［0，1］的局部和全局显著图；是调制函数.接着进一步对最终显著图s′（x，y）进行基于注意力焦点（focus of attention，FOA）的优化处理：

以0.8为阈值，将最终显著图中显著性值大于0.8的像素点视为 FOA，即式（4）中的（x′，y′）表示位置处，其显著性值为s′（x′，y′），越靠近FOA的区域显著性越强，越远离的区域显著性越弱.于是将图像中的非FOA点的显著值按式（4）进行更新，其中的dFOA（x，y）表示位置（x，y）与其最近的 FOA 点（x′，y′）之间的距离.这样处理后的显著图更符合人的视觉机制.

步骤7 采用自适应阈值的显著区域分割方法，对生成的显著图S（I）进行二值化，划分为显著区域Os和非显著区域ONS，分别代表目标区域和背景区域.

步骤8 形态学处理消除一些小区域或弱连通部分，提取出一个输入图像I的初始轮廓C.

这样得到初始轮廓C一般都是包含目标物体的不规则曲线，是符合人类的视觉感知特性的，具有一定的自适应性.所以也可以作为输入图像的先验形状约束，引入到几何活动轮廓模型中，很好地利用了图像的视觉特性.

图1分别是褐斑病叶片基于小波变换的视觉显著性检测所得到局部显著图，全局显著图和最终显著图.图2是通过对图1生成的最终显著图进行二值化分割后得到的初始轮廓.

2 向量值图像的边界检测

用来做病害检测的植物叶片图像大都是采集的RGB彩色图像，这也是一种向量值图像.植物叶片彩色图像分割的存在难题在于：难以检测目标的真实边界，背景的复杂性，目标和背景的颜色是不同质的，通常基于梯度的边缘描述子并不能表示彩色图像的边界特征.

在后面构造的基于水平集的活动轮廓演化模型中，除了利用到区域信息外，边缘信息也加入到能量函数的构建中，可以减少分割结果对光照，阴影和噪音的敏感性.而要检测向量值图像的边缘，往往需要计算图像的梯度值.在当前的边缘检测方法中，形态梯度法通常是一个简单有效的方法，该方法通过膨胀和腐蚀等形态运算来计算像素邻域上最大梯度值与最小梯度值之间的差值.但形态学运算要求的是标量或灰度图像的输入，这就成了向量值图像处理的一个障碍.如果按照传统做法，将向量图像转换成标量图像后进行分割，或对每个通道进行分割然后再将结果进行融合.这样做既不能充分利用图像信息，又会损失可利用的信息，更忽略了通道间的相关性，导致无法检测完整目标，对分割精度也会受到影响.

图1 褐斑病叶片的显著性检测Fig.1 Saliency detection of brown spot leaf

图2 初始轮廓Fig.2 Initial contour

而Cumani［8］考虑到人眼的视觉感知特点，将一维的像素值空间扩展成多维向量，从而克服了以上的问题.一个m通道的图像可以用函数f：R2→Rm来表示，图层上的点P（x1，x2）映射为一个m维的向量f=（f1（x1，x2），f2（x1，x2），…，fm（x1，x2）），单色图像是m=1时的特例.若彩色图像用RGB颜色模型表示，那么该映射可以转化为，f：（x1，x2）→（fR（x1，x2），fG（x1，x2），fB（x1，x2）），其中fR（x1，x2），fG（x1，x2），fB（x1，x2）分别是对应的点（x1，x2）处的RGB分量.考虑图像中两个临近点P和Q处值的差分：

当Q-P是一个无穷小的位移dP=（dx1，dx2）时，此差分变成为微分：

再次求微分得：

这里的”·”表示Rm中的向量内积.令x1=x，x2=y，则可以表示如下：

其中的M是结构张量矩阵，式（5）的极值对应M的特征值.并且通过求M的特征值和特征向量，可以得到图像边缘幅度和相位.矩阵M的特征值λ±和特征向量n±分别如下：

λ+和λ-分别描述图像中的最大和最小变化率.在灰度图像中，边缘的强度是最大梯度的大小.但是在向量值图像中，该强度是关于λ+和λ-间的一个比较函数，向量值图的边缘图一般可用：fedge=f（λ+，λ-）=λ+-λ-来描述.

然后再考虑运用基于多尺度形态学梯度的方法来重建图像边缘信息，计算图像的梯度值▽B（fedge），从而得到所要的边缘截止函数.

多尺度形态学梯度应用于图像分割，可以对阶跃边缘和模糊边缘进行有效处理，能够在大尺度下抑制噪声，可靠地识别边缘，在小尺度下定位，再由粗到细跟踪边缘，得到边缘的位置.在这定义多尺度梯度为

式中⊕和Θ分别表示膨胀和腐蚀运算，B为结构元素.形态学梯度的性能关键是B的选取，而多尺度形态学梯度可以结合利用大结构元素和小结构元素的各自优点.这里设Bi（0≤i≤n）为一组正方形的结构元素，Bi的大小为（2i+1）×（2i-1）像素.

在以上得到的结果中每个像素的值表现了它与其邻域的对比强度，这对于强调边缘上的梯度峰值很重要.同时，为了增强显著区域像素的效应，加速轮廓朝目标边缘的演化，于是可设计如下的边缘停止函数：

式中的S（I）为得到的显著图.图3是式（7）生成的褐斑病叶片多尺度形态学梯度图.

图3 多尺度形态学梯度图Fig.3 Multi-scales morphology gradient map

3 基于水平集的活动轮廓演化

3.1 基于显著区域的活动轮廓模型

在第1节显著性检测的基础上，这里提出一种基于区域显著性的活动轮廓模型，首先根据显著性检测算法得到目标边界的形状近似，即初始轮廓C，然后将这种先验形状信息加入活动轮廓模型中，构建一个图像分割能量泛函，再采用水平集方法实现曲线的演化和图像分割.

水平集［9］方法的核心思想是将移动界面C（t）作为零水平集嵌入到更高维的闭超水平集函数中，记为ø（x，y，t），由闭超曲面的演化过程可以得到水平集函数的演化方程，而嵌入的闭超曲面总是其零水平集，最终只需要确定水平集函数的零水平集即可确定移动界面的演化结果.通常记ø=0为零水平集，ø（x，y，）为水平集函数.初始水平集函数通常设定为初始闭合曲线C0生成的符号距离函数，即ø（x，y，0）=±d（x，y），其中d（x，y）表示点（x，y）到曲线C0的距离，点（x，y）位于闭合曲线内部取值为正，否则为负.在演化方程的控制下不断更新曲线在固定坐标系下的位置，最终使其逼近目标物体的边缘.

对于一个给定彩色图像u0，被闭合的初始轮廓C划分为显著区域ΩS和非显著性区域ΩNS两个区域，它们可分别看作是目标区域和背景区域的先验形状近似.这两个区域的通道图像在Ω内的像素均值用ci+（或ci-）表示，i表示彩色图像的R、G或B三个通道.于是定义先验形状适应能量函数为

式（9）中的两个能量项分别是初始轮廓C内部和外部区域的灰度值与c+i和c-i的平方误差，即实际图像与假定的“先验形状近似”图像之间的偏离.当初始轮廓经过演化运动后的曲线与目标边界不符合时，式（9）中的能量函数达不到最小，只有当曲线运动到目标边界时，能量函数才能达到最小值.此外，λi+和λi-

对应不同通道的系数，可以调节演化曲线内部和目标平均灰度差别，以及曲线外部和背景平均灰度差别对曲线形变的影响程度.根据先前的研究［10］，在任意时刻，水平集曲线C将整幅图像分成内外两个部分，分别计算C内部区域灰度熵Ein和外部区域灰度熵Eout：

式中：N为通道图像的灰度级，Pi（øn）为曲线演化第n次迭代中第i级灰度出现的概率.式中定义的熵从某种程度上反映了图像灰度的多样性，灰度相对丰富时熵较大，而灰度相对单一时熵较小.曲线内外部图像灰度熵的变化将有助于引导曲线的演化过程，使之更快地收敛于理想目标区域.实际上，可以通过实时计算曲线内外部区域图像灰度熵以获得拟合参数和）.

当然，除了可以通过熵值来确定λ+i和λ-i的合理比例，还可以根据熵在形变过程中的变化来判断曲线的演化是否趋于稳定.

3.2 改进的距离正则化水平集演化模型

距离正则化水平集演化（distance regularized level set evolution，DRLSE）模型是一种边缘信息的变分水平集分割方法.该方法在传统水平集方法中添加了一个距离规则化项，解决了重新初始化的问题，提高了水平集的演化速度.在文献［11］中，其能量泛函如下：

令ø：Ω→R是定义在Ω域上的水平集函数（level set function，LSF）.能量泛函ε（ø）：

式中：Rp（ø）是水平集规则化项；μ＞0是常数；εext（ø）为外部能量项.Rp（ø）定义为：

式中：dX是一个二重积分，p是势能（或能量密度）函数p：［0，∞）→R.加入规则化项的目的不仅能平滑水平集ø，而且在零水平集附近能保持符号距离属性|▽ø|=1.当p=p1（s）≜ø（s-1）2，其中s=1，p有唯一的最小值点.根据梯度下降流的公式：∂ø／∂t=μdiv（dp（|▽ø|）▽ø）.其中dp（·）表示为函数p求导.dp（p）=p'（s）／s=1-（1／s），该方程是一个扩散方程，扩散率为D=ξdp（▽ø），ξ为系数.当|▽ø|→0时，扩散率的值会-∞，影响计算的精确性.为解决该问题，Li等［11］引入了一个具有2个最小值点的函数p2（s）：

式（13）能保持在零水平集附近符号距离属性|▽ø|=1，而在远离零水平集附近保持LSF为常数，同时|▽ø|=0.为了保持LSF的轮廓，势能函数在s=0和s=1处有最小值.这样的函数也叫双阱势能函数.

式（11）中函数Rp（ø）的 Gateaux导数为∂Rp／∂ø=-div（dp（|▽ø|）▽ø），其中div（·）为散度算子，dp是一个函数，定义为.ξR（ø）的p

式中，ξ，λ，α分别为各能量项的系数，且λ＞0.δ为Dirac函数，H为 Heaviside函数，δ（ø）=H'（ø）.在梯度下降流为

于是，一个包含了外部能量项εext（ø）的总能量为：实际应用中，H（ø）和δ（ø）由平滑函数Hε（ø）和δε（ø）代替.同时，本文使用式（8）中定义的边缘图gcolor来代替通常基于图像梯度的边界停止函数g.

但如果仅仅是这样的话，DRLSE方法的演化曲线对初始轮廓的位置会很敏感，并且也没有包含任何的区域信息，会给分割的效果带来较大的影响.为了让式（14）中的DRLSE能更好地控制轮廓的演化，同时包含区域信息，再将前面基于显著性检测的活动轮廓模型如式（9）所示，也整合进来，那么对于一个水平集函数ø，定义一个新的能量泛函为

用梯度下降法求解得到的对应变分水平集公式为

其中的平滑函数Hε（ø）和δε（ø）分别为

参数ε根据实际情况设置.

式（15）中c+i和c-i由如下公式迭代更新：

图4所示的即为本文方法采用式（15）所得到的褐斑病叶片的分割结果.

4 实验结果分析

本文选用Matlab2014a编程实现本文提出的图像分割算法，运行环境是在 MacBook Pro OS X 10.9.5的虚拟机Paralles配置上（Intel（R）Core（TM）i5-4258UCPU＠2.40GHz，内存4.00GB）安装的Windows7 64位操作系统.

图4 本文方法的褐斑病叶片分割Fig.4 Segmentation of proposal method

实验数据库来源于某农作物病害防治研究机构，以实际采集到的叶稻瘟病，黄瓜霜霉病两种病害叶片为测试对象进行实验验证.其中每种病害各选取50幅不同的叶片进行分割实验，下面每个示例图像就是选自这50幅之一，且都具有不同程度的复杂背景.从这些图中可以看出，不同病害叶片的病斑所表现的症状不同，叶稻瘟病早期显褐点型病斑，黄瓜霜霉病呈黄色不规则形霉斑.这些差异特征是利用叶片识别病害的主要依据.

图5 叶稻瘟病叶片分割Fig.5 Segmentation of blast leaf

用来进行比较的算法是k-Means聚类彩色图像分割算法和C-V的向量图像分割算法［12］.从图5—6可以看出，对叶稻瘟病进行分割时，直接用C-V的向量图像分割算法时，贴近目标和覆盖目标时的初始轮廓，都不能有效提取病斑，且算法的迭代次数要达500～800次才趋于稳定；即使用k-Means算法，提取出来的病斑也不是太理想；而利用本文设计的方法分割效果相对很好，基于显著性检测的改进距离正则化水平集活动轮廓模型的分割迭代次数只需要25次.图7—8对黄瓜霜霉病图像进行分割时，三种方法效果相差不大，其中k-Means算法的效果最好，C-V的向量分割算法其次，但对初始轮廓的设置很敏感，且运算迭代次数要300次以上，如图7—8所示.而本文的算法虽然这时效果较逊些，但从表1中，可知算法的运行时间最少，效率最高，运算迭代次数也只需要30次.

图6 叶稻瘟病叶片的C-V向量值图像分割Fig.6 C-V segmentation of vector blast leaf image

表1是在本实验所采用的100幅抽样病害图像上运行算法得到的，充分反映了基于显著性自适应初始轮廓分割算法的高效率.

图7 黄瓜霜霉病叶片分割Fig.7 Segmentation of downy mildew cucumber leaf

图8 黄瓜霜霉病叶片的C-V向量值图像分割Fig.8 C-V segmentation of downy mildew cucumber leaf

表1 运行时间比较Tab.1 Comparison of Runtime

5 结论

本文在距离正则化的水平集演化过程结合了图像的显著性分析，这样，通过除了引入基于多尺度形态学梯度边缘检测算子，以保证DRLSE能更好地控制轮廓的演化外，同时还包含了区域信息，保证了水平集函数至少在零水平集附近理想的符号距离函数的性质，解决了DRLSE方法的演化曲线对初始轮廓位置的敏感问题.仿真实验表明，本文分割算法具有演化速度快，能够减少某些情况下复杂背景信息对目标的干扰，且具有一定的自适应性，提高了传统活动轮廓模型的效果和效率.

［1］ Zhang S，Zhang C.Orthogonal locally discriminant projection for classification of plant leaf diseases［C］／／Computational Intelligence and Security.Leshan：IEEE，2013：241-245.

［2］ Camargo A，Smith J S.An image-processing based algorithm to automatically identify plant disease visual symptoms［J］.Biosystems Engineering，2009，102（1）：9.

［3］ Kurniawati N N，Abdullah S N H S，Abdullah S.Investigation on image processing techniques for diagnosing paddy diseases［C］／／Soft Computing and Pattern Recognition.Malacca：IEEE，2009：272-277.

［4］ 胡秋霞，田杰，何东健，等.基于改进型C-V模型的植物病斑图像分割［J］.农业机械学报，2012，43（5）：157.HU QiuXia，TIAN Jie，HE Dongjian，etal.Segmentation of Plant Lesion Image Using Improved C-V model ［J］.Transactions of the Chinese Society for Agricultural Machinery，2012，43（5）：157.

［5］ Wang X，Zhang M，Zhu J，etal.Spectral prediction of Phytophthora infestans infection on tomatoes using artificial neural network［J］.International Journal of Remote Sensing，2008，29（6）：1693.

［6］ Camargo A，Smith J.An image processing based algorithm to automatically identify plant disease visual symptoms［J］.Biosystems Engineering，2009，102（1）：9.

［7］ Imamoglu N，Lin W，Fang Y.A saliency detection model using low-level features based on wavelet transform［J］.Multimedia，IEEE Transactions on，2013，15（1）：96.

［8］ Cumani A.Edge detection in multi-spectral images［J］.Graphical Models and Image Processing，1991，53（1）：40.

［9］ Osher S，Sethian J.Fronts propagating with curvaturedependent speed： Algorithms based on Hamilton Jacobi formulations［J］.Journal of Computational Physics，1988，79：12.

［10］ 陈宇飞，吴启迪，赵卫东，等.基于图像熵的快速Chan-Vese模型分割算法［J］.同济大学学报：自然科学版，2011，39（5）：738.CHEN Yufei，WU Qidi，ZHAO Weidong，etal.Fast chan-vese segmentation algorithm based on image entropy［J］.Journal of Tongji University：Natural Science，2011，39（5）：738.

［11］ LI Chunming，XU Chenyang，GUI Changfeng，etal.Distance regularized level set evolution and its application to image segmentation［J］.IEEE Transactions on Image Processing，2010，19（12）：3243.

［12］ Chan T F，Sandberg B Y，Vese L A.Active contours without edges for vector-valued images ［J］.Journal of Visual Communication and Image Representation，2000，11（2）：130.