CN101751674A - 基于Graph-cut和广义高斯模型的遥感影像变化检测方法 - Google Patents

Abstract

本发明公开了一种基于Graph-cut和广义高斯模型的遥感影像变化检测方法，主要解决了现有差异图像分类时阈值难以确定以及比值法构造的差异图像难以分析的问题。其实现过程是：(1)采用对数比值法构造差异图像；(2)利用Graph-cut算法对差异图像初分类；(3)对初分类的结果采用FCM算法聚类；(4)利用EM算法估计广义高斯模型类别参数；(5)根据贝叶斯决策判断像素所属类别，得出变化检测结果。实验表明本发明具有检测精度高，虚检较少，适用性强的优点，可用于对遥感图像的灾情评估和土地利用。

Description

基于Graph-cut和广义高斯模型的遥感影像变化检测方法

技术领域

本发明属于图像处理技术领域，涉及遥感影像变化检测，具体地说是一种基于Graph-cut和冠以高斯模型的遥感影像变化检测方法。该方法可用于对遥感影像变化检测中的差异图像分类。

背景技术

遥感技术是从远距离感知目标反射或自身辐射的电磁波、可见光、红外线对目标进行探测和识别的技术，是获取地面信息的快速有效的手段之一。随着遥感技术的不断发展，利用同一地区多时相的遥感影像进行变化检测已经成为遥感技术的一个重要应用领域。遥感变化检测通常包含了四个方面的内容：判断是否发生了变化；确定发生变化的区域；鉴别变化的性质；评估变化的时间和空间分布模式。前两个方面是变化检测所要达到的基本目标，也是变化检测研究需要首先解决的问题。

现阶段国内外变化检测的研究蓬勃发展。变化检测广泛地应用于环境监测，土地利用/覆盖变化研究，森林植被变化分析，灾害评估，农业调查，城市规划布局分析，洪水监测和地图更新等领域。然而从整体上来说遥感影像的变化检测仍然存在不少困难。遥感影像的变化检测目前缺乏坚实有力的数学理论和模型支持。在进行变化检测时，一般都没有获取检测地区的变化地物的先验信息，只是根据图像上的信息来检测变化。绝大多数的变化检测方法对数据的统计分布要求比较严格，方法实施的时候需要事先知道数据的统计分布模型，同时，大多数变化检测方法对输入影像之间的辐射校正、几何配准等要求较高，对噪声等因素的影响比较敏感。

传统的遥感影像变化检测方法通常是直接比较同一位置不同时相的像元特征值来检测变化，通常采用数学变换的方式产生不同时相间的差异图像，再对差异图像进行阈值化处理，从中提取变化区域。这些方法存在着一些自身的缺陷，变化阈值的选择很大程度上影响变化检测的精度。因此对差异图像的阈值化处理提取变化区域成为变化检测的难题和关键问题。传统的代数运算类变化检测方法有差值法、比值法、图像回归法和变化矢量法。传统的差值法在构造差异图像时简单、直观，结果比较容易解译。但是不能提供变化类型信息；当两次成像条件(季节、太阳高度角、地表湿度等)不同时，也可能造成灰度差异，但它不一定能代表目标发生了变化。图像回归法由于在建立线性回归关系时考虑了不同时相图像中由于大气条件、季节、太阳角度等因素引起的灰度差异，所以可以消除这些因素对变化检测的影响。但是建立高精度的回归关系往往比较困难，计算量大。变化矢量法可以认为是图像差值法的扩展，可以提供变化类型信息，但是当两次成像条件不同时，也可能造成灰度差异。传统的比值法构造差异图像时通过除法运算可以消除一些由于太阳高度角、阴影和地形引起的乘性误差。但是生成的结果图像往往不服从正态分布，给差异图像的分析造成一定的困难。因此在比值法构造的差异图像中选择合适的模型分析差异图像是基于比值法构造差异图像需要解决的问题。

发明内容

本发明的目的在于：针对变化检测中运用传统方法构造的差异图像分类时阈值难以确定的问题，提出了一种基于Graph-cut和广义高斯模型的遥感影像变化检测方法。该方法用于差异图像分类，不需要阈值，解决了阈值选择的问题。并且选择广义高斯分布近似对数比值法构造的差异图像的类别概率分布函数。解决了比值法差异图像分析困难的问题。

本发明的技术方案是：将对数比值法构造的差异图像运用Graph-cut算法进行初始分类，将分类结果运用FCM算法聚为两类，然后采用EM算法估计广义高斯分布参数，优化分类结果。具体实现步骤如下：

1、一种基于Graph-cut和广义高斯模型的遥感影像变化检测方法，包括如下步骤：

(1)选取两个不同时间采集到的待检测遥感影像，分别表示为X₁和X₂，采用对数比值法构造差异图像为：

X_LR＝log(X₂/X₁)＝logX₂-logX₁；

(2)利用Graph-cut方法将差异图像分为N类：G₁，G₂，…，G_N，N＞2；

(3)运用模糊C均值FCM聚类算法将上述分类结果合并为变化类和非变化类，分别表示为ω_c和ω_u，该ω_c和ω_u服从广义高斯分布，即

$p\left(x_{i,j}\right|{\mathrm{\ω}}_i)=a_i{e}^{-{\left[b_i\right|x_{i,j}-m_i\left|\right]}^{{\mathrm{\β}}_i}},i\∈\{c,u\}$

其中x_i，j为差异图像的像素值，m_i为该类别的均值，c表示变化类，u表示非变化类， $a_i=\frac{b_i{\mathrm{\β}}_i}{2\mathrm{\Γ}(1/{\mathrm{\β}}_i)},$ $b_i=\frac{1}{{\mathrm{\σ}}_i}\sqrt{\frac{\mathrm{\Γ}(3/{\mathrm{\β}}_i)}{\mathrm{\Γ}(1/{\mathrm{\β}}_i)}},$ β_i为形状参数；

(4)采用最大数学期望EM算法估计ω_c和ω_u广义高斯模型参数m_u、σ_u ²、β_u、p(ω_u)、m_c、σ_c ²、β_c和p(ω_c)，其中mx为非变化类的均值，σ_u ²为非变化类的方差，β_u为非变化类的形状参数、p(ω_u)为非变概率密度函数，其中m_c为非变化类的均值，σ_c ²为非变化类的方差，β_c为非变化类的形状参数、p(ω_c)为非变概率密度函数；

(5)根据贝叶斯最小误差准则，针对图像中的每个像素点，分别计算p(ω_c)*p(x_i，j|ω_c)和p(ω_u)*p(x_i，j|ω_u)，根据计算结果判断像素点的变化类，如果p(ω_c)*p(x_i，j|ω_c)大于p(ω_u)*p(x_i，j|ω_u)，则该像素点判断为变化类，否则判为非变化类，最终得到变化检测结果。

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

1、本发明由于采用Graph-cut分类方法，不需要设置阈值，因而避免了由于阈值选择对变化检测精度的影响；

2、本发明由于采用广义高斯模型，相比于传统的高斯模型，能够更加精确地近似差异图像的类别概率分布函数，提高了分类精度。

附图说明

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

图2是本发明用于实验的数据集原始图像及变化参考图像；

图3是本发明利用对数比值法构造的差异图像；

图4是本发明采用Graph-cut算法对差异图像初分类得到的结果图；

图5是本发明将初分类结果采用FCM算法聚为两类的结果图；

图6是本发明采用广义高斯模型分类变化检测结果图；

图7是本发明中采用广义高斯模型分类结果和基于高斯模型的变化检测结果比较图；

图8是广义高斯分布函数与高斯分布函数的比较图；

图9是本发明与现有的二种方法得到的变化检测结果图。

具体实施方式

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

步骤1，设两幅配准的遥感图像X₁＝{X₁(i，j)，1≤i≤I，1≤j≤J}和X₂＝{X₂(i，j)，1≤i≤I，1≤j≤J}分别是t1和t2时间同一地区采集到的数据，t1时间图像如图2(a)所示，t2时间图像如图2(b)所示，设X_R是通过比值法逐像元比较得到的差异图像，对比值差异图像取对数变换，得到对数比值差异图像，如图3所示，

$X_{\mathrm{LR}}=\log X_R=\log \frac{X_2}{X_1}=\log X_2-\log X_1.$

步骤2，将差异图像X_LR看作一个无向图G：(V，E)，其中V表示图中所有顶点的集合，E表示图中顶点之间边缘的集合，对于图像中任意空间位置为(i，j)的像素点有一个相应的顶点v_i，j∈V，它与其相邻的空间位置为(r，s)的像素点之间有边缘e_{i，j；r，s}∈E，每个边缘e_{i，j；r，s}∈E有一个衡量基于像素点(i，j)和(r，s)之间相似性的权值w(e_{i，j；r，s})，按w(e_{i，j；r，s})值大小，采用Graph-cut方法将图像分成N类：G₁，G₂，…，G_N，N＞2，如图4所示。

步骤3，将上述Graph-cut方法得到的分类结果，采用模糊C均值聚类的方法，分为变化类ω_c和非变化类ω_u，如图5所示。

步骤4，采用以下步骤估计形状参数β_i：

4a)将β_i的取值范围限定为1.0～10.0，每隔为0.1计算一次广义高斯比率函数 $r\left({\mathrm{\β}}_i\right)=\frac{\mathrm{\Γ}(1/{\mathrm{\β}}_i)\mathrm{\Γ}(3/{\mathrm{\β}}_i)}{{\mathrm{\Γ}}^2(2/{\mathrm{\β}}_i)},$ 并建立r(β_i)与β_i一一对应的一个数据表；

4b)计算每个像素点像素值与均值差的数学期望，计算公式为：

$E\left(\right|x\left|\right|{\mathrm{\ω}}_i)=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}|x_i-m_i|$

其中x_i为像素值，m_i为均值，M为像素点的总数目；

4c)计算该类别所有像素点方差σ² _i与上述数学期望E(|x||ω_i)的比值

${\mathrm{\ρ}}_i=\frac{{{\mathrm{\σ}}^2}_i}{E\left(\right|x\left|\right|{\mathrm{\ω}}_i)};$

4d)在步骤中4a)建立的数据表中搜索，找到与ρ_i对应的β_i的值，作为形状参数β_i的估计值。

步骤5，假设差异图像的变化类ω_c和非变化类ω_u都服从广义高斯分布，则将ω_c和ω_u的概率分布函数表示为：

$p\left(x_{i,j}\right|{\mathrm{\ω}}_i)=a_i{e}^{-{\left[b_i\right|x_{i,j}-m_i\left|\right]}^{{\mathrm{\β}}_i}},i\∈\{c,u\}$

其中x_i，j为差异图像的像素值，m_i为该类别的均值，c表示变化类，u表示非变化类， $a_i=\frac{b_i{\mathrm{\β}}_i}{2\mathrm{\Γ}(1/{\mathrm{\β}}_i)},$ $b_i=\frac{1}{{\mathrm{\σ}}_i}\sqrt{\frac{\mathrm{\Γ}(3/{\mathrm{\β}}_i)}{\mathrm{\Γ}(1/{\mathrm{\β}}_i)}},$

根据形状参数β_i的估计值，可以计算得到ω_c和ω_u的概率分布函数。

步骤6，采用最大数学期望EM算法估计变化类ω_c和非变化类ω_u广义高斯模型参数p(ω_c)和p(ω_u)，其中p(ω_c)为ω_c概率密度函数，p(ω_u)为ω_u概率密度函数。

步骤7，利用步骤5中计算得到p(x_i，j|ω_c)和p(x_i，j|ω_u)，利用步骤6中估计参数得到p(ω_c)和p(ω_u)，根据贝叶斯最小误差准则，针对图像中的每个像素点，分别计算p(ω_c)*p(x_i，j|ω_c)和p(ω_u)*p(x_i，j|ω_u)，根据计算结果判断像素点的变化类，如果p(ω_c)*p(x_i，j|ω_c)大于p(ω_u)*p(x_i，j|ω_u)，则该像素点判断为变化类，否则判为非变化类，最终得到变化检测结果，如图6所示。

本发明的效果可以通过以下内容进行说明：

1.实验数据

图2(a)、2(b)是实验数据原图，图2(c)为参考变化图，图中白色区域表示变化的区域，其中，变化的像元数为4236，未变化像元数为153214。

2.对比实验

本发明中采用两种对比方案，第一种对比方案将本发明方法变化检测结果与采用高斯模型的变化检测结果比较，第二种对比方案将本发明方法变化检测结果与传统的阈值法变化检测结果和传统的FCM算法变化检测结果进行比较。

3.实验结果与分析

图7(a)是基于广义高斯模型GGM的变化检测结果，图7(b)是基于高斯模型GM的变化检测结果，可以看出广义高斯模型优于高斯模型的检测结果，GGM得到的检测结果虚检明显少于GM得到的检测结果，这两种模型的检测率比较见表1。

图8(a)是变化类概率密度函数，图8(b)是非变化类概率密度函数，从图8中可以看出，对于广义高斯密度函数非常逼近非变化类别概率密度函数，而高斯函数与非变化类别概率密度函数相距较远，因此相对于基于高斯模型的变化检测结果，基于广义高斯模型的变化检测结果较好。

图9(a)是阈值法检测结果，图9(b)是FCM方法检测结果，图9(c)是本发明方法检测结果，从图9中可以看出，相对于差异图像阈值分类方法和差异图像FCM分类方法，基于Graph-cut和广义高斯模型的变化检测方法减少了许多虚检，检测结果较好，这三种方法的检测率比较见表2。

综上，本发明提出的基于Graph-cut和广义高斯模型的变化检测方法能够达到较高的检测精度，比差异图像阈值分类方法、差异图像FCM分类方法和基于Graph-cut和高斯模型的变化检测方法具有更好的检测精度，虚检较少，得到的反映检测结果内部一致性的kappa系数较高。

表1基于广义高斯模型和高斯模型的检测结果比较

模型	检测数	虚检数	总错误数	检测率	总错误率	Kappa系数
							GGM	4015	122	343	94.78％	8.17％	95.79％
GM	4183	4570	4623	98.75％	53.46％	63.01％

表2本发明方法与阈值法和FCM算法分类差异图像的检测结果比较

方法	检测数	虚检数	总错误数	检测率	总错误率	Kappa系数
							阈值法	3960	395	671	93.48％	15.59％	91.97％
FCM	4017	633	852	94.83％	18.78％	90.38％
							本发明方法	4015	122	343	94.78％	8.17％	95.79％

Claims (2)

1.一种基于Graph-cut和广义高斯模型的遥感影像变化检测方法，包括如下步骤：

(1)选取两个不同时间采集到的待检测遥感影像，分别表示为X₁和X₂，采用对数比值法构造差异图像为：

X_LR＝log(X₂/X₁)＝logX₂-logX₁；

(2)利用Graph-cut方法将差异图像分为N类：G₁，G₂，…，G_N，N＞2；

(3)运用模糊C均值FCM聚类算法将上述分类结果合并为变化类和非变化类，分别表示为ω_c和ω_u，该ω_c和ω_u服从广义高斯分布，即

$p\left(x_{i,j}\right|{\mathrm{\ω}}_i)=a_i{e}^{-{\left[b_i\right|x_{i,j}-m_i\left|\right]}^{{\mathrm{\β}}_i}}i\∈\{c,u\}$

其中x_i，j为差异图像的像素值，m_i为该类别的均值，c表示变化类，u表示非变化类， $a_i=\frac{b_i{\mathrm{\β}}_i}{2\mathrm{\Γ}(1/{\mathrm{\β}}_i)},$ $b_i=\frac{1}{{\mathrm{\σ}}_i}\sqrt{\frac{\mathrm{\Γ}(3/{\mathrm{\β}}_i)}{\mathrm{\Γ}(1/{\mathrm{\β}}_i)}},$ β_i为形状参数；

(4)采用最大数学期望EM算法估计ω_c和ω_u广义高斯模型参数m_u、σ_u ²、β_u、p(ω_u)、m_c、σ_c ²、β_c和p(ω_c)，其中m_u为非变化类的均值，σ_u ²为非变化类的方差，β_u为非变化类的形状参数、p(ω_u)为非变概率密度函数，其中m_c为非变化类的均值，σ_c ²为非变化类的方差，β_c为非变化类的形状参数、p(ω_c)为非变概率密度函数；

(5)根据贝叶斯最小误差准则，针对图像中的每个像素点，分别计算p(ω_c)*p(x_i，j|ω_c)和p(ω_u)*p(x_i，j|ω_u)，根据计算结果判断像素点的变化类，如果p(ω_c)*p(x_i,j|ω_c)大于p(ω_u)*p(x_i，j|ω_u)，则该像素点判断为变化类，否则判为非变化类，最终得到变化检测结果。

2.根据权利要求1所述的基于Graph-cut和广义高斯模型的遥感影像变化检测方法，其中步骤(3)所述的形状参数β_i，按如下步骤估计：

2a)将β_i的取值范围限定为1.0～10.0，每隔为0.1计算一次广义高斯比率函数 $r\left({\mathrm{\β}}_i\right)=\frac{\mathrm{\Γ}(1/{\mathrm{\β}}_i)\mathrm{\Γ}(3/{\mathrm{\β}}_i)}{{\mathrm{\Γ}}^2(2/{\mathrm{\β}}_i)},$ 并建立r(β_i)与β_i一一对应的一个数据表；

2b)计算每个像素点像素值与均值差的数学期望，计算公式为：

$E\left(\right|x\left|\right|{\mathrm{\ω}}_i)=\frac{1}{M}\underset{i=1}{\overset{M}{\mathrm{\Σ}}}|x_i-m_i|$

其中x_i为像素值，m_i为均值，M为像素点的总数目；

2c)计算该类别所有像素点方差σ² _i与上述数学期望E(|x||ω_i)的比值

${\mathrm{\ρ}}_i=\frac{{\mathrm{\σ}}_i^2}{E\left(\right|x\left|\right|{\mathrm{\ω}}_i)};$

2d)在步骤中2a)建立的数据表中搜索，找到与ρ_i对应的β_i的值，作为形状参数β_i的估计值。

Images