遥感图像变化检测是通过对比同一地区不同时相的遥感图像数据，提取和描述感兴趣物体或现象随时间变化的特征，并定量分析和确定其变化的理论和技术^[1]。随着航天、遥感和信息技术的发展，变化检测已成为遥感领域研究的热点之一，广泛应用于国民经济和国防建设领域，如地理空间信息更新、城市用地变化、资源与环境检测、地震洪水森林大火等自然灾害监测等^[2]。

遥感图像变化检测方法有很多种，其中的无监督方法不需要施加任何先验知识，可以直接进行比较，是一种应用较为广泛的变化检测方法^[3]。而模糊C均值(fuzzy c-means, FCM)是应用最为广泛的无监督遥感图像变化检测方法之一，但该方法也存在着对孤立点敏感、抗噪性能极差等缺点^[4]。针对FCM方法的局限性，国内外很多学者提出了多种改进方法。赵磊等^[5]提出利用邻域像素到聚类中心的距离来调整隶属度，然而其邻域像素加权采用的是均匀的权值，忽略了像素之间位置的差异。Krishna等人^[6]提出了混合遗传FCM的变化检测方法，利用遗传算法来最小化FCM的目标函数，但是参数的选择采用的是经验值。文献[7]中的多时相遥感图像变化检测方法，采用FCM进行变化区域检测时并未考虑邻域像素的影响，易受噪声的影响。Bovolo等^[8]提出了一种基于多尺度分解和可靠尺度融合的SAR图像变化检测方法，采用对数比值法构造差异图，再对差异图进行二维的平稳小波分解，并通过可靠尺度融合的方法得到变化检测结果图，但它只对小波分解中各个尺度的低频信息进行处理，忽略了高频信息的影响。Li等人^[9]提出了基于Gabor小波变换、马尔科夫随机场与FCM相结合的遥感图像变化检测方法，具有较好的检测精度但易受噪声的影响且算法执行效率较低。Harshada等人^[10]提出了DWT和K-means相结合的变化检测方法，该方法的DWT融合规则过于简单，并不能很好地突出变化信息，此外K-means属于硬划分很容易受到噪声的干扰。慕彩虹等人^[11]提出了小波融合和PCA-核模糊聚类的遥感图像变化检测方法，该方法首先采用小波融合方式获取差异图像，再求取出不同小块图像的特征向量，最后采用基于核的模糊C均值提取变化区域，但需要反复调节核函数的参数才能取得较好的检测结果。

针对已有方法的不足，本文提出了一种基于离散小波变换和邻域信息FCM的遥感图像变化检测方法。该方法主要包括2个部分：差异图像获取和变化区域检测。首先，分别计算出已配准和校正后的2幅图像的差值图像和比值图像，把区域强度和区域能量信息作为融合规则对小波变换后的差值图像和比值图像的低频小波系数和高频小波系数进行融合，再利用小波逆变换获取差异图像；采用结合邻域信息和邻域灰度差值信息的改进FCM算法从差异图像中提取变化信息，在最大限度提高算法抗噪性能的同时得到精确的遥感图像变化结果。

1 差异图像获取

传统差异图的构造方法主要有图像差值法和图像比值法。图像差值法计算简单，容易受到图像成像质量、光谱特性的影响而产生“伪变化”信息^[12]。比如，1幅遥感图像的灰度值从0~50与205~255的变化，虽然其差值相等，但反应的真实地面变化却不相同。图像比值法可以在增强变化信息的同时减少外界条件对图像造成的影响，但是有时候会过于夸大变化的部分^[13]。比如，图像灰度值从200~20与从20~10的变化，比值法就不能很好地区分出变化的信息。

为了得到既能增强变化区域信息特征又能抑制背景信息的差异图像，提出了基于离散小波变换的差异图像构造方法，通过对低频小波系数和高频小波系数采用基于区域强度信息和区域能量信息的融合标准，得到最优的差异图像。图 1所示为本文所提出的基于离散小波变换获取差异图像的流程图。

1.1 差值与比值图像获取

采用公式(1)和(2)对已校正并配准的两幅遥感图像分别做差值和比值运算：

(1)

(2)

式中, F_t1(x)和F_t2(x)分别为t₁和t₂时刻的2幅影像, D_A(x)和D_B(x)分别表示最终所得的差值图像和比值图像。

1.2 基于离散小波变换的差异图像获取

采用公式(1)和(2)分别求出2幅遥感图像的差值图和比值图后; 再对所获取的差值图和比值图像进行归一化处理, 并对其进行1层离散小波变换, 得到1个低频带和3个高频带^[14-15]; 然后选用不同的融合规则对低频带系数和高频带系数进行融合。

1.2.1 低频带小波系数融合规则

文中对低频带小波系数采用区域强度作为融合的标准, 其定义为:

(3)

式中, I_{j, F}(x, y)表示图像F在2^-j分辨率下, 以(x, y)为中心位置的局部区域强度, w为窗口大小为N*N的掩膜算子(本文选择的窗口为w=[0, 1, 0;1, 4, 1;0, 1, 0]/8), A_{j, F}(x, y)为图像在该位置的低频带小波系数。

根据(3)式分别求取差值图像和比值图像的局部区域强度I_{j, A}(x, y)和I_{j, B}(x, y), 定义区域强度比R_{j, AB}(x, y)作为融合测度, 其表达式为:

(4)

为了增强图像细节信息, 选择强度大的系数作为融合系数, 反之采用算术平均的形式求出融合系数。因此, 差异图像的低频小波系数的选择因子为:

(5)

式中, T∈(0, 1), 一般取值区间为(0.5, 1), 本文选择的阈值为0.75, 则最终所得到的融合后差异图像的低频带小波系数为:

(6)

2.2.2 高频带小波系数融合规则

高频分量表征图像的细节信息, 而区域能量作为图像的特征之一, 反映了1幅图像的细节信息的丰富程度, 其表达式为:

(7)

式中, ε=LH, HL, HH, LH, HL, HH分别对应垂直方向, 水平方向和对角线方向的高频带, w(m, n)=[1, 2, 1;2, 3, 2;1, 2, 1]。

最终所得的高频带小波融合系数为:

(8)

对融合后的低频系数和高频系数进行离散小波逆变换^[16], 即可得到最终的差异图像。

2 变化区域检测

获取差异图像之后, 就要把差异图像中的像素点分为变化类和非变化类, 以获得最终的变化区域。本文提出了综合考虑邻域信息的FCM的变化检测方法, 把邻域与中心像素的距离信息和灰度差值信息引入到FCM的目标函数中, 使邻域像素对中心像素的影响程度随着邻域像素与中心像素的空间欧氏距离的增大而变小, 随着邻域内像素点与中心像素点的灰度值差异增大而增大, 以提高算法抑制噪声的能力, 避免误分类进而提高检测精度。

1) 传统FCM算法

模糊C均值聚类用隶属度来确定每个样本属于某种类别的程度的一种算法, 并通过迭代计算的方式求出其目标函数最小的1组样本隶属度^[17]。FCM算法的目标函数可以表示为:

(9)

式中, n表示待分类的样本数目, c表示类别数目, 隶属度值u_ik表示第k个样本在第i类上的隶属度, d_ik²=‖x_k-v_i‖²表示第i个聚类中心与第k个样本之间的欧氏距离, m∈(1, ∞)为模糊加权指数。

传统的FCM算法应用于遥感图像变化中存在着以下不足^[17], 文中做出如下改进:①根据小波融合获取的差异图像特点, 灰度值最高和最低的点, 是有最大概率发生变化和非变化的点, 因此可以选此类点为初始的聚类中心; ②对于噪声敏感的问题, 文中提出了结合空间距离信息和邻域灰度差值信息的FCM方法。

2) 邻域FCM方法

在遥感图像变化检测中, 确定一个点是否变化, 与它所处的空间位置有极大的相关性, 如果一个孤立点邻域全是确定变化的, 则该点发生变化的概率很大, 而且一个孤立点与其邻域像素的灰度差值较小, 则该点发生变化的概率也很大。本文通过将均匀的邻域信息加权值变为综合考虑空间邻域信息和邻域灰度差值信息的加权值, 使得邻域像素对中心像素的影响程度随着邻域像素与中心像素的空间欧氏距离的增大而变小, 随着邻域内像素点与中心像素点的灰度值差异的增大而增大。则结合邻域信息和邻域灰度差值信息的欧氏距离为:

(10)

(11)

式中，d_kl^N_k表示落入到所选取的N*N窗口(本文选择5*5大小的窗口, 中心像素点为k)的像素点l到中心像素点的欧氏距离, s_kl表示像素点l与中心像素点的灰度差值, 由(11)式可以看出, 本文所选择的邻域欧氏距离信息和邻域灰度差值信息的权重各占一半。

因此, 新的聚类中心v_i^N和隶属度值u_ik^N可以表示为:

(12)

(13)

利用该算法处理遥感图像时, 当遥感图像中的噪声点是邻域中心像素时, 由于它受周围邻域大量正常像素点的影响, 通过自动调整其隶属度, 可以较好地避免误分类; 此外该算法还可以最大限度地降低邻域噪声点对正常像素分类造成的影响, 提高算法抑制噪声的能力。

本文提出的邻域信息FCM算法的目标函数可以表示为:

(14)

算法步骤为:

1) 确定聚类数c, 加权指数m, 终止误差ε, 循环次数N_loop。选择差异图像中灰度值最大和最小的2个点作为聚类的2个初始聚类中心。

2) 初始化隶属度函数U=[u_ik^N]_c×n。

3) 开始循环, 计算目标函数J^N(U, V)。如果它相对于上次目标函数值的该变量小于ε或者达到最大循环次数N_loop, 则算法停止。

4) 计算新的聚类中心v_i^N, 并返回第2)步。

实验中, 选择的聚类个数为2, 最大循环次数为50, 中止误差为0.01, m的取值为2.5。

3 实验结果与分析

对于遥感图像变化检测的结果, 采用检测概率(detection probability, DP)、虚警率(false probability, FP)、漏警率(missing probability, MP)3个标准进行评价^[18]。检测概率为检测正确的像素点在整个图像中的百分比; 虚警率为实际非变化像素点, 却被检测为变化的像素点在整个图像中的百分比; 漏警率为实际变化像素点, 却没有被检测出来的像素点在整个图像中的百分比。将本文方法与文献[10]方法、文献[11]方法进行比较分析。

3.1 城市面积变化检测

图 2所示为Landsat7第四波段的关于某地区城市建设的用地面积发生变化遥感图像, 分别拍摄于2000年和2002年。2幅实验数据大小为512×512。

由图 2和表 1的实验结果可知, 采用文献[10]所提出的方法, 受到噪声影响使得虚警率较高, 而且没有考虑邻域信息的影响，一些变化的信息没有被检测处理, 因此最终的准确率为94.38%;采用文献[11]所提出的方法, 检测概率达到了98.53%;本文所提出方法, 虽然最终的检测概率略低于文献[11]的方法, 但是由于本算法考虑到邻域信息的影响, 其漏检率始终优于文献[11]的方法。

3.2 森林面积变化检测

图 3中实验采用的数据为www.land cover.org网站上的森林面积变化遥感图像, 图 3a)和图 3b)分别拍摄于2000年和2005年, 文中选用的实验数据大小为512×512。

结合图 3和表3的实验结果可知, 采用文献[10]所提出的方法, 由于噪声的影响出现了很多杂点, 原本变化的信息由于采用硬划分方法无法很好地检测出变化区域, 其漏警率已达到了20.97%, 使得最终的检测概率只有73.35%;采用文献[11]所提出的方法, 由于没有考虑到邻域信息对变化区域的影响, 本属于变化的区域仍然没有被检测出来, 其检测概率也只有82.48%;采用本文所提出的基于小波变换和邻域信息FCM方法, 很好地消除了噪声的影响, 使虚警率和漏警率较其他2种方法都降低了很多, 使检测概率提高到了87.5%。

对比2组数据的实验结果，对于变化信息较为明显的图像(如城市面积变化)，3种方法的最终检测的结果相差不是很明显，其检测概率都在92%以上，但对于变化区域不明显的图像，尤其是对于变化区域不连续的图像(如森林面积变化)，本文所提出的方法的检测结果明显优于另外2种方法。

4 结论

本文首先采用离散小波变换方法对差值图像和比值图像进行融合得到最终的差异图像，又提出了基于空间邻域信息的FCM变化检测方法，把空间距离信息和邻域灰度差值信息引入到目标函数中。为了验证算法的可行性和有效性，实验中分别采用城市用地面积和森林面积的变化情况，并与其他方法进行对比分析。对比实验结果可以看出，本文提出的算法对噪声具有较强的抑制能力，其检测结果更理想更有效。离散小波变换的融合规则不同，其所构造差异图的效果就不同，这将直接影响到后期的分割处理，如何得到更优的图像融合准则是本文进一步研究的方向^[19]。此外，由于分割算法采用欧氏距离，而遥感数据像元的分布不服从各向同性或球体分布，因此实际使用时也可采用更适合遥感数据特点的距离度量方式，可以更好地发挥FCM在遥感数据分类方面的优点，这也是本文要进一步努力的另一个方向。