图像变化检测方法

‘壹’ 基于典型相关分析的多元变化检测

( 一) 基于典型相关分析的多元变化检测方法

为了解决传统方法应用于多源遥感影像变化检测中的不足，丹麦学者 Nielsen 等( 1998) 率先提出了多元变化检测 ( Multivariate Alteration Detection，MAD) 的概念和方法。这种方法的基本思想是基于典型相关分析方法，对两幅不同时相的分别具有 p 和 q 个通道的影像 ，不妨假设 p≤q，且 X 和 Y 都具有零均值，即 E{ X} = E{ Y} = 0，分别构造两个线性组合:

退化废弃地遥感信息提取研究

然后将它们之间的差值a^TX－b^TY作为衡量变化信息的测度。

理论上，任意一对向量a和b对应于一种变化信息测度，因此需要根据一定的判断准则来找出最优的a和b，使之最有利于变化分析。MAD方法将差值变量a^TX－b^TY的方差最大化作为其判断准则。考虑到a和b同时乘上常数c，a^TX－b^TY的方差会变成原来的c²倍，因此要对a和b进行限制，最常见的就是限定变量a^TX和b^TY均具有单位方差，即Var{a^TX}=Var{b^TY}=1，这样，MAD问题可以表述为以下数学形式:

退化废弃地遥感信息提取研究

在单位方差约束条件下，根据变量和的方差计算公式，可得:

退化废弃地遥感信息提取研究

MAD方法要求上式中的相关系数Corr{a^TX，b^TY}≥0，也就是线性组合变量a^TX和b^TY之间成正相关，因为只有这样的a^TX－b^TY才能表示差异信息。从上式可以看出，欲使a^TX－b^TY的方差最大，只需找出使a^TX和b^TY之间的相关性最小的满足单位方差约束条件的向量a和b即可。这说明MAD方法的核心思想在于寻找相关性最小的针对原始变量的线性组合对，它们之间的差值就包含了最大的差异信息，通过分析这些差值就可以找出时相之间的变化。

要解决这一问题，可以从典型相关分析原理描述中寻求思路。典型变换中得到的每一对典型变量都是两组原始变量的一对线性组合，它们之间的相关系数依次递减。而在所有可能的针对原始变量的、满足单位方差约束条件的线性组合对中，第一对典型变量之间的相关系数是最大的，而最后一对即第p对典型变量之间的相关系数是最小的。因此，MAD就转化为了典型变换求解问题，只需把典型变换得到的一系列典型变量对按典型相关系数从小到大重新排序即可。这样，就把MAD定义为如下变换:

退化废弃地遥感信息提取研究

式中:a_i和b_i就是典型变换所得到的各典型变量U_i和V_i对应的线性组合系数向量，这一变换被称为MAD变换，所得到的差值变量称之为MAD变量，其中第i个MAD变量为第p+1－i对典型变量之间的差值，即:

退化废弃地遥感信息提取研究

(二)MAD变化检测方法的基本流程

1.数据预处理

数据预处理包括影像数据的几何纠正，两时相影像之间的几何配准和辐射归一化，以及影像剪裁。从两时相影像中剪裁出两幅大小合适、位置范围完全相同的研究区域影像块，作为后面进行计算的基础数据。

2.计算协方差和互协方差矩阵

按以下公式来计算两时相影像X和Y之间的协方差和互协方差矩阵:

退化废弃地遥感信息提取研究

式中:n为样本数即像元个数;X_i和Y_i分别为X和Y上第i个像元的光谱向量值; 分别为X和Y的均值向量。

3.求解广义特征方程

构造典型相关的广义特征方程，求解从大到小排序的特征值ρ²_j和相应的满足单位方差约束条件的特征向量a_j和b_j(j=1，2，…，P)。

4.计算典型变量

用上一步得到的特征向量a_j和b_j作为系数向量，分别对X和Y做线性组合，得到典型变量U_j和V_j(j=1，2，…，p)。

5.计算MAD变量

根据式(6－6)计算典型变量之间的差值，得到MAD变量，即构成MAD变化检测的结果影像。

(三)结果与分析

对经过预处理的两幅影像做典型变换，求解出典型相关特征方程的特征值和特征向量，如表6－4所示。值得注意的是，两个方程的特征向量维数不同，对应于两个时相原始影像的波段数不相同。特征向量矩阵各列从左到右依次为对应特征值从大到小的特征向量。

表 6 －4 典型相关广义特征方程的特征值和特征向量

典型变换得到的4对典型变量影像如图6－6所示，图中左右两列分别为2001年和2005年的典型变量影像，从上到下依次为对应特征值从大到小的典型变量对，记第i对典型变量为2001CV_i和2005CV_i，它们之间的相关系数即典型相关系数，以及典型相关系数的平方即方差重叠度如表6－5所示。

表6－5 典型相关系数和典型变量方差重叠度

从表6－5可以看出，两时相遥感影像的四对典型变量中，第一对典型变量2001CV₁和2005CV₁之间的相关系数达到0.4939，方差重叠度为0.2439，反映了两时相之间大部分为不变信息。第二、三对典型变量之间的相关系数相对较低，分别只有0.2429和0.2720，既包含了部分的不变信息，也包含了部分变化信息。第四对典型变量间的相关系数最低，方差重叠度也仅有0.0019，即前者的全部变差中只有0.19%能够用后者的变差来解释，这意味着二者之间存在着很大的差异，可能包含了大量的时相间差异或变化信息。对图6－6进行目视判读，可以看出这四对典型变量影像的相似性的确是依次降低的，即所包含的不变信息依次减少，这与理论上的结论是一致的。但是两时相的变化信息并没有只集中表现在第四对变量影像中，可以看出，变化信息是分布于各对典型变量中的。

表6－6列出了各对典型变量与原始变量之间的相关系数，反映了它们之间的相关关系。

图6－6 四对典型变量影像

表6－6 典型变量与原始变量之间的相关系数

将4对典型变量依次相减后逆序排列，就得到了4幅MAD变量影像，如图6－7所示，图中从上到下、从左到右依次为MAD₁、MAD₂、MAD₃和MAD₄变量影像。这些MAD变量影像反映了两时相原始影像之间的差异，MAD影像中像元值越接近于0，说明对应位置上的不同时相之间的差异越小，越有可能没有发生变化;反之，像元值越远离0，即绝对值越大，表现在影像中很亮或很暗的区域，对应位置的不同时相之间的差异越大，标志着越有可能发生了变化。

图6－7MAD变量影像

计算得到各MAD变量的方差、方差贡献率和累计方差贡献率，如表6－7所示。

表6－7 各MAD变量的方差、方差贡献率、累计方差贡献率

从表6－7中可以看出，MAD₃和MAD₄的方差贡献率分别为39.9135%和31.6852%，从理论上讲，这两个变量集中了两个时相的大部分变化信息，但是图6－7对这些变化信息表现得并不是很清晰，MAD₂的方差贡献率虽然仅有8.7439%，但是也清晰地反映出了部分变化信息。这表明对于MAD变换得到的各MAD变量影像，并非方差越大的MAD变量包含的变化信息就越多。

表6－8列出了4个MAD变量与10个原始影像波段变量之间的相关系数。

表6－8 MAD变量与原始变量之间的相关系数

‘贰’ 边缘检测

图像边缘信息主要集中在高频段，通常说图像锐化或检测边缘，实质就是高频滤波。我们知道微分运算是求信号的变化率，具有加强高频分量的作用。
在空域运算中来说，对图像的锐化就是计算微分。
由于数字图像的离散信号，微分运算就变成计算差分或梯度。
图像处理中有多种边缘检测（梯度）算子，常用的包括普通一阶差分，Robert算子（交叉差分），Sobel算子等等，是基于寻找梯度强度。拉普拉斯算子（二阶差分）是基于过零点检测。通过计算梯度，设置阀值，得到边缘图像。

边缘检测 是图像处理和计算机视觉中的基本问题。
边缘检测的目的 是标识数字图像中亮度变化明显的点。
图像属性中的显着变化通常反映了属性的重要事件和变化。

边缘检测是图像处理和计算机视觉中，尤其是特征提取中的一个研究领域。

除非场景中的物体非常简单并且照明条件得到了很好的控制，否则确定一个用来判断两个相邻点之间有多大的亮度变化才算是有边界的阈值，并不是一件容易的事。实际上，这也是为什么边缘检测不是一个微不足道问题的原因之一。

图像边缘检测大幅度地减少了数据量，并且剔除了可以认为不相关的信息，保留了图像重要的结构属性。
有许多方法用于边缘检测，它们的绝大部分可以划分为两类：

滤波 做为边缘检测的预处理通常是必要的，通常采用高斯滤波。
已发表的边缘检测方法应用计算边界强度的度量，这与平滑滤波有本质的不同。 正如许多边缘检测方法依赖于图像梯度的计算，他们用不同种类的滤波器来估计x-方向和y-方向的梯度。

其它一些边缘检测操作是 基于亮度的二阶导数 。这实质上是亮度梯度的变化率。

在理想的连续变化情况下，在二阶导数中检测过零点将得到梯度中的局部最大值。另一方面，二阶导数中的峰值检测是边线检测，只要图像操作使用一个合适的尺度表示。

如上所述，边线是双重边缘，这样我们就可以在边线的一边看到一个亮度梯度，而在另一边看到相反的梯度。这样如果图像中有边线出现的话我们就能在亮度梯度上看到非常大的变化。

①滤波： 边缘检测算法主要是基于图像强度的一阶和二阶导数，但导数的计算对噪声很敏感，因此必须使用滤波器来改善与噪声有关的边缘检测器的性能。需要指出，大多数滤波器在降低噪声的同时也导致了边缘强度的损失，因此，增强边缘和降低噪声之间需要折中。
②增强: 增强边缘的基础是确定图像各点邻域强度的变化值。增强算法可以将邻域(或局部)强度值有显着变化的点突显出来。边缘增强一般是通过计算梯度幅值来完成的。
③检测： 在图像中有许多点的梯度幅值比较大，而这些点在特定的应用领域中并不都是边缘，所以应该用某种方法来确定哪些点是边缘点。最简单的边缘检测判据是梯度幅值阈值判据。
④定位： 如果某一应用场合要求确定边缘位置，则边缘的位置可在子像素分辨率上来估计，边缘的方位也可以被估计出来。
在边缘检测算法中，前三个步骤用得十分普遍。这是因为大多数场合下，仅仅需要边缘检测器指出边缘出现在图像某一像素点的附近，而没有必要指出边缘的精确位置或方向。

边缘检测的实质是采用某种算法来提取出图像中对象与背景间的交界线。我们将 边缘定义 为图像中灰度发生急剧变化的区域边界。
边缘检测方法
图像灰度的变化情况可以用图像灰度分布的梯度来反映，因此我们可以用 局部图像微分技术 来获得边缘检测算子。 经典的边缘检测方法，是通过对原始图像中像素的某小邻域构造边缘检测算子来达到检测边缘这一目的的。

‘叁’ 图像观察法

图像观察法是珠宝玉石材料无损检测中使用最广泛的方法。常用的仪器是各类放大镜,宝石显微镜(倍率为45~75,带暗场、明场等照明功能),偏光显微镜以及其他低倍率、功能单一的显微镜,超景深显微镜等。

(一)宝石显微镜

宝石显微镜是用来观测宝石表面及内部特征最常用的仪器,对彩色宝石产地特征的观察和信息采集具有最重要意义。

宝石显微镜主要由以下几个部件组成:

(1)光学系统(透镜系统),包括目镜、物镜等。

(2)照明系统,包括底光源、顶光源、光量强度调节按钮等。

(3)机械系统,包括支架、宝石夹、焦距调节旋钮等。

图2-1为直立式宝石显微镜的原理及结构图。

图2-8 DiamondView^TM下可见红宝石内的生长条带和裂隙(10×)

DiamondView^TM在使用时可将具光滑表面的样品置于紫外光下,拍摄并记录其紫外荧光图样。不同品种的宝石具有不同的发光特征,不同产地的同一品种宝石其荧光特性(如荧光强度)也存在差异。对于钻石而言,可将其荧光图像直接与DiamondView^TM软件中存有的各种天然和合成钻石的紫外荧光图样进行对比,从而得出结论。

‘肆’ 遥感变化信息检测的主要方法

本节所研究的主要是基于像元级的遥感变化信息检测方法。当今国内外常用的遥感变化信息检测方法主要有分类后比较法、代数运算法、光谱特征变异法、主成分分析法等。

( 一) 光谱特征变异法

光谱特征变异法是使用最为广泛的一种遥感变化信息检测方法，其原理是将两时相遥感影像的相关波段进行融合、组合，如果地物类型发生了变化，该区域的光谱就会发生变异或突变，与周围地物失去协调性，使得能通过肉眼识别出来。该方法容易受到研究区域地物光谱特征的影响，可能会丢失细小的变化图斑，但是在一般情况下总体效果良好 ( 冯德俊，2004) 。

图 4 －8 HIS 融合法结果

把研究区震前 IKONOS 的全色波段影像和震后QuickBird 的多光谱影像运用光谱特征变异法中的 HIS进行了融合，结果见图 4 －8。从图中可以看出，没有发生变化的区域光谱特征和分辨率都得到了加强 ( 空间分辨率变为 1m) ，其中水体、河流为蓝色，山地为褐色，植被信息为绿色，白色区域明显与四周的地物和背景信息光谱不协调，这就是发生变化的区域。

( 二) 主成分分析法

主成分分析法在遥感变化信息检测中使用很多，一般能够取得很好的检测效果，它能够在一定程度上消除影像内部各波段间的相关性，提高变化信息检测的效率和精度。

图 4 － 9 为 IKONOS 融合后影像与 QuickBird 融合后影像求差并取绝对值后进行主成分变换的第一个主成分图像。由主成分变换的特性知道，变换结果的第一分量集中了影像绝大部分的变化信息，而其他分量则主要反映了波段之间的差异性。第一分量就代表着变化信息。图 4 －9 中白色区域为发生变化的区域。

基于主成分分析的遥感变化信息检测方法仍然存在着一些缺陷: ① 由于主成分变换所得到的各个分量往往失去了原有的物理意义，所以还需要参考其他数据来分析地物类型变化与否及其因果关系。② 主成分分析 ( PCA) 是基于二阶统计的方法，只有在信号的统计分布满足高斯分布的条件下，才能完全消除信号间的相关性，而对于非高斯信号则只能去除信号间的二阶相关性 ( 钟家强、王润生，2006) 。在多时相遥感影像中，各种地类的光谱特性几乎都不能满足高斯分布，因此经过主成分变换后的各成分图像间仍然存在高阶相关性，而这些相关信息会直接影响到变化信息的检测和提取。由于这个原因，在做主成分分析时，常常导致把这些高阶相关信息转变为噪声，如图 4 －9 所示，白色的 “斑点”遍布整个分量影像，又和变化信息掺杂在一起 ( 变化的区域也为高亮的白色) ，这使得在提取真正变化信息的时候遇到困难，如何有效地消除多时相图像间的高阶相关信息，避免这些 “伪变化”的噪声，对于变化信息的检测和提取具有非常重要的意义。