1 引言

图像分割是图像处理到图像分析的一个关键步骤, 其目的为将一幅图像划分为若干个具有不同特性且有意义的区域, 每个区域具有相似的特征. 图像分割算法大致可分为基于阈值的分割^[1]、基于聚类的分割^[2]、基于区域的分割^[3]和基于图论的分割^[4]等. 其中, 基于聚类的分割方法以其快速、高效的特点, 越来越广泛的用于遥感图像、医学图像和自然图像的分割.

模糊C均值聚类(FCM)^[5]及其改进算法是各种聚类技术中使用最广泛的方法. FCM算法为每一个样本的归属引入了模糊性, 保留了原始图像更多的信息. 但是该算法存在3个缺点: (1)只考虑像素灰度忽视了像素间的空间信息, 算法缺乏抗噪性和鲁棒性; (2)没有充分利用不同特征的优势; (3)使用欧氏距离作为距离度量.

为了提高传统FCM算法的抗噪性和鲁棒性, Ahmed等^[6]在FCM目标函数中增加了空间信息约束, 提出了FCM_S算法, 由于每次迭代时都要计算邻域, 因此该算法时间复杂度较高. 为了降低FCM_S算法的复杂度, Chen等^[7]提出FCM_S1和FCM_S2算法, 由于均值/中值滤波图像可以提前计算, 因此在算法迭代前就获得了邻域信息, 降低了计算成本. Zhong等^[8]利用熵的概念提出了自适应结合空间信息的AFCM_S1算法, 同时解决了手动调参问题. Lei等^[9]提出了结合形态学重建与隶属度滤波的FRFCM分割算法, 该算法利用形态学重建平滑原图像, 提高了算法的抗噪性和细节保护能力.

为了充分提取图像特征, 一些学者提出了融合图像多特征的算法. Belongie等^[10]提出了一种基于颜色和纹理特征的图像分割算法, 该算法将三种颜色特征和三种纹理特征嵌入到单个特征向量中. Yu等^[11]提出了融合多特征的图像分割算法, 该算法通过亮度相似度、纹理相似度和边缘相似度来定义区域差异性. 尽管这些算法得到了可靠的分割结果, 但不同特征之间的权重调节仍存在问题. Rajaby等^[12]提出的WHIFCM算法使用图像的色相和强度分量, 并通过自适应调整的权重将其组合到模糊C均值目标函数中.

很多分割算法大都是使用欧氏距离计算图像中目标的相似性, 该距离是一种线性度量. 近年, 研究者发现使用非线性度量可以更多地发现图像中的分布信息. Gong等^[13]提出了KWFLICM算法. 该算法使用像素的空间距离和灰度值来重新定义权重因子, 并引入核诱导距离进行距离度量, 实现了更好的分割效果. Shang等^[14]提出CKS_FCM算法, 该算法将欧氏距离度量替换为核诱导距离度量, 提高了图像分割的精确度.

基于上述3个问题, 本文提出了结合形态学重建和超像素的多特征FCM分割算法(SRMFCM). 首先使用形态学闭合重建和Mean-Shift超像素分割方法^[15]预处理原始图像; 然后提取重建图像各像素的颜色、纹理和梯度特征, 利用平均策略定义各超像素的颜色、纹理和梯度特征. 最后, 运用区域代表像素点, 用核诱导距离进行距离度量, 实现多特征加权聚类, 合并具有相同标签的区域来获得最终分割结果.

2 相关工作 2.1 形态学重建

形态学重建能在不知道噪声类型的情况下较好的去除噪声并保留物体轮廓^[9]. 最基本的形态学腐蚀重建 ${R^\varepsilon }$ 和形态学膨胀重建 ${R^\delta }$ 如式(1)所示.

$\left\{ {\begin{array}{*{20}{c}} {R_f^\varepsilon \left( g \right) = \varepsilon _f^{\left( i \right)}\left( g \right),g \ge f}\\ {R_f^\delta \left( g \right) = \delta _f^{\left( i \right)}\left( g \right),g \le f} \end{array}} \right. $

(1)

其中, $\varepsilon $ 表示腐蚀操作, $\delta $ 表示膨胀操作. $f$ 表示原始图像, $g$ 表示标记图像. 式(1)满足的条件包括 $\varepsilon _f^{\left( 1 \right)}\left( g \right) = \varepsilon \left( g \right) \vee f$ , $\varepsilon _g^{\left( i \right)}\left( f \right) = \varepsilon \left( {{\varepsilon ^{\left( {i - 1} \right)}}\left( g \right)} \right) \vee f$ , $\delta _f^{\left( 1 \right)}\left( g \right) = \delta \left( g \right) \wedge f$ , $\delta _g^{\left( i \right)}\left( f \right) = \delta \left( {{\delta ^{\left( {i - 1} \right)}}\left( g \right)} \right) \wedge f$ , 其中 $ \vee $ 和 $ \wedge $ 分别代表逐点的最大值和最小值. 在实际运用中, $g = \varepsilon \left( f \right)$ 满足膨胀重建中的条件 $ g\le f $ , $g = \delta \left( f \right)$ 满足腐蚀重建中的条件 $ g\ge f $ . 因此, 为了简洁有效, 总是利用 $g = \varepsilon \left( f \right)$ 和 $g = \delta \left( f \right)$ 来形成标记图像.

${R^\varepsilon }$ 和 ${R^\delta }$ 可以结合成有较强滤波能力的形态学开启重建算子 ${R^O}$ 和形态学闭合重建算子 ${R^C}$ . 因为 ${R^C}$ 更适于平滑图像纹理细节, 本文选择形态学闭合重建 ${R^C}$ 来处理原始图像^[9]. ${R^C}$ 如式(2)所示, 其中B表示腐蚀和膨胀操作所需的结构元.

${R^C}\left( f \right) = R_{R_f^\delta \left( {\varepsilon \left( f \right)} \right)}^\varepsilon \left( {\delta_B \left( {R_f^\delta \left( {\varepsilon_B \left( f \right)} \right)} \right)} \right)$

(2)

2.2 超像素分割

超像素分割利用图像的局部相似性, 将图像分割成若干个具有相似颜色、亮度和纹理等特征的图像块. 相比于以像素为单位的图像分割方法, 以超像素为单位的图像分割方法更有利于提取图像局部特征来获取更有效的图像信息, 同时降低计算复杂度. 超像素分割的基本原理是将一幅大小为 $M \times N$ 的图像P划分为n个非空区域 $\left( {{R_1},{R_2}, \cdots ,{R_n}} \right)$ , 这些区域满足 ${R_1} \cup {R_2} \cup \cdots \cup {R_n} = P$ 且当 $i \ne j$ 时, ${R_i} \cap {R_j} = \emptyset $ . 现有的超像素分割方法有基于图论的方法和基于梯度下降的方法. 本文采用的基于梯度下降的Mean-Shift算法, 该算法生成的超像素块不规则, 边缘贴合度好, 具有很好的抗噪性.

2.3 EWFCM

EWFCM^[16]是最大熵正则化的加权模糊C-均值算法, 该算法通过最小化类内离散度同时最大化属性权重熵来获得最佳的聚类结果. 给定样本集合 $X = \left\{ {{x_1},} \right.$ $\left. {x_2}, \cdots ,{{x_n}} \right\} \subset {R^m}$ , m是样本空间的维数, n是样本个数. EWFCM算法通过最小化式(3)将样本X分为c个类.

$\begin{split} {F_{\rm EWFCM}} = &\sum\limits_{k = 1}^c {\sum\limits_{i = 1}^n {u_{ki}^\alpha } } \sum\limits_{j = 1}^m {{w_{kj}}} {\left( {{x_{ij}} - {c_{kj}}} \right)^2} \\ &+\lambda \sum\limits_{k = 1}^c {\sum\limits_{j = 1}^m {{w_{kj}}} } \log\left( {{w_{kj}}} \right) \end{split}$

(3)

其中, ${u_{ki}}$ 表示第i个样本属于第k类的程度, 满足条件 $\displaystyle\sum\nolimits_{k = 1}^c {{u_{ki}}} = 1$ , $0 \le {u_{ki}} \le 1$ . $\alpha $ 表示模糊因子, 一般取2. ${c_{kj}}$ 表示第k类中第j维特征的聚类中心. ${w_{kj}}$ 表示第k类中第j维特征所占的权重, 满足条件 $\displaystyle\sum\nolimits_{j = 1}^m {{w_{kj}}} = 1$ , $0 \le {w_{kj}} \le 1$ . ${\left( {{x_{ij}} - {c_{kj}}} \right)^2}$ 表示样本点 ${x_{ij}}$ 与中心点 ${c_{kj}}$ 的欧氏距离. 式(3)第一部分加权的欧氏距离度量用来控制集群的形状和大小, 促进了集群的聚集. 第二部分是属性权值的负熵, 它根据现有数据对所有属性权值的最优分布进行了正则化. $\lambda \left( {\lambda > 0} \right)$ 是正则化可调参数. 只要选择合适的 $\lambda $ , 就能平衡式(3)的两部分, 得到最优稳定解.

3 本文算法 3.1 图像预处理

首先, 利用形态学闭合重建处理原始图像P, 得到重建图像 $P'$ ; 然后, 按照Mean-Shift算法将重建图像 $P'$ 预分割成n块区域, 写成集合形式为 $R = \left\{ {{R_1},} \right.\left. {{R_2}, \cdots ,{R_n}} \right\}$ , 每个区域 ${R_i}\left( {1 \le i \le n} \right)$ 的总像素数表示为 $\left| {{R_i}} \right|$ .

3.2 特征提取

HSV颜色模型是一种与人类视觉感知很相似的颜色模型, 本文分别从RGB空间和HSV空间中提取图像 $P'$ 的颜色特征向量. 在RGB空间中, 图像 $P'$ 的每个像素 ${P_{x,y}}$ 都会有一个3维颜色特征向量 $\left[ {P_{x,y}^r,P_{x,y}^g,P_{x,y}^b} \right]$ . 同理, 在HSV空间中, 图像 $P'$ 每个像素 ${P_{x,y}}$ 也有一个3维颜色特征向量 $\left[ {P_{x,y}^h,P_{x,y}^s,P_{x,y}^v} \right]$ . 每个区域中所有像素在各颜色通道上的均值作为本区域的颜色特征. 例如, 区域R_i在R颜色通道上的颜色特征如式(4)所示.

$R_i^r = \frac{{\displaystyle\mathop \sum \nolimits_{\left( {x,y} \right) \in {R_i}} P_{x,y}^r}}{{\left| {{R_i}} \right|}}$

(4)

同理可得该区域在其他颜色通道上的颜色特征向量.

Gabor是用于边缘提取的线性滤波器, 它具有良好的方向选择和尺度选择特性, 对光照变化不敏感, 因此十分适合纹理分析. 本文选择1个尺度和8个方向的Gabor滤波器来提取图像的纹理特征. 首先将图像 $P'$ 转为灰度图像 $P_{\rm gray}^{'}$ , 然后利用Gabor滤波器对 $P_{\rm gray}^{'}$ 进行卷积, 则每个像素 ${P_{x,y}}$ 都得到8维的纹理特征向量 $\left[ {P_{x,y}^{t1},P_{x,y}^{t2},P_{x,y}^{t3},P_{x,y}^{t4},P_{x,y}^{t5},P_{x,y}^{t6},P_{x,y}^{t7},P_{x,y}^{t8}} \right]$ , 每个区域中所有像素在各方向上的纹理特征均值作为本区域的纹理特征, 例如, 区域R_i在方向t1上的纹理特征如式(5)所示.

$R_i^{t1} = \frac{{\displaystyle\mathop \sum \nolimits_{\left( {x,y} \right) \in {R_i}} P_{x,y}^{t1}}}{{\left| {{R_i}} \right|}}$

(5)

同理可得该区域在其它方向上的纹理特征向量.

经典的图像梯度算法是考虑图像的每个像素的某个邻域内的灰度变化, 利用边缘临近的一阶或二阶导数变化规律, 对原始图像中像素某个邻域设置梯度算子, 通常运用小区域模板进行卷积来计算. 通常有Sobel算子、Robinson算子、Laplace算子等. 本文利用Sobel算子生成梯度图像, Sobel算子有两个, 一个用于检测水平边缘, 另一个用于检测垂直边缘. 如式(6)和式(7)所示. ${G_x}$ 和 ${G_y}$ 分别代表图像 $P_{\rm gray}^{'}$ 经横向及纵向边缘检测的图像灰度值.

${G_x} = \left[ {\begin{array}{*{20}{c}} { - 1}&0&1 \\ { - 2}&0&2 \\ { - 1}&0&1 \end{array}} \right]*P_{\rm gray}^{'}$

(6)

${G_y} = \left[ {\begin{array}{*{20}{c}} { - 1}&{ - 2}&{ - 1} \\ 0&0&0 \\ 1&2&1 \end{array}} \right]*P_{\rm gray}^{'}$

(7)

图像 $P_{\rm gray}^{'}$ 的梯度图像 ${P^G}$ 如式(8)所示.

${P^G} = \sqrt {G_x^2 + G_y^2} $

(8)

区域R_i中所有像素的梯度均值作为本区域的梯度特征向量, 如式(9)所示.

$R_i^G = \frac{{\displaystyle\mathop \sum \nolimits_{\left( {x,y} \right) \in {R_i}} P_{x,y}^G}}{{\left| {{R_i}} \right|}}$

(9)

3.3 目标函数

原图像经过预处理后得到n块区域, 写成集合形式为 $R = \left\{ {{R_1},} \right.\left. {{R_2}, \cdots ,{R_n}} \right\}$ . 每个区域R_i都有15维的特征向量, 分别是6维颜色特征向量、8维纹理特征向量和1维梯度特征向量, 所以整体的数据集规模为 $n \times 15$ . 本文最初的目标函数如式(10)所示.

$ \begin{split} F =& \sum\limits_{k = 1}^c {\sum\limits_{i = 1}^n {u_{ki}^\alpha } } \sum\limits_{j = 1}^m {{w_{kj}}} {\left\| {\Phi \left( {{x_{ij}}} \right) - \Phi \left( {{c_{kj}}} \right)} \right\|^{\rm{2}}} \\ &+\lambda \sum\limits_{k = 1}^c {\sum\limits_{j = 1}^m {{w_{kj}}} } log\left( {{w_{kj}}} \right) \end{split}$

(10)

其中, ${\left\| {\Phi \left( {{x_{ij}}} \right) - \Phi \left( {{c_{kj}}} \right)} \right\|^{\rm{2}}}$ 是利用核诱导距离表示的归一化量度. $ \Phi (·)$ 表示非线性映射, 在特征空间的内积表示为 $\Phi {\left( {{x_{ij}}} \right)^{\rm{T}}}\Phi \left( {{c_{kj}}} \right) = K\left( {{x_{ij}},{c_{kj}}} \right)$ . 核替换的过程如式(11)所示.

$ \begin{split} &{\left\| {\Phi \left( {{x_{ij}}} \right) - \Phi \left( {{c_{kj}}} \right)} \right\|^{\rm{2}}}\\ = & {\left( {\Phi \left( {{x_{ij}}} \right) - \Phi \left( {{c_{kj}}} \right)} \right)^{\rm{T}}}\left( {\Phi \left( {{x_{ij}}} \right) - \Phi \left( {{c_{kj}}} \right)} \right)\\ = & \Phi {\left( {{x_{ij}}} \right)^{\rm{T}}}\Phi \left( {{x_{ij}}} \right) - \Phi {\left( {{x_{ij}}} \right)^{\rm{T}}}\Phi \left( {{c_{kj}}} \right)\\ &- \Phi {\left( {{c_{kj}}} \right)^{\rm{T}}}\Phi \left( {{x_{ij}}} \right) + \Phi {\left( {{c_{kj}}} \right)^{\rm{T}}}\Phi \left( {{c_{kj}}} \right)\\ =& K\left( {{x_{ij}},{x_{ij}}} \right) - 2K\left( {{x_{ij}},{c_{kj}}} \right) + K\left( {{c_{kj}},{c_{kj}}} \right)\\ =& 2 - 2K\left( {{x_{ij}},{c_{kj}}} \right) \end{split}$

(11)

最终的目标函数如式(12)所示.

$F \!\!=\!\! \mathop \sum \limits_{k = 1}^c \mathop \sum \limits_{i = 1}^n u_{ki}^\alpha \mathop \sum \limits_{j = 1}^m {w_{kj}}\left( {1 - K\left( {{x_{ij}},{c_{kj}}} \right)} \right)\!\! +\!\! \lambda \mathop \sum \limits_{k = 1}^c \mathop \sum \limits_{j = 1}^m {w_{kj}}\log\left( {{w_{kj}}} \right)$

(12)

其中,

$K\left( {{x_{ij}},{c_{kj}}} \right) = \exp\left( {\frac{{ - \|{x_{ij}} - {c_{kj}}\|^2}}{\sigma }} \right)$

(13)

式(13)表示高斯核函数, σ表示函数的宽带参数, 且满足k(x,x)=1. 使用拉格朗日乘子法最小化式(12), 整理得 ${u_{ki}},{c_{kj}},{w_{kj}}$ 的迭代公式如式(14)~式(16)所示.

${u_{ki}} = \frac{1}{{\displaystyle\mathop \sum \nolimits_{h = 1}^c {{\left( {\frac{{\displaystyle\mathop \sum \nolimits_{j = 1}^m {w_{kj}}\left( {1 - K\left( {{x_{ij}},{c_{kj}}} \right)} \right)}}{{\displaystyle\mathop \sum \nolimits_{j = 1}^m {w_{hj}}\left( {1 - K\left( {{x_{ij}},{c_{hj}}} \right)} \right)}}} \right)}^{\textstyle\frac{1}{{\alpha - 1}}}}}}$

(14)

${c_{kj}} = \frac{{\displaystyle\mathop \sum \nolimits_{i = 1}^n \left( {u_{ki}^\alpha K\left( {{x_{ij}},{c_{kj}}} \right){x_{ij}}} \right)}}{{\displaystyle\mathop \sum \nolimits_{i = 1}^n \left( {u_{ki}^\alpha K\left( {{x_{ij}},{c_{kj}}} \right)} \right)}}$

(15)

${w_{kj}} = \frac{{\exp\left( { - {\lambda ^{ - 1}}\displaystyle\mathop \sum \nolimits_{i = 1}^n u_{ki}^\alpha \left( {1 - K\left( {{x_{ij}},{c_{kj}}} \right)} \right)} \right)}}{{\displaystyle\mathop \sum \nolimits_{s = 1}^m \exp\left( { - {\lambda ^{ - 1}}\displaystyle\mathop \sum \nolimits_{i = 1}^n u_{ki}^\alpha \left( {1 - K\left( {{x_{is}},{c_{ks}}} \right)} \right)} \right)}}$

(16)

3.4 算法流程

本文提出的SRMFCM算法流程如算法1.

算法1. SRMFCM算法

输入: 原始图像P, 聚类数目C, 模糊因子m, 参数 $\scriptstyle \lambda$ , 最大迭代次数max_iter和高斯核kernel.

输出: 图像分割结果Segm_P.

Step 1.根据式(1)和式(2)对图像P进行形态学闭合重建, 生成重建图像P'.

Step 2. 利用Mean-Shift方法将图像P'分割成N块区域, 区域代表像素点完成后续聚类.

Step 3. 根据式(4)提取各区域的颜色特征. 根据式(5)提取各区域的纹理特征. 根据式(6)~式(9) 提取各区域的梯度特征. 将所有特征向量整合为一个特征向量矩阵.

Step 4. 初始化隶属度矩阵 $\scriptstyle {U^0}$ , 设置循环计算器t=0.

Step 5. 用式(15)更新聚类中心; 用式(14)更新隶属度矩阵; 用式(16)更新权重矩阵; 用式(12)计算目标函数F; t++.

Step 6. 如果t>max_iter, 转向Step 7, 否则返回Step 5.

Step 7. 根据隶属度矩阵得到最终分割结果Segm_P.

4 实验结果与分析

为了验证本文所提SRMFCM算法的有效性, 本节测试了BSDS300数据集的6幅自然图像. 本文选择的对比方法有4种, 第1种结合形态学重建和隶属度滤波的FRFCM算法; 第2种是基于像素点的用熵来调节类内紧凑度和像素空间信息的AFCM_S1算法; 第3种是基于核诱导距离度量的KFCM算法; 第4种是基于区域的融合图像颜色和纹理特征的FCM算法, 简写为FCM(RC+RT).

图1中Image1-Image6表示6幅测试图像, 每幅图像的尺寸大小都为321×481像素. Image1-Image6的主要分割目标分别为海星, 马, 老虎, 鹰, 蜻蜓和草坪. 本实验过程中, 只需要分割出上述主要目标物体即可, 其他部分都作为背景区域对待, 因此聚类数目都设为2.

图2–图7展示了使用上述5种算法对这6幅自然图像的分割结果, 其中GT表示真值图像.

首先分析图2, 由于FRFCM算法、AFCM_S1算法和KFCM算法都是以像素为单位, 所以海星缺失了很多像素, 虽然FCM(RC+RT)算法和SRMFCM算法的分割结果优于前三者, 但是FCM(RC+RT)算法分割的海星中仍缺少一些像素, 使用SRMFCM算法分割的结果中仍存在部分噪声.

对于图3, 草地影响了AFCM_S1算法、KFCM算法和FCM(RC+RT)算法对马的分割效果. 虽然FRFCM算法中草地的影响较小, 但是分割的目标缺少部分像素. 相比于前4种算法, SRMFCM算法分割结果比较精确.

对于图4, FRFCM算法的分割结果较差. 其它4种算法都较准确的分割出了老虎, 但使用SRMFCM算法分割的结果噪声最少.

对于图5, 由于原图中的天空颜色由内向外逐渐加深, 所以影响了AFCM_S1算法和KFCM算法的分割效果. FRFCM算法、FCM(RC+RT)算法和SRMFCM算法的分割结果较好, 但是相比于FRFCM算法和FCM(RC+RT)算法, SRMFCM算法对于鹰的尾巴有更好的细节保留.

对于图6, 算法都可以将蜻蜓的主体分割出来, 对蜻蜓的足部等细节也有很好的保留, 但使用SRMFCM算法得到的分割结果最好.

对于图7, AFCM_S1算法和KFCM算法分割的草坪中有很多噪声, 而使用FRFCM算法、FCM(RC+RT)算法和SRMFCM算法得到的结果优于前两种算法, 且SRMFCM算法得到的分割结果最好.

其次, 利用信息检索指标(F-measure, F)和错误率(Error Rate, ER)这两种经典的测试指标评价不同算法的分割结果. 这两种指标如式(17), 式(18)所示.

$ ER=\frac{FP+FN}{TP+FP+FN+TN} $

(17)

$F = \frac{{\left( {1 + a} \right) \times Precision \times Recall}}{{a \times Precision + Recall}}$

(18)

其中, 系数a一般1, 精确率 $Precision = TP/\left( {TP + FP} \right)$ , 召回率 ${{Re}} call = TP/\left( {TP + FN} \right)$ . TP表示划分正确的目标像素数, FP表示误分到目标中的像素数, FN表示目标区中的点错分到背景区的像素数, TN表示划分正确的背景像素点数. 由式(17)和式(18)可知, F值(ER值)越高(越低), 图像分割效果越好.

表1和表2分别为使用5种算法对6幅图像分割结果的F值和ER值. 由表中数据可知, SRMFCM算法在分割精度上优于FGFCM算法、AFCM_S1算法、KFCM算法和FCM(RC+RT)算法. 因此, 本文提出的结合形态学重建和超像素的多特征FCM分割算法比其它基于单特征或者多特征FCM算法有更好的性能.

5 结论与展望

为了解决FCM算法存在的3个问题, 本文提出了结合形态学重建和超像素的多特征FCM算法(SRMFCM). 该算法实现了图像多特征的有机结合, 提高了图像分割的精度, 弥补了以往FCM算法的不足. 利用BSDS300数据集中的6幅自然图像进行实验对比, 结果证明该算法有更高的分割精度. 本文的不足之处在于对比实验不够完善. 因此, 优化对比实验是接下来研究的重点.