CN103700109B - 基于多目标优化moea/d和模糊聚类的sar图像变化检测方法 - Google Patents

Abstract

本发明公开一种基于多目标进化算法MOEA/D和模糊聚类的SAR图像变化检测方法，主要是通过使用多目标优化的方法来解决SAR图像变化检测中细节保持和噪声去除这两个目标之间的权衡问题。其步骤为：(1)对待检测图像使用对数比方法生成差异图；(2)对差异图进行滤波，得到去噪后的差异图；(3)根据细节保持和噪声去除两个目标确定出两个目标函数，并组合成多目标优化问题；(4)使用MOEA/D算法求得多目标问题的Pareto前端和对应的结果图；(5)根据需要，从所有结果中选择合适的变化检测结果图。本发明最大的优势在于，相对于其它变化检测算法只得到一个解的情况，本发明得到的是一个最优解集，用户可以根据自己对细节保持和噪声去除的偏重程度，从中选择更加合适的解。

Description

基于多目标优化MOEA/D和模糊聚类的SAR图像变化检测方法

技术领域

本发明属于图像处理和进化计算技术领域，涉及多目标优化算法MOEA/D和SAR图像变化检测，特别是基于多目标优化MOEA/D和模糊聚类的SAR图像变化检测方法，可用于环境监测、农业调查、城市研究、救灾工作等SAR图像变化检测相关领域中。

背景技术

近年来，合成孔径雷达技术得以迅速发展，星载合成孔径雷达系统（SyntheticApertureRadar，SAR）已经对地球表面观测了数年时间，获得了大量的多时相地面观测数据。很多遥感研究都试图开发出能够很好利用这些信息的技术，包括目标提取、地物分类、边缘检测、干涉测量、变化检测等，其中研究最广泛的就是变化检测技术。

SAR图像变化检测通过对不同时期SAR图像的比较分析，根据图像之间的差异来得到所需要的地物变化信息。近年来，由于SAR图像的变化检测可以广泛用于环境监测、农作物测量、城市研究、森林监测等领域，人们对变化检测技术的兴趣大大增加。因为合成孔径雷达（SAR）具有高分辨率、全天候，全天时的特点，是很好的变化检测信息源，所以研究SAR图像变化检测技术有着非常广阔的应用前景。

随着SAR成像技术的日益成熟、成像质量的逐渐提高和SAR图像的分辨率不断增强，基于SAR图像的变化检测研究受到了越来越多的关注，各种方法相继被提出，用于提高SAR图像变化检测的性能。在遥感应用中有许多普遍使用的变化检测技术，现行研究应用的变化检测技术主要有以下五种，图像差值法、图像比值法、图像回归法、主分量分析法、变化向量法、分类比较法以及统计测试法等。从基本上，SAR图像变化检测方法也可以分为图像阈值法和图像分类法两大类。在通过对比SAR图像来辨识图像的“变化”和“未变化”区域时，人们通常将变化检测问题转化为图像的二值分类问题。

然而，由于SAR图像本身存在大量的相干斑噪声，许多传统的变化检测方法又对噪声过于敏感，导致误检和漏检的现象比较严重，结果表现为离散的检测点比较多，从而影响了变化检测的精度，如广泛使用的基于FCM聚类的方法，它虽能较好的保留细节信息，但却对噪声的存在非常敏感。虽然可以通过相干斑抑制算法来降低相干斑噪声的影响，如多视处理、空间滤波算法等，但这些处理都将会导致图像细节信息的丢失。

针对FCM算法的噪声敏感问题，Ahmed等人提出的FCM_S算法，该算法在FCM算法的目标函数的基础上加入了空间域信息，从而实现在一定程度上抑制噪声的目标，参见M.Ahmed,S.Yamany,N.Mohamed,A.Farag,andT.Moriarty,“AmodifiedfuzzyC-meansalgorithmforbiasfieldestimationandsegmentationofMRIdata,”IEEETrans.Med.Imag.,vol.21,pp.193–199,2002。但是该算法存在一个很大的弊端，就是它需要人为选择一个参数来平衡图像的细节保持能力和噪声去除能力，这在没有对所存在噪声的先验知识的情况下是非常困难的。

正是由于上述这些问题的存在，使得在SAR图像变化检测问题中细节保持和噪声去除之间的均衡成为了研究中亟待解决的问题。因此，研究一种行之有效的、能够根据不同情况兼顾图像细节保持和噪声去除的变化检测算法成为当务之急。

发明内容

本发明的目的在于克服上述已有技术的不足，提出一种基于多目标优化MOEA/D和模糊聚类的SAR图像变化检测方法，通过将细节保持和噪声去除作为两个单独的目标进行多目标优化，从而解决SAR图像变化检测中细节保持和噪声去除这两个本质目标之间的权衡问题。

本发明的技术方案是，基于多目标优化MOEA/D和模糊聚类的SAR图像变化检测方法，其特征是：包括如下步骤：

步骤101：对经过配准、几何校正等预处理后得到的同一地区不同时刻的两幅图像X₁、X₂使用对数比方法生成差异图X_l，该差异图是根据如下对数比方法的计算公式得到的：

$X_l=\left|\log \frac{X_2}{X_1}\right|=|\log X_2-\log X_1|$

其中X_l为计算得到的差异图，X₁和X₂分别为同一地区不同时刻的两幅预处理后的图像，log为自然对数运算符；

步骤102：使用均值滤波方法对步骤101中产生的差异图X_l进行处理，得到去噪后的差异图该差异图按如下公式计算得到：

$\stackrel{\‾}{x_i}=\frac{1}{N_R}\underset{i\∈N_i}{\mathrm{\Σ}}{x}_i$

其中为得到的去噪后的差异图的第i个像素的灰度值，xi为差异图的第i个像素的灰度值，N_R为第i个像素的邻域所包含的像素个数，本发明中选用的是3x3的邻域窗口，即N_R为9，N_i为落在第i个像素邻域内的像素集合；

步骤103：根据细节保持和噪声去除这两个目标确定出两个目标函数f₁、f₂，并将其组合成多目标优化问题；

步骤104：使用多目标优化算法MOEA/D求解由步骤103得到的多目标问题，得到符合条件的Pareto前端，然后根据每一个Pareto最优解对应的隶属度矩阵求得对应的变化检测结果图；这一步骤中使用的MOEA/D算法的思想是把多目标优化问题分解成若干个单目标优化子问题，每个子问题在当前的种群中都能找到一个局部最优解；用权重向量的欧氏距离来描述子目标中邻居之间的远近程度，而这个权重向量就是所有子目标的聚合系数；然后通过进化过程，同时求解这些子问题，由于两个相邻子问题的优化解理论上非常相似，所以在MOEA/D算法中每个子问题均可借助于其相邻子问题的优化信息，从而可以得到更加可靠的Pareto前端；

步骤105：根据需要，从Pareto前端中选择合适的解，得到对应的SAR图像变化检测的结果图，其选择原则是：若重点在于图像细节的保持，则选择使第一个目标函数取得最小值的解；若重点在于图像中噪声的去除，则选择使第二个目标函数取得最小值的解，否则，根据对两个目标的偏重程度差异选择一个折衷解。

所述的步骤103，包括如下步骤：

步骤201：从图像细节保持的角度出发，选用如下函数作为第一个目标函数f₁：

$f_1(v_1,v_2)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\underset{j=1}{\overset{c}{\mathrm{\Σ}}}{u}_{\mathrm{ji}}^m{\left|\right|x_i-v_j\left|\right|}^2---\left(1\right)$

其中N为图像总的像素个数，c为聚类数目，u_ji∈(0,1)为第i个像素相对于第j类的隶属度，m为大于1的加权指数，x_i为差异图的第i个像素的灰度值，v_j为第j个聚类中心，由随机初始化得到，然后根据如下公式得到隶属度矩阵：

$u_{\mathrm{ij}}=\frac{1}{\underset{k=1}{\overset{c}{\mathrm{\Σ}}}{\left(\frac{\left|\right|x_i-v_j\left|\right|}{\left|\right|x_i-v_k\left|\right|}\right)}^2}$

因为本发明是用于解决变化检测问题，所以聚类个数选为2，即分为变化类和未变化类，加权指数m也取为2；

步骤202：从图像噪声去除的角度出发，选用如下函数作为第二个目标函数f₂：

$f_2(v_1,v_2)=\underset{i=11}{\overset{N}{\mathrm{\Σ}}}\underset{j=1}{\overset{c}{\mathrm{\Σ}}}{u}_{\mathrm{ji}}^m{\left|\right|\stackrel{\‾}{x_i}-v_j\left|\right|}^2---\left(2\right)$

它与第一个目标函数(1)的区别在于为对差异图进行均值滤波后得到的图像第i个像素的灰度值；

步骤203：由以上两个目标函数(1)和(2)组合成为一个具有两个决策变量和两个目标变量的多目标优化问题：

$\left\{\begin{array}{c}\min F(v_1,v_2)={(f_1,f_2)}^T\\ s.t.{(v_1,v_2)}^T\∈\mathrm{\Ω}\end{array}\right.---\left(5\right)$

其中变化类和未变化类的聚类中心(v₁,v₂)分别为两个决策变量，它们组成一个决策矢量v＝(v₁,v₂)^T，Ω为决策空间，这样就把问题转化成了对这个多目标问题(3)的求解。

本发明的优点是：1）本发明通过将图像细节的保持和噪声的去除作为多目标优化问题中的两个目标进行优化，解决了现有变化检测算法中普遍存在的图像细节保持和噪声去除间的均衡问题；2）在本发明中，进化算法作为一种并行的搜索技术被加入到聚类过程中，克服了传统聚类方法用于变化检测时存在的对初始聚类中心敏感的缺点，提高了其收敛到全局最优解的概率，提升了变化检测的稳定性；3）本发明由于采用基于分解的MOEA/D多目标进化算法进行求解，将多目标优化问题分解成若干个单目标优化子问题，降低了算法的时间复杂度，并且相比于NSGA-II以及MOGLS等多目标算法，MOEA/D算法的到的解集更加逼近真实的Pareto前端。

附图说明

图1是本发明实现步骤的主流程框图；

图2是所述步骤103的流程图；

图3是第一组仿真实验图，其中图3a）和图3b）分别拍摄于1999.04和1999.05，大小均为301×301，图3c）为参考图；

图4是由步骤101中的对数比方法产生的对应于第一组实验的两幅差异图，其中图4a)是未经过滤波的，4b)是经过均值滤波的差异图；

图5是从所有结果中随机选择的对应于第一组实验图的三幅变化检测结果图；

图6第二组仿真实验图，其中图6a)和图6b）分别拍摄于1997.05和1997.08，大小均为290×350，图6c）为参考图；

图7是由步骤101中的对数比方法产生的对应于第二组实验的两幅差异图，其中图7a)是未经过滤波的，7b)是经过均值滤波的差异图；

图8是从所有结果中随机选择的对应于第二组实验图的三幅变化检测结果图；

图9第三组仿真实验图，其中图9a)和图9b）分别拍摄于2008年和2009年，大小均为306×291，图9c）为参考图；

图10是由步骤101中的对数比方法产生的对应于第三组实验的两幅差异图，其中图10a)是未经过滤波的，10b)是经过均值滤波的差异图；

图11是从所有结果中随机选择的对应于第三组实验图的三幅变化检测结果图；

图12是从Pareto前端选出的在细节保持和噪声去除能力方面具有不同侧重的一部分变化检测结果图。

具体实施方式

如图1所示，

主流程图步骤特征是：

步骤101：对经过配准、几何校正等预处理后得到的同一地区不同时刻的两幅图像X₁、X₂使用对数比方法生成差异图X_l，该差异图是根据如下对数比方法的计算公式得到的：

$X_l=\left|\log \frac{X_2}{X_1}\right|=|\log X_2-\log X_1|$

其中X_l为计算得到的差异图，X₁和X₂分别为同一地区不同时刻的两幅预处理后的图像，log为自然对数运算符；

步骤102：使用均值滤波方法对步骤101中产生的差异图X_l进行处理，得到去噪后的差异图该差异图按如下公式计算得到：

$\stackrel{\‾}{x_i}=\frac{1}{N_R}\underset{i\∈N_i}{\mathrm{\Σ}}{x}_i$

其中为得到的去噪后的差异图的第i个像素的灰度值，xi为差异图的第i个像素的灰度值，N_R为第i个像素的邻域所包含的像素个数，本发明中选用的是3x3的邻域窗口，即N_R为9，N_i为落在第i个像素邻域内的像素集合；

步骤103：根据细节保持和噪声去除这两个目标确定出两个目标函数f₁、f₂，并将其组合成多目标优化问题；

步骤104：使用多目标优化算法MOEA/D求解由步骤103得到的多目标问题，得到符合条件的Pareto前端，然后根据每一个Pareto最优解对应的隶属度矩阵求得对应的变化检测结果图；这一步骤中使用的MOEA/D算法的思想是把多目标优化问题分解成若干个单目标优化子问题，每个子问题在当前的种群中都能找到一个局部最优解；用权重向量的欧氏距离来描述子目标中邻居之间的远近程度，而这个权重向量就是所有子目标的聚合系数；然后通过进化过程，同时求解这些子问题，由于两个相邻子问题的优化解理论上非常相似，所以在MOEA/D算法中每个子问题均可借助于其相邻子问题的优化信息，从而可以得到更加可靠的Pareto前端；

步骤105：根据需要，从Pareto前端中选择合适的解，得到对应的SAR图像变化检测的结果图，其选择原则是：若重点在于图像细节的保持，则选择使第一个目标函数取得最小值的解；若重点在于图像中噪声的去除，则选择使第二个目标函数取得最小值的解，否则，根据对两个目标的偏重程度差异选择一个折衷解。

如图2所示，

所述的步骤103，包括如下步骤：

步骤201：从图像细节保持的角度出发，选用如下函数作为第一个目标函数f₁：

$f_1(v_1,v_2)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\underset{j=1}{\overset{c}{\mathrm{\Σ}}}{u}_{\mathrm{ji}}^m{\left|\right|x_i-v_j\left|\right|}^2---\left(1\right)$

其中N为图像总的像素个数，c为聚类数目，u_ji∈(0,1)为第i个像素相对于第j类的隶属度，m为大于1的加权指数，x_i为差异图的第i个像素的灰度值，v_j为第j个聚类中心，由随机初始化得到，然后根据如下公式得到隶属度矩阵：

$u_{\mathrm{ij}}=\frac{1}{\underset{k=1}{\overset{c}{\mathrm{\Σ}}}{\left(\frac{\left|\right|x_i-v_j\left|\right|}{\left|\right|x_i-v_k\left|\right|}\right)}^2}$

因为本发明是用于解决变化检测问题，所以聚类个数选为2，即分为变化类和未变化类，加权指数m也取为2；

步骤202：从图像噪声去除的角度出发，选用如下函数作为第二个目标函数f₂：

$f_2(v_1,v_2)=\underset{i=11}{\overset{N}{\mathrm{\Σ}}}\underset{j=1}{\overset{c}{\mathrm{\Σ}}}{u}_{\mathrm{ji}}^m{\left|\right|\stackrel{\‾}{x_i}-v_j\left|\right|}^2---\left(2\right)$

它与第一个目标函数(1)的区别在于为对差异图进行均值滤波后得到的图像第i个像素的灰度值；

步骤203：由以上两个目标函数(1)和(2)组合成为一个具有两个决策变量和两个目标变量的多目标优化问题：

$\left\{\begin{array}{c}\min F(v_1,v_2)={(f_1,f_2)}^T\\ s.t.{(v_1,v_2)}^T\∈\mathrm{\Ω}\end{array}\right.---\left(5\right)$

其中变化类和未变化类的聚类中心(v₁,v₂)分别为两个决策变量，它们组成一个决策矢量v＝(v₁,v₂)^T，Ω为决策空间，这样就把问题转化成了对这个多目标问题(3)的求解。

本发明的效果可通过以下仿真进一步说明：

1.仿真条件

本实例是在Intel(R)Core(TM)2CPU1.86GHz1.97G内存WindowsXP系统下，Matlabr2012a运行平台上，完成的仿真实验。

2.仿真实验内容

A.将本发明应用在如图3所示两幅Bern城市的不同时刻的SAR图像上，其大小均为301×301，图3a)和图3b)的拍摄时间分别为1999.04和1999.05，用步骤101中的对数比方法产生差异图，如图4a)所示，并用步骤102中的均值滤波方法产生去噪后的差异图，如图4b)所示；

B.将本发明应用在如图6所示的Ottawa地区遭遇水灾前后的两幅SAR图像上，其大小均为290×350，图6a)和图6b)的拍摄时间分别为1997.05和1997.08，用步骤101中的对数比方法产生差异图，如图7a)所示，并用步骤102中的均值滤波方法产生去噪后的差异图，如图7b)所示；

C.将本发明应用在如图9所示的黄河部分区域的两幅不同时刻的SAR图像上，其大小均为306×291，由于原始黄河区域图像过于庞大，因此本发明中的这两幅图9a)和9b)是分别从2008和2009年拍摄的两幅黄河区域的大图上截下来的相同部分，用步骤101中的对数比方法产生差异图，如图10a)所示，并用步骤102中的均值滤波方法产生去噪后的差异图，如图10b)所示。

3.仿真实验结果及分析

从三组实验对应的变化检测结果图5、8、11可以看出，本发明得到的仿真实验结果具有较好的主观视觉效果，并且为进一步验证本发明的良好效果，在实验过程中还对所有Pareto占优解产生出的结果进行了变化检测指标评价。对于第一组实验图像的正确率在99.63%到99.68%之间，kapaa系数在0.8673到0.8708之间；对于第二组实验图像的正确率在97.64到98.99%之间，kappa系数在0.9488到0.9618之间；对于第三组实验图像的正确率在97.43%到98.30%之间，而kappa系数介于0.8455和0.8645之间。通过这些指标可以定量的看出本发明应用于SAR图像变化检测产生了较好的效果。

从图12可以看出，本发明对应于不同的Pareto占优解所产生的一组结果中在图像细节保持和噪声去除方面均有各自的侧重，或者是细节保持效果更好一些，或者是噪声去除效果更好一些，或者是两者效果均衡并重。

由以上的仿真实验可以说明，针对SAR图像变化检测，本发明存在一定的优势，克服了现有方法在SAR图像变化检测上的细节保持能力和噪声去除能力间的权衡问题。它避免用于平衡图像细节保持能力和噪声去除能力的人工参数的选择，在完全没有牺牲时间复杂度的前提下，对于没有对所存在噪声的先验知识的情况，可以得到一组检测结果均较好、同时又在细节保持和噪声去除效果上有所侧重的变化检测结果，从而使得用户可以根据不同情况不同需要选择更加合适的结果。

综上所述，本发明能够更加有效的应用于SAR图像变化检测中。

Claims (2)

1.基于多目标优化MOEA/D和模糊聚类的SAR图像变化检测方法，其特征是：包括如下步骤：

步骤101：对经过配准、几何校正预处理后得到的同一地区不同时刻的两幅图像X₁、X₂使用对数比方法生成差异图X_l，该差异图是根据如下对数比方法的计算公式得到的：

$X_l=\left|log\frac{X_2}{X_1}\right|=\left|\log X_2-\log X_1\right|$

其中X_l为计算得到的差异图，X₁和X₂分别为同一地区不同时刻的两幅预处理后的图像，log为自然对数运算符；

步骤102：使用均值滤波方法对步骤101中产生的差异图X_l进行处理，得到去噪后的差异图该差异图按如下公式计算得到：

$\stackrel{\‾}{x_i}=\frac{1}{N_R}\underset{i\∈N_i}{\mathrm{\Σ}}{x}_i$

其中为得到的去噪后的差异图的第i个像素的灰度值，x_i为差异图的第i个像素的灰度值，N_R为第i个像素的邻域所包含的像素个数，本发明中选用的是3x3的邻域窗口，即N_R为9，N_i为落在第i个像素邻域内的像素集合；

步骤103：根据细节保持和噪声去除这两个目标确定出两个目标函数f₁、f₂，并将其组合成多目标优化问题；

步骤104：使用多目标优化算法MOEA/D求解由步骤103得到的多目标问题，得到符合条件的Pareto前端，然后根据每一个Pareto最优解对应的隶属度矩阵求得对应的变化检测结果图；这一步骤中使用的MOEA/D算法的思想是把多目标优化问题分解成若干个单目标优化子问题，每个子问题在当前的种群中都能找到一个局部最优解；用权重向量的欧氏距离来描述子目标中邻居之间的远近程度，而这个权重向量就是所有子目标的聚合系数；然后通过进化过程，同时求解这些子问题，在MOEA/D算法中每个子问题均可借助于其相邻子问题的优化信息，从而可以得到更加可靠的Pareto前端；

步骤105：根据需要，从Pareto前端中选择合适的解，得到对应的SAR图像变化检测的结果图，其选择原则是：若重点在于图像细节的保持，则选择使第一个目标函数取得最小值的解；若重点在于图像中噪声的去除，则选择使第二个目标函数取得最小值的解，否则，根据对两个目标的偏重程度差异选择一个折衷解。

2.根据权利要求1所述的基于多目标优化MOEA/D和模糊聚类的SAR图像变化检测方法，其特征是：所述的步骤103，包括如下步骤：

步骤201：从图像细节保持的角度出发，选用如下函数作为第一个目标函数f₁：

$f_1(v_1,v_2)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\underset{j=1}{\overset{c}{\mathrm{\Σ}}}{u}_{ji}^m\left|\right|x_i-v_j|{|}^2---\left(1\right)$

其中N为图像总的像素个数，c为聚类数目，u_ji∈(0,1)为第i个像素相对于第j类的隶属度，m为大于1的加权指数，x_i为差异图的第i个像素的灰度值，v_j为第j个聚类中心，由随机初始化得到，然后根据如下公式得到隶属度矩阵：

$u_{ji}=\frac{1}{\underset{k=1}{\overset{c}{\mathrm{\Σ}}}{\left(\frac{\left|\right|x_i-v_j\left|\right|}{\left|\right|x_i-v_k\left|\right|}\right)}^2}$

因为本发明是用于解决变化检测问题，所以聚类个数选为2，即分为变化类和未变化类，加权指数m也取为2；

步骤202：从图像噪声去除的角度出发，选用如下函数作为第二个目标函数f₂：

$f_2(v_1,v_2)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\underset{j=1}{\overset{c}{\mathrm{\Σ}}}{u}_{ji}^m\left|\right|\stackrel{\‾}{x_i}-v_j|{|}^2---\left(2\right)$

它与第一个目标函数(1)的区别在于为对差异图进行均值滤波后得到的图像第i个像素的灰度值；

步骤203：由以上两个目标函数(1)和(2)组合成为一个具有两个决策变量和两个目标变量的多目标优化问题：

$\left\{\begin{array}{c}\min F(v_1,v_2)={(f_1,f_2)}^T\\ s.t.{(v_1,v_2)}^T\∈\mathrm{\Ω}\end{array}\right.---\left(3\right)$

其中变化类和未变化类的聚类中心(v₁,v₂)分别为两个决策变量，它们组成一个决策矢量v＝(v₁,v₂)^T，Ω为决策空间，这样就把问题转化成了对这个多目标问题(3)的求解。