CN102360496B - 基于gmm参数迁移聚类的sar图像分割方法 - Google Patents

Abstract

本发明公开了一种基于GMM参数迁移聚类的SAR图像分割技术方法，主要解决现有技术分割结果不稳定、精度和区域一致性不理想的问题。其实现过程是：1)输入图像，确定分割类数；2)提取特征；3)设置初始参数；4)对样本聚类7次，得到7组聚类结果；5)求聚类一致性值；6)根据聚类一致性值把样本划分成源域和目标域；7)用EM算法对源域进行参数估计；8)寻找目标域样本在源域样本中的K近邻点并求出这些点的聚类一致性值；9)根据K近邻点的一致性值和源域参数，求目标域的新参数；10)根据新参数求目标域每个样本的概率值，得到图像最终分割结果。本发明具有分割效果稳定良好、区域一致性好的优点，可用于雷达目标检测和识别。

Description

基于GMM参数迁移聚类的SAR图像分割方法

技术领域

本发明属于图像处理技术领域，涉及SAR图像分割，可用于雷达目标检测和目标识别。

背景技术

合成孔径雷达SAR具有全天时、全天候的探测与侦察能力。它利用脉冲压缩技术获得高的距离分辨率，利用合成孔径原理提高方位分辨率，从而相比真实孔径雷达在遥感领域具有独特的优势。对SAR图像的理解和解译属于图像处理范畴，还涉及信号处理，模式识别及机器学习等众多学科。由于SAR具有的独特作用，使得SAR图像的理解与解译在国防和民用领域正受到越来越广泛的关注，SAR图像分割作为SAR图像后续解译处理的关键环节之一，也就显得愈加重要。现有的SAR图像分割方法大致可以分为基于区域分割方法和基于边缘的分割方法，基于区域的分割方法有聚类、阈值分割、区域生长法等，基于边缘的分割方法有边界跟踪法和图像滤波法等。其中基于聚类的SAR图像分割方法，是将SAR图像中具有某方面相似特征的区域划分成一类。聚类中常用的技术是机器学习，并已取得了显著的成功，但是许多机器学习方法由于都是基于训练数据和测试数据来自相同的分布和同一个特征空间的假设，因而当数据分布改变时，大多数的机器学习方法需要从头开始学习，要求使用者重新收集大量的训练数据来训练模型。而在真实世界中重新收集数据然后标记需要花费很昂贵的代价。此外，经典机器学习方法EM算法因简单、容易实现而被广泛使用，但它本身还存在对初始值敏感，容易陷入局部最优的缺陷，造成图像分割的不稳定以及错分现象。

发明内容

本发明的目的在于克服上述已有技术的不足，提出了一种基于GMM参数迁移聚类的SAR图像分割方法，以合理利用服从不同分布下的数据，提高SAR图像分割的稳定性和分割效果。

实现本发明的技术思路是：从待分割SAR图像中提取小波特征为待划分的样本，初始聚类几次，通过聚类一致性值把样本划分成易分与不易分两类，即源域和目标域，然后在源域中学习参数以指导目标域的划分，最终得到图像分割结果。具体实现步骤包括如下：

(1)输入待分割SAR图像，确定需要将图像划分的类别数C；

(2)对SAR图像提取特征：对待分割的SAR图像进行3层平稳小波变换，图像像素点总个数为n，对每个像素点提取出10维子带能量特征，构成大小为n×10的输入数据样本集X＝{x_i|i＝1，2，...，n}：

$x_i=\frac{1}{M\×N}\underset{e=1}{\overset{M}{\mathrm{\Σ}}}\underset{u=1}{\overset{N}{\mathrm{\Σ}}}\left|\mathrm{coef}(e,u)\right|$

其中，M×N为利用滑动窗口确定的子带大小，这里取值为16×16，coef(e，u)为平稳小波子带中第e行第u列的系数值；

(3)设置初始参数：样本近邻数K＝7，用k-means聚类算法获得样本初始标签L_re，聚类一致性值的阈值thre＝0.6；

(4)对数据样本集X分别用EM算法聚类4次，用k-means聚类3次，得到7组聚类结果；

(5)对7组聚类结果求每个样本的聚类一致性值CI；

(6)根据CI和thre，把数据样本集X划分成源域和目标域：CI≥0.6的样本为源域X^S，CI＜0.6为目标域X^T；

(7)用基于GMM的EM算法对源域X^S参数估计，得到GMM的均值μ_j和协方差σ_j两个类参数，j表示类别；

(8)用K近邻算法寻找目标域X^T的每个样本在源域中的K个近邻点(K＝6)，并根据步骤(4)中获得的标签求6个近邻点的聚类一致性值CI′；

(9)根据CI′并结合步骤(7)中获得的源域参数，计算目标域的均值μ_jnew和协方差σ_jnew两个新参数，即将源域参数迁移至目标域；

(10)将所述两个新参数代入高斯分布的公式求得目标域每个样本的概率值，得到图像的最终划分结果。

本发明与现有技术相比具有以下优点：

(1)本发明通过寻找目标域样本在源域样本中的K近邻点，由近邻点的聚类一致性值来影响目标域参数的变化，实际上是加入了样本的局部空间信息，从而能有效提高聚类算法的性能。

(2)本发明通过迁移学习将GMM的源域参数迁移至目标域，使得目标域参数更加的接近真实参数，从而分割效果稳定良好，分割图像区域一致性好，能够更加准确的识别地物。

仿真实验表明本发明在细节划分、区域一致性上都优于传统EM算法。

附图说明

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

图2是本发明与现有方法对包含山地和平坦地区的两类SAR图像的仿真分割对比图；

图3是本发明与现有方法对包含森林和田地两类SAR图像的仿真分割对比图。

具体实施方式

参照图1，本发明的具体实施过程如下：

步骤1.输入待分割SAR图像，根据图像内容，判断需要识别的主要目标及背景，确定分割类别数C，本实例中的C取值为2。

步骤2.提取待分割SAR图像的特征。

SAR图像不仅数据量大，在成像过程中不同的地物有不同的后向发射和散射特性，从而具有丰富的幅度、相位、极化和纹理信息，并且，图像固有的相干斑噪声对于分割性能直接产生影响，因此，有必要在图像分割前对SAR图像进行纹理分析，提取有效的纹理特征进行聚类。

在以上分析基础上，对待分割的SAR图像进行3层平稳小波变换，图像像素点总个数为n，对每个像素点提取出10维子带能量特征，构成大小为n×10的输入数据样本集X＝{x_i|i＝1，2，...，n}：

$x_i=\frac{1}{M\×N}\underset{e=1}{\overset{M}{\mathrm{\Σ}}}\underset{u=1}{\overset{N}{\mathrm{\Σ}}}\left|\mathrm{coef}(e,u)\right|$

其中，M×N为利用滑动窗口确定的子带大小，这里取值为16×16，coef(e，u)为平稳小波子带中第e行第u列的系数值，这样，对于本实例中大小为256×256的SAR图像，就构成65536×10的矩阵作为GMM参数迁移聚类方法的输入数据。

步骤3.设置初始参数：样本近邻数K为7，聚类一致性值的阈值thre＝0.6。

步骤4.对数据样本集X分别用基于GMM的EM算法聚类4次、k-means聚类算法聚类3次。

(4a)对数据样本集X＝{x_i|i＝1，2，...，n}随机生成一个归属矩阵z_ij

$z_{\mathrm{ij}}=\left\{\begin{array}{cc}1& x_i\∈j\\ 0& x_i\∉j\end{array}\right.,i=\mathrm{1,2},...,n;j=\mathrm{1,2},...,C$

其中样本x_i＝{x_i1，x_i2，...，x_im}，m为维数，n为数据个数，C为类别数，j为类别，若z_ij＝1表示样本x_i属于第j类，z_ij＝0表示样本x_i不属于第j类；

(4b)利用如下公式求GMM初始类混合权重α_j、类均值μ_j和类方差σ_j

α_j＝1/C

${\mathrm{\μ}}_j=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{z}_{\mathrm{ij}}{x}_{\mathrm{il}}/\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{z}_{\mathrm{ij}},l=\mathrm{1,2},...,m$

${\mathrm{\σ}}_j=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{z}_{\mathrm{ij}}(x_i-{\mathrm{\μ}}_j){(x_i-{\mathrm{\μ}}_j)}^T/\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{z}_{\mathrm{ij}};$

(4c)根据上述类参数求得任意一项样本x_i的概率为：

$p\left(x_i\right)=\underset{j=1}{\overset{C}{\mathrm{\Σ}}}{\mathrm{\α}}_{j}G(x_i;{\mathrm{\μ}}_j,{\mathrm{\σ}}_j),i=\mathrm{1,2},...,n$

其中 $G(x_i,{\mathrm{\μ}}_j,{\mathrm{\σ}}_j)=\frac{1}{{\left(2\mathrm{\π}\right)}^{\frac{m}{2}}{\left|{\mathrm{\σ}}_j\right|}^{\frac{1}{2}}}\exp \{-\frac{1}{2}{(x_i-{\mathrm{\μ}}_j)}^T{{\mathrm{\σ}}_j}^{-1}(x_i-{\mathrm{\μ}}_j)\},$ 表示样本x_i在第j类的高斯分布概率，m为维数；

(4d)利用上述所求得的概率计算步骤(4a)所生成归属矩阵z_ij的期望值E；

E＝E(z_ij|x_i；α_j，μ_j，σ_j)＝α_jG(x_i；μ_j，σ_j)/p(x_i)

其中E＝E(z_ij|z_i；α_j，μ_j，σ_j)表示期望值E是在类混合权重α_j、类均值μ_j和类方差σ_j条件下求得的样本x_i的归属矩阵z_ij的期望值；

(4e)用产生的期望值E替换类均值μ_j和类方差σ_j公式中的z_ij，得到新的类均值μ_j′和类方差σ_j′，并用下式更新类混合权重得到α_j′：

${{\mathrm{\α}}_j}^{\′}=\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}E;$

n是图像像素点总个数；

(4f)将归属矩阵z_ij的期望值E最大化：

首先，用所述更新后的类均值μ_j′、类方差σ_j′和类混合权重α_j′重新进行步骤(4c)-(4e)，得到归属矩阵z_ij的新期望E′；

然后，反复迭代步骤(4c)-(4e)直至|E_t+1′-E_t′|≤0.000001循环结束，E_t′表示第t次循环得到的z_ij的期望，E_t+1′表示第t+1次循环得到的z_ij的期望；

(4g)利用z_ij的期望E′更新归属矩阵z_ij，得到样本的最终划分：

z_ij＝1当且仅当 $j=\arg \underset{j=1}{\overset{C}{\max }}\left(E^{\′}\right)$

z_ij＝0当且仅当 $j\≠\arg \underset{j=1}{\overset{C}{\max }}\left(E^{\′}\right);$

(4h)重复(4a)-(4g)四次，即可以得到用EM算法聚类的4组结果；

(4i)对数据样本集X用k-means算法聚类3次，得到另外的3组聚类结果。

步骤5.根据下式计算7组聚类结果的聚类一致性值CI：

$\mathrm{CI}\left(x_i\right)=\frac{1}{H}\max {\left\{\mathrm{\Σ\δ}(\mathrm{\π}\left(x_i\right),L)\right\}}_{L\∈\{\mathrm{1,2},...,C\}}$

其中，H表示样本x_i的聚类次数，此处H＝7，π(x_i)表示样本x_i的聚类结果，L为每一个类别标号，∑δ(π(x_i)，L)表示类别L在π(x_i)中的总个数，C是图像划分的类别数。

步骤6.根据聚类一致性值CI和聚类一致性值的阈值thre，把样本集X划分成源域X^S和目标域X^T：CI≥0.6的样本为源域X^S，CI＜0.6为目标域X^T；

步骤7.重复步骤(4a)-(4f)对源域X^S进行参数估计，得到源域X^S的均值μ_j和协方差σ_j两个类参数，j表示类别。

步骤8.用K近邻算法寻找目标域X^T的每个样本在源域中的6个近邻点，并根据步骤(4)中获得的聚类结果求这6个近邻点的聚类一致性值CI′，其计算公式与步骤(5)中的聚类一致性值CI的公式相同。

步骤9.根据所述聚类一致性值CI′并结合步骤(7)中获得的源域的均值μ_j和协方差σ_j，计算目标域的均值μ_jnew和协方差σ_jnew这两个新参数：

μ_jnew＝CI′μ_j

σ_jnew＝CI′σ_j

步骤10.根据所计算的两个新参数求目标域X^T每个样本的概率值，得到图像的最终划分。

(10a)将目标域的均值μ_jnew和协方差σ_jnew代入如下高斯分布公式，计算目标域每个样本x_i属于第j类的概率P_ij：

$P_{\mathrm{ij}}=\frac{1}{{\left(2\mathrm{\π}\right)}^{\frac{m}{2}}{\left|{\mathrm{\σ}}_{\mathrm{jnew}}\right|}^{\frac{1}{2}}}\exp \{-\frac{1}{2}{(x_i-{\mathrm{\μ}}_{\mathrm{jmew}})}^T{{\mathrm{\σ}}_{\mathrm{jnew}}}^{-1}(x_i-{\mathrm{\μ}}_{\mathrm{jnew}})\}$

m为每个样本x_i的维数；

(10b)通过如下公式计算归属矩阵z_ij，获得图像的最终分割结果：

z_ij＝1当且仅当 $j=\arg \underset{j=1}{\overset{C}{\max }}\left(P_{\mathrm{ij}}\right)$

z_ij＝0当且仅当 $j\≠\arg \underset{j=1}{\overset{C}{\max }}\left(P_{\mathrm{ij}}\right).$

本发明的效果可以通过如下仿真实验进行验证：

(一)实验条件设置

选择2幅SAR图像进行了验证，图像大小为256×256，划分类别数C取值为2。

(二)仿真内容及结果

A)用现有的EM聚类算法和本发明方法对图2(a)包含山地和平坦地区的两类SAR图像进行仿真实验，分割结果如图2所示。其中图2(b)为EM算法分割结果，图2(c)为本发明分割结果。

从图2(b)可以看出，EM分割结果未分割出向纵深延伸的山脊，同时山地部分也存在很多误分。

从图2(c)可以看出，本发明方法由于考虑到图像的局部信息，分割出了一些向纵深延伸的山脊，分割结果明显地呈现山地和平坦地区两部分，这也是符合地质勘探上对地形整体把握需要的。

B)用现有的EM聚类算法以及本发明方法对图3(a)包含森林和田地两类SAR图像进行仿真实验，分割结果如图3所示。其中图3(b)为EM算法分割结果，图3(c)为本发明方法分割结果。

从图3(b)可以看出，EM聚类的结果中不仅田地中混杂了森林区域，而且轮廓划分也较为粗糙。

从图3(c)可以看出，本发明方法的分割结果在森林和田地的一致性有很大的改善，同时右下角的细节部分保留较好。

Claims (4)

1.一种基于GMM参数迁移聚类的SAR图像分割方法，包括以下步骤： 

（1）输入待分割SAR图像，确定需要将图像划分的类别数C； 

（2）对SAR图像提取特征：对待分割的SAR图像进行3层平稳小波变换，图像像素点总个数为n，对每个像素点提取出10维子带能量特征，构成大小为n×10的输入数据样本集X＝{x_i|i＝1,2,...,n}： 

其中，M×N为利用滑动窗口确定的子带大小，这里取值为16×16，coef(e,u)为平稳小波子带中第e行第u列的系数值； 

（3）设置初始参数：样本近邻数K＝7，用k-means聚类算法获得样本初始标签L_re，聚类一致性值的阈值thre＝0.6； 

（4）对数据样本集X分别用EM算法聚类4次，用k-means聚类3次，得到7组聚类类别标签； 

（5）对7组聚类结果求每个样本的聚类一致性值CI； 

（6）根据CI和thre，把数据样本集X划分成源域和目标域：CI≥0.6的样本为源域X^S，CI＜0.6为目标域X^T； 

（7）用基于GMM的EM算法对源域X^S进行参数估计，得到GMM的均值μ_js和协方差σ_js两个类参数，j表示类别； 

（8）用K近邻算法寻找目标域X^T的每个样本在源域中的K个近邻点，K＝7，并根据步骤（4）中获得的标签求7个近邻点的聚类一致性值CI′； 

（9）根据CI′并结合步骤（7）中获得的源域参数，计算目标域的均值μ_jnew和协方差σ_jnew两个新参数，即将源域参数迁移至目标域； 

（10）将所述两个新参数代入高斯分布的公式求得目标域每个样本的概率值，得到图像的最终划分结果。 

2.根据权利要求1所述的SAR图像分割方法，其中步骤（4）所述的对数据样本集X分别用EM算法聚类4次、用k-means聚类3次，得到7组聚类标签，按照如下步骤进行： 

（4a）数据集为X＝|{x_i|i＝1,2,...,n}，其中n为数据个数，x_i＝{x_i1,x_i2,...,x_im}^T，m为维数，首先随机生成一个归属矩阵z_ij： 

若z_ij＝1表示样本x_i属于第j类，z_ij＝0表示样本x_i不属于第j类； 

（4b）利用如下公式求GMM初始类混合权重α_j、类均值μ_j和类方差σ_j： 

α_j＝1C 

（4c）根据上述类参数求得任意一项样本x_i的概率为： 

其中表示样本x_i在第j类的高斯分布概率，m为维数； 

（4d）利用所求得的概率计算归属矩阵z_ij的期望值E_ij,即E-step； 

E_ij＝E(z_ij|x_i;α_j,μ_j,σ_j)＝α_jG(x_i;μ_j,σ_j)/p(x_i) 

其中E_ij＝E(z_ij|x_i;α_j,μ_j,σ_j)表示期望值，E_ij是在类混合权重α_j、类均值μ_j和类方差σ_j条件下求得的样本x_i的归属矩阵z_ij的期望值； 

（4e）用产生的期望值E_ij替换类均值μ_j和类方差σ_j公式中的z_ij，得到新的类均值μ_j′和类方差σ_j′，并用下式更新类混合权重得到α_j′： 

（4f）将归属矩阵z_ij的期望值E_ij最大化，即M-step： 

首先，用所述更新后的类均值μ_j′、类方差σ_j′和类混合权重α_j′重新进行步骤(4c)-(4e)得到归属矩阵z_ij的新期望E′_ij； 

然后，反复迭代步骤(4c)-(4e)直至|E_t+1′-E_t′|≤0.000001循环结束，E_t′表示第t次循环得到的z_ij的期望，E_t+1′表示第t+1次循环得到的z_ij的期望； 

（4g）利用z_ij的期望E′_ij更新归属矩阵z_ij，得到样本的最终划分： 

z_ij＝1当且仅当

z_ij＝0当且仅当

（4h）重复(4a)-(4g)四次，即可以得到用EM算法聚类的4组结果； 

（4i）对数据样本集X用k-means算法聚类3次，得到另外的3组聚类结果。 

3.根据权利要求1所述的SAR图像分割方法，其中步骤（5）所述的对7组聚类结果求样本聚类一致性值CI，根据下式计算： 

其中，H表示样本x_i的聚类次数，此处H＝7，π(x_i)表示样本x_i的聚类结果，L为每一个类别标号，Σδ(π(x_i),L)表示类别L在π(x_i)中的总个数，C是图像划分的类别数。 

4.根据权利要求1所述的SAR图像分割方法，其中步骤（9）所述计算目标域的均值μ_jnew和协方差σ_jnew两个新参数，按照如下公式进行， 

μ_jnew＝CI′μ_js

σ_jnew＝CI′σ_js

CI′是每个目标域样本在其源域样本中寻得的7个近邻点的聚类一致性值，μ_js是源域的类均值，σ_js是源域的类方差。 

Images