CN104021563B - 基于多目标模糊聚类和反对学习的噪声图像分割方法 - Google Patents

Abstract

本发明公开了一种基于多目标模糊聚类和反对学习的噪声图像分割方法，主要解决已有技术分割噪声图像效果差的问题。其实现步骤为：(1)提取噪声图像的邻域信息；(2)确定相邻像素间的约束关系；(3)生成初始种群，并进行初始设定；(4)构建初始种群中所有个体的图像细节保持和去除噪声两个目标函数及隶属度矩阵；(5)计算种群中所有个体的子问题目标函数；(6)对种群进行反对学习、交叉和变异操作，并进行更新；(7)对更新后的种群的Pareto前端，选取其拐点对应的个体的隶属度矩阵去模糊，输出分割结果。本发明既能保持图像细节，又能去除噪声，同时提高了图像分割效率，可用于提取后续图像处理所需的目标。

Description

基于多目标模糊聚类和反对学习的噪声图像分割方法

技术领域

本发明属于图像处理技术领域，特别涉及一种图像的分割方法，可用于提取后续图像处理如图像分类、图像语义识别等操作所需要的目标。

背景技术

图像分割在图像处理技术中占有很重要的地位，它是图像处理中的关键步骤。图像分割的目的是把图像分成若干个特定的、具有特殊性质的区域并提取出感兴趣的目标，这些被提取的目标可在后续图像处理中运用于图像语义识别等领域。

目前，常用图像分割方法主要有：阈值分割法，区域分割法，边缘分割法和基于特定理论分割法等。其中，基于特定理论分割法能够将各学科中一些特定的理论及方法和图像分割问题有效结合，因此在近几年得到广泛的应用和关注。在基于特定理论分割法中，模糊聚类法是一种将模糊集理论和聚类分析相结合的有效的图像分割方法。其中模糊集理论具有描述事物不确定性的能力，能够有效处理事物内在的不确定性，因此非常适用于图像分割问题；聚类分析则是通过将待分割图像像素映射到特征空间进行聚集，再将聚集结果映射回原图像空间而获得最终分割结果，因此能够有效分析图像像素之间的相似性。虽然传统的模糊聚类法在普通图像分割能够获得良好的分割结果，但其对噪声和灰度不均匀的情况敏感。因此在分割含有噪声的图像时，传统的模糊聚类法不能同时完成图像细节的保持和噪声的去除，所以不能获得有效的分割结果。

发明内容

本发明的目的在于针对上述已有技术在分割噪声图像时的不足，提出一种基于多目标模糊聚类和反对学习的噪声图像分割方法，以获得既能保持图像细节又能去除图像噪声的分割结果。

实现本发明目的的技术方案是：针对待分割噪声图像提取图像的邻域信息，构建能保持图像细节和去除噪声的两个目标函数，并采用和反对学习相结合的基于分解的多目标优化算法进行求解，从求得的Pareto前端中自动选取出最适合的解作为最终分割结果。其具体步骤包括如下：

(1)输入一幅含像素数为N的噪声图像，采用3×3的矩形窗提取出图像的邻域空间信息和邻域灰度信息；

(2)根据步骤(1)中的邻域空间信息和邻域灰度信息，确定图像相邻像素之间的约束关系；

(3)根据步骤(1)输入的噪声图像生成多目标模糊聚类的初始种群，并进行初始设定：

3a)设定聚类类数为c，当前迭代次数t＝0，最大迭代次数为t_max，其中c和t_max均大于0；

3b)将待分割图像随机生成100组聚类中心作为初始种群中的个体，共同构成多目标模糊聚类的初始种群；

3c)对初始种群中的所有个体定义其权值向量：λ^(l)＝[λ_l,1-λ_l]^T，其中λ_l＝l/100，l表示个体在初始种群中的序号，l的取值范围是[1,100]；

(4)根据权值向量λ^(l)和步骤(2)中的相邻像素之间的约束关系，构建初始种群中所有个体的图像细节保持目标函数f₁ ^(l)和去除噪声目标函数f₂ ^(l)，以及隶属度矩阵U^(l)；

(5)根据所述的两个目标函数f₁ ^(l)、f₂ ^(l)和权值向量λ^(l)，利用基于分解的多目标优化算法，计算种群中所有个体对应的子问题目标函数g^(l)，其计算公式为g^(l)＝λ_lf₁ ^(l)+(1-λ_l)f₂ ^(l)，其中λ_l和1-λ_l是权值向量λ^(l)中的元素；

(6)根据上述子问题目标函数g^(l)，对步骤(3)中的初始种群进行反对学习，生成当前种群；

(7)对当前种群所有个体采用传统差分策略进行交叉操作，采用高斯变异策略进行变异操作，生成对应的候选个体，并执行如下操作：

7a)按照步骤(4)构建所有候选个体的两个目标函数f₁ ^(l)'、f₂ ^(l)'和隶属度矩阵U^(l)'；

7b)按照步骤(5)计算所有候选个体的子问题目标函数g^(l)'；

(8)根据候选个体的子问题目标函数g^(l)'更新步骤(6)中的当前种群，生成新种群；

(9)判断当前迭代次数t是否等于最大迭代次数t_max，如果当前迭代次数t小于最大迭代次数t_max，则令t＝t+1，对步骤(8)中的新种群按照步骤(6)进行反对学习，循环执行步骤(6)～步骤(8)，直到满足循环结束条件；若当前迭代次数t等于最大迭代次数t_max，则停止循环，将步骤(8)中的新种群作为最终种群，执行步骤(10)；

(10)对最终种群进行非支配排序得到Pareto前端，利用基于角度的拐点确定法选取Pareto前端的拐点，将拐点对应的个体的隶属度矩阵进行去模糊，输出最终分割结果。

本发明具有如下优点：

(A)本发明将图像的邻域空间信息和邻域灰度信息，运用于目标函数的构建，有助于去除噪声。

(B)本发明构建了能保持图像细节和去除噪声的两个目标函数，并采用基于分解的多目标优化算法进行同时优化，在图像分割的过程中实现了图像细节的保持和噪声的去除。

(C)本发明采用的反对学习，加快了多目标优化算法的搜索速度，提高了图像分割的效率。

(D)本发明求得的Pareto前端上的每个点均对应着一个可行的分割结果，其中选取出的拐点对应着最适合的图像分割结果。

实验证明，本发明对含有噪声的图像能够进行既保持图像细节又去除噪声的分割。

附图说明

图1是本发明的整体实现流程图；

图2是本发明的反对学习子流程图；

图3是用本发明与现有NNcut和FCM方法，对flower自然图像的分割结果对比图；

图4是用本发明与现有NNcut和FCM方法，对coins自然图像的分割结果对比图；

图5是用本发明与现有NNcut和FCM方法，对MR医学图像的分割结果对比图。

具体实施方式

参照图1，本发明的具体实现步骤如下：

步骤1，输入一幅噪声图像并提取邻域信息，确定图像中像素之间的约束关系。

在本发明实施例中，输入一幅含像素数为N的噪声图像，采用3×3矩形窗对该噪声图像提取邻域信息，其中包括图像的邻域空间信息和邻域灰度信息。

步骤2，根据邻域信息，确定待分割图像中相邻像素间的约束关系。

2a)将待分割图像的所有像素进行转换，令x_r为图像中某个像素的灰度值，其转换公式如下：

$C_r=\frac{\mathrm{var}\left(x\right)}{{\left[\mathrm{mean}\left(x\right)\right]}^2},$

其中，C_r表示像素x_r转换后的灰度值，x是像素x_r作为中心像素的矩形窗中所有像素组成的集合；

2b)根据经过上述转换后的灰度值，确定相邻像素之间的灰度差异约束：

令x_i为矩形窗内的中心像素，x_j为中心像素x_i的某个邻居像素，C_j为邻居像素x_j经过步骤2a)转换后的灰度值，按照以下公式计算中心像素x_i与邻居像素x_j之间的灰度差异ξ_ij：

${\mathrm{\ξ}}_{\mathrm{ij}}=\exp [-(C_j-\frac{1}{n_i}\underset{\underset{k\≠i}{k\∈N_i}}{\mathrm{\Σ}}{C}_k)],$

式中，N_i是由中心像素x_i的所有邻居像素组成的集合，n_i是N_i所含元素的个数，根据上述灰度差异ξ_ij确定中心像素x_i与邻居像素x_j之间的灰度差异约束w_gc：

$w_{\mathrm{gc}}=\left\{\begin{array}{cc}2+{\mathrm{\&xi;}}_{\mathrm{ij}}/\underset{\underset{k\&NotEqual;i}{k\&Element;N_i}}{\mathrm{\&Sigma;}}{\mathrm{\&xi;}}_{\mathrm{ik}}& C_j<\stackrel{\&OverBar;}{C}\\ 2-{\mathrm{\&xi;}}_{\mathrm{ij}}\underset{\underset{k\&NotEqual;i}{k\&Element;N_i}}{\mathrm{\&Sigma;}}{\mathrm{\&xi;}}_{\mathrm{ik}}& C_j\&GreaterEqual;\stackrel{\&OverBar;}{C}\end{array}\right.,$

式中，表示矩形窗内所有像素通过步骤2a)转换后的灰度平均值；

2c)确定相邻像素之间的邻域空间约束：

$w_{\mathrm{sc}}=\frac{1}{d_{\mathrm{ij}}+1},$

其中，w_sc表示步骤2b)中的中心像素x_i与邻居像素x_j之间的邻域空间约束，d_ij表示中心像素x_i与邻居像素x_j在图像中的空间位置之间的欧式距离；

2d)根据步骤2b)中的灰度差异约束w_gc和步骤2c)中的邻域空间约束w_sc确定相邻像素之间的约束关系为：

w_ij＝w_sc·w_gc，

其中，w_ij表示中心像素x_i与邻居像素x_j之间的约束关系。

步骤3,根据步骤1输入的噪声图像生成多目标模糊聚类的初始种群，进行初始设定。

3a)设定聚类数目为c，当前迭代次数为t＝0，最大迭代次数为t_max，其中c和t_max取值均大于0，在本发明实施例中设定t_max＝30；

3b)将待分割图像随机生成100组聚类中心，每组聚类中心作为初始种群中的一个个体，共同构成初始种群；

3c)对初始种群中的所有个体定义其权值向量λ^(l)＝[λ_l,1-λ_l]^T，其中λ_l＝l/100，l表示个体在初始种群中的序号，l的取值范围是[1,100]，权值向量λ^(l)和种群中的个体是一一对应的关系，因此初始种群中100个个体对应100个权值向量；

3d)从上述100个权值向量中选取与λ^(l)欧式距离最近的20个权值向量，共同组成λ^(l)的邻居集合B^(l)。

步骤4，根据权值向量λ^(l)和步骤2中相邻像素间的约束关系，构建初始种群中所有个体的图像细节保持目标函数f₁ ^(l)和去除噪声目标函数f₂ ^(l)，以及隶属度矩阵U^(l)。

4a)令z_l表示初始种群中第l个个体，l的取值范围是[1,100]，z_l,p表示初始种群个体z_l的第p个聚类中心，计算中心像素x_i与聚类中心z_l,p之间的距离D²(x_i,z_l,p)，以及邻居像素x_j与聚类中心z_p之间的距离D²(x_j,z_l,p)：

$D^2(x_i,z_{l,p})=1-\exp (-{\left|\right|x_i-z_{l,p}\left|\right|}^2/\mathrm{\&sigma;}),$

$D^2(x_j,z_{l,p})=1-\exp (-{\left|\right|x_j-z_{l,p}\left|\right|}^2/\mathrm{\&sigma;}),$

其中，σ表示带宽， $\mathrm{\&sigma;}=\sqrt{\frac{1}{N-1}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{({\left|\right|x_i-\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{x}_i\left|\right|}^2-\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\left|\right|x_i-\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{x}_i\left|\right|}^2)}^2}.$

式中，N是待分割图像所含的像素数；

4b)根据步骤4a)所述距离D²(x_i,z_l,p)和D²(x_j,z_l,p)、权值向量λ^(l)以及步骤2中相邻像素间的约束关系，构建隶属度矩阵U^(l)：

$U^{\left(l\right)}=\left\{u_{\mathrm{ip}}^{\left(l\right)}\right|i\&Element;[1,N],p\&Element;[1,c]\},$

其中表示中心像素x_i在初始种群个体z_l的第p个聚类中心z_l,p的隶属度，其计算公式如下：

$u_{\mathrm{ip}}^{\left(l\right)}=\frac{{\mathrm{\&lambda;}}_l{D}^2(x_i,z_{l,p})+(1-{\mathrm{\&lambda;}}_l)\underset{\underset{j\&NotEqual;i}{j\&Element;N_i}}{\mathrm{\&Sigma;}}{(1-{\hat{u}}_{\mathrm{jp}}^{\left(l\right)})}^2{w}_{\mathrm{ij}}{D}^2(x_j,z_{l,p})}{\underset{q=1}{\overset{c}{\mathrm{\&Sigma;}}}({\mathrm{\&lambda;}}_l{D}^2(x_i,z_{l,q})+(1-{\mathrm{\&lambda;}}_l)\underset{\underset{j\&NotEqual;i}{j\&Element;N_i}}{\mathrm{\&Sigma;}}{(1-{\hat{u}}_{\mathrm{jp}}^{\left(l\right)})}^2{w}_{\mathrm{ij}}{D}^2(x_j,z_{l,q}))},$

式中，c是聚类类数，λ_l和1-λ_l是权值向量λ^(l)中的元素，w_ij是中心像素x_i与邻居像素x_j之间的约束关系，N_i是中心像素x_i的所有邻居像素组成的集合，表示在上一次迭代结束时邻居像素x_j在种群中第l个个体的第p个聚类中心的隶属度，当迭代次数t＝0时，取值为0；

4c)根据步骤4a)所述距离D²(x_i,z_l,p)、D²(x_j,z_l,p)和隶属度矩阵U^(l)，分别按如下公式构建初始种群个体z_l的图像细节保持目标函数f₁ ^(l)和去除噪声目标函数f₂ ^(l)：

${f_1}^{\left(l\right)}=\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}\underset{p=1}{\overset{c}{\mathrm{\&Sigma;}}}{\left(u_{\mathrm{ip}}^{\left(l\right)}\right)}^2{D}^2(x_i,z_{l,p}),$

${f_2}^{\left(l\right)}=\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}\underset{p=1}{\overset{c}{\mathrm{\&Sigma;}}}{\left(u_{\mathrm{ip}}^{\left(l\right)}\right)}^2\underset{\underset{j\&NotEqual;i}{j\&Element;N_i}}{\mathrm{\&Sigma;}}{(1-u_{\mathrm{jp}}^{\left(l\right)})}^2{w}_{\mathrm{ij}}{D}^2(x_j,z_{l.p}),$

式中，表示邻居像素x_j在初始种群个体z_l的第p个聚类中心z_l,p的隶属度。

步骤5，根据步骤4所述两个目标函数f₁ ^(l)、f₂ ^(l)和权值向量λ^(l)，利用基于分解的多目标优化算法，计算初始种群中所有个体对应的子问题目标函数：g^(l)＝λ_lf₁ ^(l)+(1-λ_l)f₂ ^(l)，其中λ_l和1-λ_l是权值向量λ^(l)中的两个元素。

步骤6，对步骤3中的初始种群进行反对学习，生成当前种群。

反对学习是一种已有的机器学习方法，其主要优势在于能够加快算法的搜索速度，使算法更快搜索到最优解。

参照图2，反对学习的具体实现如下：

6a)对初始种群中所有个体按如下公式求反：

其中，z_l,p是初始种群个体z_l的第p个聚类中心，表示z_l,p的反，按上式对初始种群个体z_l包含的所有聚类中心求反，共同构成初始种群个体z_l的反MIN和MAX分别表示待分割图像中的最小和最大灰度值；

6b)构建初始种群所有个体的反的两个目标函数和隶属度矩阵

6b1)分别计算中心像素x_i和初始种群个体z_l的反的第p个聚类中心之间的距离以及邻居像素x_j和初始种群个体z_l的反的第p个聚类中心之间的距离

其中，σ是带宽， $\mathrm{\&sigma;}=\sqrt{\frac{1}{N-1}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{({\left|\right|x_i-\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{x}_i\left|\right|}^2-\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\left|\right|x_i-\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{x}_i\left|\right|}^2)}^2}.$

式中，N是待分割图像所含的像素数；

6b2)根据步骤6b1)所述距离和权值向量λ^(l)以及步骤2中相邻像素间的约束关系，构建初始种群个体z_l的反的隶属度矩阵

其中表示中心像素x_i在初始种群个体z_l的反的第p个聚类中心的隶属度，其计算公式如下：

式中，c是聚类类数，λ_l和1-λ_l是权值向量λ^(l)中的两个元素，w_ij是中心像素x_i与邻居像素x_j之间的约束关系，N_i是中心像素x_i的所有邻居像素组成的集合，表示在上一次迭代结束时邻居像素x_j在种群中第l个个体的第p个聚类中心的隶属度，当迭代次数t＝0时，取值为0；

6b3)根据步骤6b1)所述距离和隶属度矩阵分别按如下公式构建初始种群个体z_l的反的图像细节保持目标函数和去除噪声目标函数

式中，表示邻居像素x_j在初始种群个体z_l的反的第p个聚类中心的隶属度；

6c)根据步骤6b3)所述两个目标函数和权值向量λ^(l)，计算初始种群个体z_l的反对应的子问题目标函数：其中λ_l和1-λ_l是权值向量λ^(l)中的两个元素；

6d)比较初始种群中所有个体的子问题目标函数g^(l)与其反的子问题目标函数若则将初始种群个体z_l的反作为当前种群的第l个个体初始种群个体z_l的反的子问题目标函数作为当前种群个体的子问题目标函数初始种群个体z_l的反的隶属度矩阵作为当前种群个体的隶属度矩阵若则将初始种群个体z_l作为当前种群的第l个个体初将始种群个体z_l的子问题目标函数g^(l)作为当前种群个体的子问题目标函数初始种群个体z_l的隶属度矩阵U^(l)作为当前种群个体的隶属度矩阵在确定当前种群个体后，将λ^(l)作为当前种群个体的权值向量。

步骤7，对当前种群所有个体进行交叉和变异操作，生成对应的候选个体。

7a)对当前种群所有个体进行交叉操作，生成对应的中间个体：

在多目标优化算法中，常用的交叉策略有单点交叉、多点交叉和差分策略等。在本发明实施例中，采用传统差分策略DE/rand/1/bin对当前种群中所有个体进行交叉操作，其具体操作如下：

令和分别表示从权值向量λ^(l)的邻居集合B^(l)中随机抽取的两个权值向量对应的个体，利用和对当前种群个体进行交叉操作，生成当前种群个体的中间个体y_l：

$y_{l,p}={\tilde{z}}_{l,p}+F\&times;({\tilde{z}}_{l1,p}-{\tilde{z}}_{l2,p}),$

式中，是当前种群个体的第p个聚类中心，y_l,p是当前种群个体的中间个体y_l的第p个聚类中心，是当前种群个体的第p个聚类中心，是当前种群个体的第p个聚类中心，F是控制差分程度的参数，在本发明实施例中设定F＝0.5；

7b)对步骤7a)所述中间个体进行变异操作：

在多目标优化算法中，常用的变异策略有二进制变异和高斯变异等，本实施例采用高斯变异对上述当前种群个体的中间个体y_l进行变异操作，生成当前种群个体的候选个体y'_l：

$y_{l,p}^{\&prime;}=\left\{\begin{array}{cc}N(y_{l,p},(\mathrm{MAX}-\mathrm{MIN})/20)& \mathrm{rand}\&le;0.9\\ y_{l,p}& \mathrm{rand}>0.9\end{array}\right.,$

式中，y_l,p是当前种群个体的中间个体y_l的第p个聚类中心，y'_l,p是当前种群个体的候选个体y'_l的第p个聚类中心，N(y_l,p,(MAX-MIN)/20)是均值为y_l,p、方差为(MAX-MIN)/20的高斯分布，rand是随机数，取值范围在[0,1]之间；

7c)若当前种群个体的候选个体y'_l的第p个聚类中心y'_l,p的取值范围超出区间[MIN,MAX]，则对y'_l,p进行修复操作，在[MIN,MAX]内生成一个随机数作为y'_l,p的取值，再执行步骤7d)，若y'_l,p的取值范围未超出区间[MIN,MAX]，则直接执行步骤7d)；

7d)分别构建当前种群个体的候选个体y'_l的图像细节保持目标函数f₁ ^(l)'和去除噪声目标函数f₂ ^(l)'，以及隶属度矩阵U^(l)'：

7d1)分别计算中心像素x_i与当前种群个体的候选个体y'_l的第p个聚类中心y'_l,p之间的距离D²(x_i,y'_l,p)，以及邻居像素x_j与当前种群个体的候选个体y'_l的第p个聚类中心y'_l,p之间的距离D²(x_j,y'_l,p)：

$D^2(x_i,y_{l,p}^{\&prime;})=1-\exp (-{\left|\right|x_i-y_{l,p}^{\&prime;}\left|\right|}^2/\mathrm{\&sigma;}),$

$D^2(x_j,y_{l,p}^{\&prime;})=1-\exp (-{\left|\right|x_j-y_{l,p}^{\&prime;}\left|\right|}^2/\mathrm{\&sigma;}),$

其中，σ是带宽， $\mathrm{\&sigma;}=\sqrt{\frac{1}{N-1}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{({\left|\right|x_i-\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{x}_i\left|\right|}^2-\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\left|\right|x_i-\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{x}_i\left|\right|}^2)}^2}.$

式中，N是待分割图像所含的像素数；

7d2)根据步骤7d1)所得到的距离D²(x_i,y'_l,p)、D²(x_j,y'_l,p)和权值向量λ^(l)，以及步骤2中相邻像素间的约束关系，构建当前种群个体的候选个体y'_l的隶属度矩阵U^(l)'：

${U^{\left(l\right)}}^{\&prime;}=\left\{{u_{\mathrm{ip}}^{\left(l\right)}}^{\&prime;}\right|i\&Element;[1,N],p\&Element;[1,c]\},$

其中表示中心像素x_i在当前种群个体的候选个体y'_l的第p个聚类中心y'_l,p的隶属度，其计算公式如下：

${u_{\mathrm{ip}}^{\left(l\right)}}^{\&prime;}=\frac{{\mathrm{\&lambda;}}_l{D}^2(x_i,y_{l,p}^{\&prime;})+(1-{\mathrm{\&lambda;}}_l)\underset{\underset{j\&NotEqual;i}{j\&Element;N_i}}{\mathrm{\&Sigma;}}{(1-{\hat{u}}_{\mathrm{jp}}^{\left(l\right)})}^2{w}_{\mathrm{ij}}{D}^2(x_j,y_{l,p}^{\&prime;})}{\underset{q=1}{\overset{c}{\mathrm{\&Sigma;}}}({\mathrm{\&lambda;}}_l{D}^2(x_i,y_{l,p}^{\&prime;})+(1-{\mathrm{\&lambda;}}_l)\underset{\underset{j\&NotEqual;i}{j\&Element;N_i}}{\mathrm{\&Sigma;}}{(1-{\hat{u}}_{\mathrm{jp}}^{\left(l\right)})}^2{w}_{\mathrm{ij}}{D}^2(x_j,y_{l,p}^{\&prime;}))},$

式中，c是聚类类数，λ_l和1-λ_l是权值向量λ^(l)中的两个元素，w_ij是中心像素x_i与邻居像素x_j之间的约束关系，N_i是中心像素x_i的所有邻居像素组成的集合，表示在上一次迭代结束时邻居像素x_j在种群中第l个个体的第p个聚类中心的隶属度，当迭代次数t＝0时，取值为0；

7d3)根据步骤7d1)得到的距离D²(x_i,y'_l,p)、D²(x_j,y'_l,p)和隶属度矩阵U^(l)'，分别按如下公式构建当前种群个体的候选个体y'_l的图像细节保持目标函数f₁ ^(l)'和去除噪声目标函数f₂ ^(l)'：

${{f_1}^{\left(l\right)}}^{\&prime;}=\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}\underset{p=1}{\overset{c}{\mathrm{\&Sigma;}}}{\left({u_{\mathrm{ip}}^{\left(l\right)}}^{\&prime;}\right)}^2{D}^2(x_i,y_{l,p}^{\&prime;}),$

${{f_2}^{\left(l\right)}}^{\&prime;}=\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}\underset{p=1}{\overset{c}{\mathrm{\&Sigma;}}}{\left({u_{\mathrm{ip}}^{\left(l\right)}}^{\&prime;}\right)}^2\underset{\underset{j\&NotEqual;i}{j\&Element;N_i}}{\mathrm{\&Sigma;}}{(1-{u_{\mathrm{jp}}^{\left(l\right)}}^{\&prime;})}^2{w}_{\mathrm{ij}}{D}^2(x_j,y_{l,p}^{\&prime;}),$

式中，表示邻居像素x_j在当前种群个体的候选个体y'_l的第p个聚类中心y'_l,p的隶属度；

7e)根据步骤7d3)得到的两个目标函数f₁ ^(l)'、f₂ ^(l)'和权值向量λ^(l)，计算当前种群个体的候选个体y'_l对应的子问题目标函数：g^(l)'＝λ_lf₁ ^(l)'+(1-λ_l)f₂ ^(l)'，其中λ_l和1-λ_l是权值向量λ^(l)中的元素。

步骤8，根据当前种群所有个体的候选个体的子问题目标函数更新当前种群，生成新种群。

8a)从权值向量λ^(l)的邻居集合B^(l)中，随机抽取一个权值向量，其中l是当前种群中个体的序号，取值范围是[1,100]，例如随机抽取出的权值向量是λ^(l＝3)，权值向量λ^(l＝3)对应的是当前种群中第3个个体是当前种群个体的子问题目标函数，是当前种群个体的隶属度矩阵；

8b)比较当前种群个体的子问题目标函数和当前种群个体的候选个体y'_l的子问题目标函数g^(l)'：

若则将当前种群个体的候选个体y'_l作为新种群的第3个个体将当前种群个体的候选个体y'_l的子问题目标函数g^(l)'作为新种群个体的子问题目标函数将当前种群个体的候选个体y'_l的隶属度矩阵U^(l)'作为新种群个体的隶属度矩阵

若则将当前种群个体作为新种群的第3个个体将当前种群个体的子问题目标函数作为新种群个体的子问题目标函数将当前种群个体的隶属度矩阵作为新种群个体的隶属度矩阵

在确定新种群个体后，将λ^(l＝3)作为新种群个体的权值向量。

步骤9，判断当前迭代次数t是否等于最大迭代次数t_max，如果当前迭代次数t小于最大迭代次数t_max，则令t＝t+1，对步骤8得到的新种群循环执行步骤6～步骤8，直到满足循环结束条件；若当前迭代次数t等于最大迭代次数t_max，则停止循环，将步骤8得到的新种群作为最终种群，执行步骤10。在本发明实施例中，设定t_max＝30为循环结束条件。

步骤10，对最终种群进行非支配排序得到Pareto前端，利用基于角度的拐点确定法选取Pareto前端的拐点，并去模糊输出最终分割结果。

在对最终种群进行非支配排序得到的Pareto前端上，每个点均对应着一个可行的分割结果。在目前已有技术中，从Pareto前端选取最终结果的最常用方法是基于用户偏好的选择法。但是，基于用户偏好的选择法不能自动选择出最终结果，需要用户的参与才可以完成。因此本实施例利用基于角度的拐点确定法，选取Pareto前端的拐点作为最终分割结果，实现了分割结果的自动选择，避免了选择过程中的人为干涉。其中，基于角度的拐点确定法是多目标优化算法在优化两个目标函数时常用的拐点确定法，其具体步骤如下：

10a)将Pareto前端上的点和其左、右两边距离最近的两个点分别进行连线，计算两连线所夹角度，作为该点的角度；

10b)选取角度最大的点作为Pareto前端的拐点；

10c)将上述拐点对应的个体的隶属度矩阵去模糊，即提取出隶属度矩阵中每个像素在所有聚类中心中最大的隶属度，将该隶属度所属聚类中心的类别数作为对应像素的类标，得到的图像像素类标矩阵就是待分割图像的分割结果。

本发明的效果可以通过以下仿真实验具体说明：

1.实验条件

实验所用微机CPU为Intel Xeon E5-2620(2.0GHz2.0GHz)、内存64.0G，编程平台为Matlab7.10.0。

2.实验内容

实验1，输入一幅含像素数为N＝128×128，含有15％高斯噪声的flower自然图像，设置聚类类数c＝3。使用现有NNcut、FCM和本发明方法分别对噪声图像进行分割，结果如图3所示。其中，图3(a)为含有15％高斯噪声的flower自然图像；图3(b)、3(c)分别是NNcut和FCM方法得到的分割效果图；图3(d)为使用本发明方法得到的分割效果图。

实验2，输入一幅含像素数为N＝242×308，含有20％椒盐噪声的coins自然图像，设置聚类类数c＝3。使用现有NNcut、FCM和本发明方法分别对噪声图像进行分割，结果如图4所示。其中，图4(a)为含有20％椒盐噪声的coins自然图像；图4(b)、4(c)分别是NNcut和FCM方法得到的分割效果图；图4(d)为使用本发明方法得到的分割效果图。

实验3，输入一幅含像素数为N＝256×256，含有20％Rician噪声的MR医学图像，设置聚类类数c＝4。使用NNcut、FCM和本发明方法分别对噪声图像进行分割，结果如图5所示。其中，图5(a)为含有20％Rician噪声的MR医学图像；图5(b)、5(c)分别是NNcut和FCM方法得到的分割效果图；图5(d)为使用本发明方法得到的分割效果图。

3.实验结果分析

由图3(b)和3(c)可见，使用NNcut和FCM方法得到的分割结果含有很多噪声，并且边缘轮廓非常不清晰；对比图3(b)、3(c)和3(d)可见，使用本发明方法得到的分割结果不仅去除了图像噪声，而且拥有非常清晰的分割边界，因此获得了比NNcut和FCM方法更优秀的图像分割效果。

由图4(b)和4(c)可见，使用NNcut和FCM方法得到的分割结果含有很多噪声；对比图4(b)、4(c)和4(d)可见，使用本发明方法得到的分割结果既保持了图像细节，又去除了图像噪声，因此获得了比NNcut和FCM方法更优秀的图像分割效果。

由图5(b)和5(c)可见，使用NNcut和FCM方法得到的分割结果含有很多噪声；对比图5(b)、5(c)和5(d)可见，使用本发明方法得到的分割结果在去除了图像噪声的同时保持了图像的细节，并且分割边界非常清晰，因此得到了比NNcut和FCM方法更优秀的图像分割效果。

综上所述，本发明提出的基于多目标模糊聚类和反对学习的噪声图像分割方法，既能保持图像细节，又能去除图像噪声，同时提高了图像分割的效率，并从一系列可行的图像分割结果中自动选取出最适合的结果作为最终分割结果。

Claims (4)

1.一种基于多目标模糊聚类和反对学习的噪声图像分割方法，包括如下步骤：

(1)输入一幅含像素数为N的噪声图像，采用3×3的矩形窗提取出图像的邻域空间信息和邻域灰度信息；

(2)根据步骤(1)中的邻域空间信息和邻域灰度信息，确定图像相邻像素之间的约束关系；

(3)根据步骤(1)输入的噪声图像生成多目标模糊聚类的初始种群，并进行初始设定：

3a)设定聚类类数为c，当前迭代次数t＝0，最大迭代次数为t_max，其中c和t_max均大于0；

3b)将待分割图像随机生成100组聚类中心作为初始种群中的个体，共同构成多目标模糊聚类的初始种群；

3c)对初始种群中的所有个体定义其权值向量：λ^(l)＝[λ_l,1-λ_l]^T，其中λ_l＝l/100，l表示个体在初始种群中的序号，l的取值范围是[1,100]；

(4)根据权值向量λ^(l)和步骤(2)中的相邻像素之间的约束关系，构建初始种群中所有个体的图像细节保持目标函数f₁ ^(l)和去除噪声目标函数f₂ ^(l)，以及隶属度矩阵U^(l)；

(5)根据所述的两个目标函数f₁ ^(l)、f₂ ^(l)和权值向量λ^(l)，利用基于分解的多目标优化算法，计算种群中所有个体对应的子问题目标函数g^(l)，其计算公式为g^(l)＝λ_lf₁ ^(l)+(1-λ_l)f₂ ^(l)，其中λ_l和1-λ_l是权值向量λ^(l)中的元素；

(6)根据上述子问题目标函数g^(l)，对步骤(3)中的初始种群进行反对学习，生成当前种群：

6a)对种群中所有个体按照如下公式求反：

其中，z_l,p是个体z_l的第p个聚类中心，表示z_l,p的反，按上式对个体z_l包含的所有聚类中心求反，共同构成个体z_l的反MIN和MAX分别表示待分割图像中的最小和最大灰度值；

6b)按步骤(4)构建所有反的两个目标函数和隶属度矩阵

6c)按步骤(5)计算所有反的子问题目标函数

6d)比较种群中所有个体的子问题目标函数g^(l)与其反的子问题目标函数保留两者中取值较小的作为当前种群的第l个个体；

(7)对当前种群所有个体采用传统差分策略进行交叉操作，采用高斯变异策略进行变异操作，生成对应的候选个体，并执行如下操作：

7a)按照步骤(4)构建所有候选个体的两个目标函数f₁ ^(l)'、f₂ ^(l)'和隶属度矩阵U^(l)'；

7b)按照步骤(5)计算所有候选个体的子问题目标函数g^(l)'；

(8)根据候选个体的子问题目标函数g^(l)'更新步骤(6)中的当前种群，生成新种群；

(9)判断当前迭代次数t是否等于最大迭代次数t_max，如果当前迭代次数t小于最大迭代次数t_max，则令t＝t+1，对步骤(8)中的新种群按照步骤(6)进行反对学习，循环执行步骤(6)～步骤(8)，直到满足循环结束条件；若当前迭代次数t等于最大迭代次数t_max，则停止循环，将步骤(8)中的新种群作为最终种群，执行步骤(10)；

(10)对最终种群进行非支配排序得到Pareto前端，利用基于角度的拐点确定法选取Pareto前端的拐点，将拐点对应的个体的隶属度矩阵进行去模糊，输出最终分割结果。

2.根据权利要求1所述的基于多目标模糊聚类和反对学习的噪声图像分割方法，其中所述步骤(2)中确定图像相邻像素之间的约束关系，按如下步骤确定：

2a)将待分割图像的所有像素进行转换，令x_r为图像中某个像素的灰度值，其转换公式如下：

$C_r=\frac{\mathrm{var}\left(x\right)}{{\&lsqb;mean\left(x\right)\&rsqb;}^2},$

其中，C_r表示像素x_r转换后的灰度值，x是像素x_r作为中心像素的矩形窗中所有像素组成的集合；

2b)根据经过上述转换后的灰度值，确定相邻像素之间的灰度差异约束：

令x_i为矩形窗内的中心像素，x_j为中心像素x_i的某个邻居像素，C_j为邻居像素x_j经过步骤2a)转换后的灰度值，按照以下公式计算中心像素x_i与邻居像素x_j之间的灰度差异ξ_ij：

${\mathrm{\&xi;}}_{ij}=\exp \&lsqb;-(C_j-\frac{1}{n_i}\underset{k\&NotEqual;i}{\underset{k\&Element;N_i}{\&Sigma;}}{C}_k)\&rsqb;,$

式中，N_i是由中心像素x_i的所有邻居像素组成的集合，n_i是N_i所含元素的个数，根据上述灰度差异ξ_ij确定中心像素x_i与邻居像素x_j之间的灰度差异约束w_gc：

$w_{gc}=\left\{\begin{array}{cc}2+{\mathrm{\&xi;}}_{ij}/\underset{k\&NotEqual;i}{\underset{k\&Element;N_i}{\&Sigma;}}{\mathrm{\&xi;}}_{ik}& C_j<\stackrel{\&OverBar;}{C}\\ 2-{\mathrm{\&xi;}}_{ij}/\underset{k\&NotEqual;i}{\underset{k\&Element;N_i}{\&Sigma;}}{\mathrm{\&xi;}}_{ik}& C_j\&GreaterEqual;\stackrel{\&OverBar;}{C}\end{array}\right.,$

式中，表示矩形窗内所有像素通过步骤2a)转换后的灰度平均值；

2c)确定相邻像素之间的邻域空间约束：

$w_{sc}=\frac{1}{d_{ij}+1},$

其中，w_sc表示步骤2b)中的中心像素x_i与邻居像素x_j之间的邻域空间约束，d_ij表示中心像素x_i与邻居像素x_j在图像中的空间位置之间的欧式距离；

2d)根据步骤2b)中的灰度差异约束w_gc和步骤2c)中的邻域空间约束w_sc确定相邻像素之间的约束关系为：

w_ij＝w_sc·w_gc，

其中，w_ij表示中心像素x_i与邻居像素x_j之间的约束关系。

3.根据权利要求1所述的基于多目标模糊聚类和反对学习的噪声图像分割方法，其中所述步骤(4)中构建种群中所有个体的图像细节保持目标函数f₁ ^(l)和去除噪声目标函数f₂ ^(l)，以及隶属度矩阵U^(l)，按如下公式构建：

${f_1}^{\left(l\right)}=\underset{i=1}{\overset{N}{\&Sigma;}}\underset{p=1}{\overset{c}{\&Sigma;}}{\left(u_{ip}^{\left(l\right)}\right)}^2{D}^2(x_i,z_{l,p}),$

${f_2}^{\left(l\right)}=\underset{i=1}{\overset{N}{\&Sigma;}}\underset{p=1}{\overset{c}{\&Sigma;}}{\left(u_{ip}^{\left(l\right)}\right)}^2\underset{j\&NotEqual;i}{\underset{j\&Element;N_i}{\&Sigma;}}{(1-u_{jp}^{\left(1\right)})}^2{w}_{ij}{D}^2(x_j,z_{l,p}),$

$U^{\left(l\right)}=\left\{u_{ip}^{\left(l\right)}\right|i\&Element;\&lsqb;1,N\&rsqb;,p\&Element;\&lsqb;1,c\&rsqb;\},$

$u_{ip}^{\left(l\right)}=\frac{{\mathrm{\&lambda;}}_l{D}^2(x_i,z_{l,p})+(1-{\mathrm{\&lambda;}}_l)\underset{j\&NotEqual;i}{\underset{j\&Element;N_i}{\&Sigma;}}{(1-{\hat{u}}_{jp}^{\left(l\right)})}^2{w}_{ij}{D}^2(x_j,z_{l,p})}{\underset{q=1}{\overset{c}{\&Sigma;}}\left({\mathrm{\&lambda;}}_l{D}^2\right(x_i,z_{l,q})+(1-{\mathrm{\&lambda;}}_l\left)\underset{j\&NotEqual;i}{\underset{j\&Element;N_i}{\&Sigma;}}{\left(1-{\hat{u}}_{jp}^{\left(l\right)}\right)}^2{w}_{ij}{D}^2\right(x_j,z_{l,q}\left)\right)},$

其中，个体z_l表示种群中第l个个体，λ_l和1-λ_l是权值向量λ^(l)中的元素，l的取值范围是[1,100]，z_l,p表示个体z_l的第p个聚类中心，N是待分割图像所含的像素数，c是聚类类数，w_ij是中心像素x_i与邻居像素x_j之间的约束关系，N_i是中心像素x_i的所有邻居像素组成的集合，表示在上一次迭代结束时邻居像素x_j在种群中第l个个体的第p个聚类中心的隶属度，当迭代次数t＝0时，取值为0，表示中心像素x_i在个体z_l的第p个聚类中心z_l,p的隶属度，表示邻居像素x_j在个体z_l的第p个聚类中心z_l,p的隶属度，D²(x_i,z_l,p)表示中心像素x_i和聚类中心z_l,p之间的距离，D²(x_j,z_l,p)表示邻居像素x_j和聚类中心z_l,p之间的距离，其计算公式分别如下：

D²(x_i,z_l,p)＝1-exp(-||x_i-z_l,p||²/σ)，

D²(x_j,z_l,p)＝1-exp(-||x_j-z_l,p||²/σ)，

其中，σ表示带宽，

4.根据权利要求1所述的基于多目标模糊聚类和反对学习的噪声图像分割方法，其中所述步骤(10)中利用基于角度的拐点确定法选取出Pareto前端的拐点，按如下步骤进行：

10a)将Pareto前端上的点和其左、右两边距离最近的两个点分别进行连线，计算两连线所夹角度，作为该点的角度；

10b)选取角度最大的点作为Pareto前端的拐点。