1 引 言

变化检测就是通过不同时间观测确定地物或自然现象状态变化的过程^[1]。遥感影像以覆盖面积大、获取周期短的优势被广泛应用于变化检测中。变化检测在土地利用/覆盖调查、灾害监测、城市扩张监测等领域中有着重要的应用^{[2, 3, 4]}。

近年来随着传感器的不断发展，遥感影像的空间分辨率越来越高。面向对象法和多特征融合法是高分辨率遥感影像变化检测中的两种重要技术。对象可以视为一系列空间相邻、光谱相似的像元的集合，又可以称之为像斑。像斑中包含着丰富的光谱信息、纹理信息及形状信息。面向对象法可以有效避免高分辨率遥感影像中椒盐噪声的影响。文献^[5]提出了基于像斑的遥感影像分类后比较法。文献^[6]提出了一种基于像斑的变化向量分析方法。文献^[7]提出了一种基于像斑KS检验的高分辨率遥感影像变化检测方法。高分辨率遥感影像中地物细节丰富、光谱混淆现象严重^[8]，因此利用单一特征进行变化检测难以达到较高的精度。融合多特征构建稳定和适用的变化检测模型，引起了许多学者的研究。文献^[9]中利用多尺度纹理特征同光谱特征结合，提高了变化检测的精度。文献^[10]中利用地统计学提取纹理特征，并将纹理特征和光谱特征一起用于变化检测。文献^[11]中融合了光谱、纹理、空间等特征进行变化检测。上述多特征融合方法虽然融合了光谱特征和纹理特征，然而光谱特征和纹理特征均被赋予相同的权重，未能实现光谱特征和纹理特征的自适应组合。针对高分辨率遥感影像变化检测问题，本文提出一种面向对象的多特征融合的变化检测方法。

2 本文方法

本文方法首先通过eCognition软件进行多尺度分割以获取像斑，然后提取各像斑的光谱直方图和LBP纹理直方图，利用G统计量来计算不同时期像斑之间的光谱距离和纹理距离；将光谱距离和纹理距离进行自适应加权获取最终的像斑异质性；将所有像斑的异质性视为两个子高斯组成的混合高斯分布，采用EM算法获取子高斯对应的参数；最后在贝叶斯最小错误率的基础上对所有像斑进行变化和未变化划分。具体的流程图如图 1所示。

2.1 影像分割

影像分割将影像划分为若干个互不重叠的连通区域。区域内部具有较好的匀质性，相邻区域之间具有明显的异质性，所有区域的总和构成了整幅影像^[12]。
本文利用eCognition软件来获取像斑。该软件中多尺度分割方法包含3个参数：紧致度(compactness)、形状因子(shape)、尺度(scale)。紧致度越小，则分割结果的形状越零碎；形状因子越大，则分割结果的形状差异越小；尺度越小，则像斑数目越多。经过多次试验，可以得出以下结论：高分辨率影像分割时，形状因子应设置在(0,0.5］区间内。过大的形状因子会导致光谱权重过小，分割结果不能较好地反映地物的实际形状。紧致度应设置在[0.5,1)区间内，过小的紧致度会导致影像分割的结果过于零碎。尺度的选择则依据变化检测要求的最小检测面积来确定。进行变化检测时，为了保证像斑内部具有较好的匀质性，尺度设置可以适当偏小，使得影像有一定的过分割。

利用时序遥感影像进行变化检测，可以将多个时期的遥感影像组合为一个新影像，对新影像进行多尺度影像分割，然后将影像分割结果套合在各时期遥感影像上以获取像斑。

2.2 LBP纹理

局部二值模式(local binary patterns，LBP)是文献^{[13, 14]}于1996年提出的一种纹理描述算子。LBP采用统计分析方法来描述纹理，计算简单，因此广泛应用于图像分析中，如图像分割^[15]、居民地提取等。LBP计算公式如下

式中，P表示半径为R的局部窗口内像元个数；g_c、g_p分别表示中心像元、邻域像元灰度值，图 2是LBP计算的示意图。原始LBP算法权值模板固定，因此对纹理的方向具有敏感性。文献^[16]对原始LBP算法进行了改进并提出具有旋转不变特性的LBP，即不断旋转权值模板，然后将所有方向的LBP值的平均值作为旋转不变LBP值即L_P,R^ri 式中，R_L(L_P,R,i)表示权值模板顺时针旋转i次后的LBP值。

2.3 光谱和纹理特征提取

本文采用像斑的光谱直方图和LBP纹理直方图来描述像斑的光谱特征和纹理特征。直方图是一种统计特征，可以较好地反映像斑中光谱和纹理的分布情况。像斑的量化等级越大，对地物细节的分辨能力越强，然而会导致计算量增大。适当减少量化等级，可以有效提高计算效率，然而会丢失重要的细节信息，因此在变化检测应用中，需要在精度和效率之间平衡。

2.4 光谱和纹理特征距离

常用的直方图距离度量方法有多种，如KL距离、直方图交叉法、G统计量等。文中选取G统计量来度量直方图之间的距离。G统计量是一种非参数统计方法，计算简单，且不需要对直方图的分布预先作任何假设，因此广泛地应用于影像分析中^{[14, 15, 16]}。令参与变化检测的两个时期分别为t₁和t₂，多光谱影像的波段数为B，直方图的量化等级为L。波段b上的概率密度函数为f^b(1≤b≤B)，则两个时期像斑直方图在波段b上的距离为

式中，f_i^b表示i值出现的频率。对于所有像斑，其直方图的累积频率为1，则 将式(4)代入式(3)中，可得

2.5 光谱和纹理距离自适应组合

不同时期像斑在波段b上的异质性H^b，可以用光谱距离G_s^b和纹理距离G_t^b的加权平均来获取，即

式中，w_s和w_t分别表示光谱距离和纹理距离的权重，且w_s+w_t=1。光谱和纹理权重的确定是光谱和纹理距离组合的关键问题。光谱和纹理权重的确定依据以下原则：在纹理主导的区域，纹理权重较大，反之则光谱权重较大。文献^[16]采用一种自适应权重确定方法，但需要设定经验值，难以推广。文献^[17]使用固定的光谱和纹理权重，不能依据像斑属性选择合适的光谱和纹理权重。文献^[18]提出了一种自适应的权重确定方法，然而权重的确定过于烦琐。

为此，本文提出一种自适应的权重确定方法。该算法思想如下：标准差可以度量像斑内灰度分布的离散程度。对于纹理丰富的像斑，对应的标准差较大，而光谱丰富的像斑，其对应的标准差较小。当像斑的标准差小于给定阈值时，则认为像斑是光谱主导区域，反之则认为是纹理主导区域。两个时期像斑中任意一个像斑的纹理信息丰富时，则设置较大的纹理权重，否则设置较大的光谱权重。阈值的设定依赖于整幅影像，对于纹理特征丰富的影像，设定较小的阈值，反之则设定较大的阈值。令像斑在t₁和t₂两个时期标准差分别为std₁和std₂，设定的阈值分别为k₁和k₂，则对应的归一化权重获取方法如下

式中，aver(t_i)和std(t_i)分别表示t_i时期所有像斑的灰度均值和标准差；n_i介于区间[-2, 2]。

对于多光谱影像，为了充分利用所有波段包含的信息，采用各波段异质性的加权平均值作为像斑的异质性H，即

式中，w(b)表示波段b的归一化权重值；E(b)表示第b个波段的熵值。熵可以用来衡量单波段影像中所包含的信息量，信息量越大的波段对应的权重越大。

2.6 基于EM算法和贝叶斯最小错误率的变化检测

令所有像斑异质性值的集合为H={h₁,h₂,…,h_N}，其中N为像斑的个数，h_i表示第i个像斑的异质性值，H中的元素可以划分为变化和未变化两大类。在变化检测中，发生变化的像斑对应的异质性值较大，未发生变化的像斑对应的异质性值较小。因此文中假设像斑的异质性值近似满足由2个子高斯组成的混合高斯分布^[19]，即

式中，w_c、w_n表示变化类和未变化类；p(w_c)、p(w_n)为变化类和未变化类像斑所占比例，且p(w_c)+p(w_n)=1。假设变化类像斑和未变化类像斑的条件概率密度函数均服从高斯分布，则 式中，l∈{w_c，w_n}。式(9)中待求解的参数是l类对应的比例p(l)、均值μ_l和标准差δ_l。由于每个像斑的变化类别是未知的，因此这是一个带有缺失数据的参数求解问题。EM算法可以较好地解决带有缺失数据的参数求解问题^{[20, 21]}。它能在不完全数据中找到极大似然值。EM算法采用迭代的方式进行，迭代公式如下 式中，t为迭代的次数。EM算法需要初始值，可以采用k均值算法将所有像斑划分为两类，从而获取所需的初始值。利用EM算法进行参数求解的过程如下：①用k均值聚类方法将所有像斑的异质性值分为两类，从而获取p(l)、μ_l、δ_l的初始值；②利用式(11)-式(13)更新p(l)、μ_l、δ_l的值；③重复步骤②，直到连续两次迭代出的参数差值小于一定阈值。获取混合高斯分布对应的参数后，计算像斑属于变化类和未变化类的后验概率，即 基于贝叶斯最小错误率原则，且忽略和类别无关的常量p(h_k)，就可以获取像斑的变化类别：若p(h_k|w_c)p(w_c)＞p(h_k|w_n)p(w_n)，则h_k∈w_c，否则h_k∈w_n。

3 试验结果与分析

文中采用虚检率、漏检率、正确率3个指标来衡量变化检测的精度。虚检率是实际未变化、检测为变化的像元占检测变化像元的比例；漏检率是实际变化、检测为未变化的像元占检测未变化像元的比例；正确率是检测结果和实际结果一致的像元占所有像元的比例。虚检率、漏检率越低，正确率越高，则变化检测的效果越好。

不同时期的遥感影像由于光照条件、太阳高度角、季节因素等不同使得同种地物在影像上具有不同的灰度值。这种差异导致地物即便未发生变化，计算出来的异质性值也会偏大，在变化检测中产生虚检，因此需要对遥感影像进行辐射校正。基于辐射传输方程的辐射校正方法是最严密的辐射校正方法，然而该方法较复杂，且通常无法获取两张影像拍摄时期对应的大气有关参数，因此文中采用相对辐射校正方法^[22]。直方图匹配法是一种非线性的相对辐射校正方法，可以有效地减小上述因素对变化检测结果造成的影响。

3.1 试验1

本次试验数据为武汉地区2002年3月和2005年7月的QuickBird多光谱影像。影像大小为1249像素×1232像素，包含蓝、绿、红、近红外4个波段，分辨率为2.44 m。将2002年影像作为参考影像，对2005年影像利用直方图匹配法进行相对辐射校正。文中将两个时期影像组合为一个波段数为8的新影像。将新影像作为输入数据，设置shape=0.45，compactness=0.65，scale=150，利用eCognition软件进行多尺度分割获取影像分割结果。影像分割后获取像斑数目为268。图 3(a)、(b)、(c)、(d)分别为2002年、2005年利用红、绿、蓝3个波段合成的真彩色影像、相对辐射校正后的2005年影像及影像分割结果。从图 3(c)中可以看出，经过相对辐射校正，影像局部区域的亮度变大。图 4所示的是分别利用光谱特征和纹理特征进行变化检测时，量化等级同变化检测精度之间的关系。从图 4中可以看出，对于光谱特征，当量化等级为256时，正确率达到最大值，同时虚检率和漏检率达到最小值；对于纹理特征，当量化等级为128时，正确率达到最大值，同时虚检率和漏检率均达到最小值。故对于光谱和纹理特征，分别选取量化等级256、128进行变化检测。

为了验证本文方法的有效性，将本文方法同固定权值法进行了比较。图 5中黑色区域表示变化区域，白色区域表示未变化区域。图 5(a)是依据目视解译和实地调绘相结合而制作的标准变化检测结果，图 5(b)-图 5(e)是光谱特征权重依次为0、0.4、0.6、1.0的变化检测结果。图 5(f)是本文自适应权重法的变化检测结果(n₁=n₂=-1)。

图 5中5种变化检测结果的精度评定如表 1所示。

3.2 试验2

试验2的数据为天津地区2008年4月和2009年4月的SPOT-5遥感影像，影像大小为4462像素×3315像素，包含绿、红、近红外3个波段，和全色波段融合后分辨率为2.5 m。以2008年影像为参考影像，对2009年影像利用直方图匹配法进行相对辐射校正。文中将两个时期的影像组合为一个波段数为6的新影像。将新影像作为输入数据，设置shape=0.4，scale=150，compactness=0.7，利用eCognition软件进行多尺度分割获取影像分割结果。影像分割后获取像斑个数为2641。图 6(a)、(b)、(c)、(d)分别为2008年、2009年利用近红外、红、绿3个波段合成的伪彩色影像、相对辐射校正后的2009年影像及影像分割结果。

图 7表示的是量化等级与变化检测精度之间的关系。从图 7中可以看出，对于光谱特征，当量化等级为16时，正确率达到最大值，虚检率达到最小值，同时漏检率也较小；对于纹理特征，当量化等级为64时，正确率达到最大值，虚检率达到最小值，此时漏检率也较小。因此对于光谱和纹理特征，依次选择量化等级16、64进行变化检测。

图 8所示的是各种方法的变化检测结果，图中黑色代表变化区域，白色代表未变化区域。图 8(a)是标准的变化检测结果。图 8(b)-图 8(e)依次是光谱权重为0、0.2、0.8、1.0的变化检测结果。图 8(f)是本文方法的变化检测结果(n₁=n₂=1.5)。分析图 8可以得出以下结论：仅利用光谱特征进行变化检测，存在着较多的虚检和一定的漏检(图 8(e)中矩形区域代表虚检，椭圆形区域代表漏检)；仅利用纹理特征进行变化检测，则存在着严重的漏检(图 8(b)中椭圆形区域)；利用光谱特征和纹理特征加权进行变化检测，当光谱特征权重为0.2时，虽然部分图 8(e)中漏检的像斑被正确判别为变化(图 8(c)中椭圆形区域)，但仍然存在着一定的虚检(图 8(c)中矩形区域)，当光谱特征权重为0.8时，存在着较多的虚检和漏检(图 8(d)中矩形区域代表虚检，椭圆形区域代表漏检)；利用本文方法进行变化检测时，许多变化的像斑均能被正确识别为变化(图 8(f)中椭圆形区域)，未变化像斑能被正确判别(图 8(f)中矩形区域)。

图 8中5种变化检测结果的精度评定如表 2所示。

将本文方法变化检测的精度同文献^[16]和文献^[17]中的方法进行对比(表 3)。从表 3中可以看出，本文方法变化检测的精度最高。

综合试验1和试验2可以得出以下结论：

(1) 针对两组试验数据，采用本文方法变化检测的效果最优。QuickBird影像变化检测正确率达到了0.823(表 3)，SPOT-5影像变化检测正确率达到了0.767(表 3)。单独利用光谱特征进行变化检测，虚检率较高，单独利用纹理特征则漏检率较高，这主要是因为不同时期的遥感影像由于光照条件、太阳高度角、季节因素等不同使得同种地物在影像上具有不同的灰度值，因此仅利用光谱特征进行变化检测容易产生虚检，而纹理特征反映的是地物的结构信息，对上述因素不敏感，因此不容易产生虚检，但纹理特征无法作为区分不同地物的主要依据，因此会导致漏检。采用固定权值法将光谱特征和纹理特征进行组合，正确率得到一定提高，然而不能依据像斑的光谱和纹理特征选择不同的权重，而本文方法可以自动调节两者的权重，充分利用光谱特征和纹理特征，因此变化检测的精度优于前者。

(2) 分析图 5(f)和图 8(f)可以看出，本文方法将许多未变化的道路错误判别为变化。这是因为两个时期影像拍摄角度和太阳高度角的不同，导致道路的光谱差异较大，相对辐射校正法未能有效减小其差异。基准期影像上道路呈白色，检测期影像上道路呈灰色，且两个时期道路的匀质性都较好，构建异质性时光谱特征占主导，导致异质性较大；另一方面道路包含像元数较少，构建的异质性不稳定，从而容易导致虚检。

(3) 试验1和试验2的数据对应着不同的最优量化等级，这表明量化等级依赖于具体影像，但是对比图 4和图 7可以发现，当量化等级为256时，像斑的灰度信息得以全部保留，变化检测的各项精度或者达到最优，或者与最优值相差小于0.05。因此进行变化检测时，可以设置量化等级为256。

表 4给出了本文方法针对两种试验数据的运行时间。编程平台为Microsoft Visual C++6.0，操作系统为Windows XP，内存为3G，CPU为酷睿2、双核E7000。

从表 4中可以看出，光谱距离和纹理距离的计算消耗时间相对较长。这是因为特征提取针对整幅影像，距离组合、EM算法、变化判别均基于像斑，而像斑数目远小于原始影像的尺寸。EM算法是一种迭代算法，经过本文改进后，EM算法耗时较短。对于一幅大小为4462像素×3315像素的影像数据，整个变化检测流程需要消耗时间约为83 s，计算效率较高。

4 结 论

本文方法实现了光谱特征和纹理特征的自适应加权组合，并以像斑为单位进行变化检测，为融合多特征进行变化检测提供了一条有效途径。采用EM算法和贝叶斯最小错误率理论获取变化结果，实现了变化检测阈值的自动获取。试验结果表明，该方法能够充分利用像斑的光谱特征和纹理特征，有效地提高变化检测的正确率，降低虚检率和漏检率。该方法兼顾了光谱特征和纹理特征，非常适合对同物异谱和同谱异物现象严重的高分辨率遥感影像进行变化检测。对于绝大多数地物类别变化可以较好地检测，但是对于线状地物(如道路)容易造成虚检。为了进一步提高变化检测的精度，在本文工作的基础上，还需要引入颜色、形状、边缘等特征。