CN105354845A - 一种遥感影像半监督变化检测方法 - Google Patents

Abstract

本发明公开了一种遥感影像半监督变化检测方法，包括：获取两时相的原始光学遥感影像；对所述原始光学遥感影像进行影像配准；利用直方图调整方法对影像配准后的遥感影像进行辐射归一化校正；根据辐射归一化校正后的遥感影像计算光谱角信息；将所述辐射归一化校正后的遥感影像和所述光谱角信息组合，作为FPTSVM的输入；采用FPTSVM方法进行学习，不断调整SVM的分类超平面，直到达到指定的学习迭代次数；利用最终的分类超平面确定影像的变化区域和非变化区域。本发明可以提高变化检测的速度和精度。

Description

一种遥感影像半监督变化检测方法

技术领域

本发明涉及遥感影像处理技术领域，尤其涉及一种遥感影像半监督变化检测方法。

背景技术

遥感影像的变化检测是根据同一区域不同时相的遥感影像，辨识观测对象或者现象的状态变化过程。它已广泛应用于资源管理与规划、环境保护等诸多领域，并为相关部门提供科学决策的依据。当前遥感影像变化检测方法主要有：代数法、变换法、分类比较法、高级模型法、GIS集成法、视觉分析法和其它方法。

其中代数法以其操作简单、易于实现成为当前变化检测中使用最广泛的方法之一。其差异影像通常是通过对两幅遥感影像的某一波段进行差值、比值运算等来构造的，由于该方法只利用了多光谱遥感影像的一个波段，未充分利用其全部有用信息，故不能较好地反映地面真实变化情况。

代数法的核心是阈值的选择，针对此，国内外学者提出了许多有效的方法，其中具有代表性的为L.Bruzzone等提出一种基于统计最小错误率的Bayes判别准则的变化检测方法，虽然该方法可使差异影像中变化与非变化像元的错分概率最小，能较好地确定变化检测的阈值，在一些特定实例中也取得了较好的效果，然而差异影像的像元值必须服从高斯分布，这一假设在有些情况下并不成立。

为了避免阈值选择对变化检测造成的不利影响，有学者提出采用分类的思想进行变化检测。该类方法是将差异影像视为变化和非变化两类像元的集合，然后利用分类器对其进行分类，最终实现遥感影像的变化检测。由于支持向量机SVM在分类时具有较好的收敛性和分类精度，因此已有学者将其用于变化检测中。在此基础上，F.Bovolo等提出一种改进的渐进转导支持向量机(ProgressiveTransductiveSupportVectorMachine，PTSVM)的遥感影像变化检测方法。其特点是每次迭代时，从离SVM分解面一定距离内的区域内选择标记样本，另外为了选择PTSVM模型的参数，基于贝叶斯准则从待检测像元中抽取部分容易判别的像元进行标记作为种子，即伪标签样本，然后利用这些伪标签样本并利用粒度搜索方法寻找最优的模型参数。

常规的PTSVM用于遥感影像的变化检测时的不足为：(1)仅使用原始影像的多光谱信息；(2)半监督学习时，每次迭代时，侯选未样本集较大，导致学习速度较慢。

发明内容

本发明所要解决的技术问题在于，提供一种遥感影像半监督变化检测方法，该方法基于快速渐进转导支持向量机(FastProgressiveTransductiveSupportVectorMachine，FPTSVM)，提高了变化检测的速度和精度。

为了解决上述技术问题，本发明提供了一种遥感影像半监督变化检测方法，包括：

获取两时相的原始光学遥感影像；

对所述原始光学遥感影像进行影像配准；

利用直方图调整方法对影像配准后的遥感影像进行辐射归一化校正；

根据辐射归一化校正后的遥感影像计算光谱角信息；

将所述辐射归一化校正后的遥感影像和所述光谱角信息组合，作为FPTSVM的输入；

采用FPTSVM方法进行学习，不断调整SVM的分类超平面，直到达到指定的学习迭代次数或前后两次迭代后的非标记样本没有变化；

利用最终的分类超平面确定影像的变化区域和非变化区域。

实施本发明，具有如下有益效果：

1)在基于FPTSVM的变化检测中，联合X_SA和原始多光谱影像的波段信息作为输入，利用这两个特征，优于仅仅使用原始多光谱影像的波段信息的SVM方法；

2)变化检测中，FPTSVM学习过程中，仅从SVM两个分界面内的未标记样本中选择侯选样本，选择的策略是基于区域的学习策略，即从离分解面一定距离内的区域选择样本，从而进一步加快学习速度。可以加快变化检测的速度。当达到一定的学习次数后，输出变化检测图像。该检测方法通过FPTSVM可以提高多时相遥感影像差异影像的分类精度及速度。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1是本发明提供的遥感影像半监督变化检测方法的一个实施例的流程示意图；

图2是原始光学遥感影像；

图3是本发明方法与其他方法的处理结果对比图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

图1是本发明提供的遥感影像半监督变化检测方法的一个实施例的流程示意图，如图1所示，包括：

S101、获取两时相的原始光学遥感影像X₁和X₂。

其中，X₁、X₂是同一区域不同时相的两幅高分辨率光学遥感影像。

S102、对所述原始光学遥感影像X₁和X₂进行影像配准。

具体的，步骤S102具体包括步骤：

S1021、采用ENVI14.8遥感软件对光学遥感影像X₁和X₂进行几何粗校正。

几何粗校正具体操作步骤为：(1)显示基准影像和待校正影像；(2)采集地面控制点GCPs；GCPs应均匀分布在整幅图像内，GCPs的数目至少大于等于9；(3)计算误差；(4)选择多项式模型；(5)采用双线性插值进行重采样输出。

其中的双线性差值法为：若求未知函数f1在点P＝(x,y)的值，假设我们已知函数f1在Q₁₁＝(x₁,y₁),Q₁₂＝(x₁,y₂),Q₂₁＝(x₂,y₁),及Q₂₂＝(x₂,y₂)四个点的值。如果选择一个坐标系统使得这四个点的坐标分别为(0,0)、(0,1)、(1,0)和(1,1)，那么双线性插值公式就可以表示为：

f1(x,y)≈f1(0,0)(1-x)(1-y)+f1(1,0)x(1-y)+f(0,1)(1-x)y+f1(1,1)xy。

S1022、利用自动匹配与三角剖分法对几何粗校正后的遥感影像X₁和X₂进行几何精校正。

其中，三角剖分法为，采用逐点插入法构建Delaunay三角网，对每一个三角形，利用其三个顶点的行列号与其对应的基准影像同名点的地理坐标来确定该三角形内部的仿射变换模型参数，对待校正影像进行纠正，得到校正后的遥感影像。

S103利用直方图调整方法对影像配准后的遥感影像进行辐射归一化校正。

具体的，步骤S103包括步骤：将影像配准后的遥感影像的像素灰度值分别调整为从而实现辐射归一化校正，

其中，影像配准后的遥感影像各个波段服从高斯分布，μ_f、σ_f分别为参考影像的像素均值和标准差，g_r、μ_r、σ_r分别为影像配准后的遥感影像的像素灰度值、均值和标准差。

本步骤中按照公式分别对影像配准后的遥感影像X₁和X₂进行辐射归一化校正。

S104、根据辐射归一化校正后的遥感影像计算光谱角信息。

具体的，步骤S104包括步骤：

根据辐射归一化校正后的遥感影像X₁和X₂计算光谱角信息X_SA，其中，

$X_{SA}(i,j)=arccos\left(\underset{b=1}{\overset{B}{\Σ}}\right(X_{1b}\left(i,j\right)X_{2b}\left(i,j\right))/\sqrt{\underset{b=1}{\overset{B}{\Σ}}{X}_{1b}^2(i,j)}\sqrt{\underset{b=1}{\overset{B}{\Σ}}{X}_{2b}^2(i,j)}),$

式中，B表示每一个时相的遥感影像的波段数目，(i,j)是影像的坐标，X_1b表示X₁的b波段影像，X_2b表示X₂的b波段影像。

S105、将所述辐射归一化校正后的遥感影像和所述光谱角信息组合，作为FPTSVM的输入。

具体的，将辐射归一化校正后的遥感影像X₁、X₂和光谱角信息组合X_SA组合为作为FPTSVM的输入。

S106、采用FPTSVM方法进行学习，不断调整SVM的分类超平面，直到达到指定的学习迭代次数或前后两次迭代后的非标记样本没有变化。

具体的，步骤S106包括步骤：

S1061、建立拉格朗日函数为

$\begin{array}{c}L(\mathrm{\α},{\mathrm{\α}}^{*},y^{*})=\underset{l=1}{\overset{n}{\Σ}}{\mathrm{\α}}_l+\underset{u=1}{\overset{m}{\Σ}}{\mathrm{\α}}_u^{*}-\frac{1}{2}\left(\underset{l=1}{\overset{n}{\Σ}}\underset{i=1}{\overset{n}{\Σ}}k\right(x_l,x_i)y_l{y}_i{\mathrm{\α}}_l{\mathrm{\α}}_i\\ +2\underset{l=1}{\overset{n}{\Σ}}\underset{u=1}{\overset{m}{\Σ}}k(x_l,x_u^{*})y_l{y}_u^{*}{\mathrm{\α}}_l{\mathrm{\α}}_u^{*}+\underset{u=1}{\overset{m}{\Σ}}\underset{j=1}{\overset{m}{\Σ}}k(x_u^{*},x_j^{*})y_u^{*}{y}_j^{*}{\mathrm{\α}}_u^{*}{\mathrm{\α}}_j^{*})\end{array}$

$Subjectto\left\{\begin{array}{c}\begin{array}{cc}0\≤{\mathrm{\α}}_l\≤C,& \begin{array}{cc}for& 1\≤l\≤n\end{array}\end{array}\\ \begin{array}{cc}0\≤{\mathrm{\α}}_u^{*}\≤C^{*},& \begin{array}{cc}for& 1\≤l\≤m\end{array}\end{array}\\ \underset{l=1}{\overset{n}{\Σ}}{y}_l{\mathrm{\α}}_l+\underset{u=1}{\overset{n}{\Σ}}{y}_u^{*}{\mathrm{\α}}_u^{*}=0\end{array}\right.$

式中，x_i为标记矢量，标记y_i∈{±1}，i＝1,2,…,n，n为标记矢量个数，-1表示没有变化类，+1表示变化类，为转导矢量，相应的伪标记为j＝1,2,…,m，m为转导矢量个数，k(·,·)表示核函数，C和C^*分别表示对有标记样本和转导样本的惩罚参数，α＝(α₁,α₂,…,α_n)和分别表示两个拉格朗日乘数；

S1062、建立分类决策函数为 $f\left(x\right)=\underset{l=1}{\overset{n}{\Σ}}{\mathrm{\α}}_l{y}_{l}k(x,x_l)+\underset{u=1}{\overset{m}{\Σ}}{\mathrm{\α}}_u^{*}{y}_{l}k(x,x_u^{*})+d;$ 式中d是一个常数，表示分类超平面的偏移量；

S1063、进行FPTSVM的归纳学习阶段：利用已知标记的样本，获取初始的超平面h₀和两个分类边界h₁、h₂，此时候选的转导样本集合式中，U₀是非标记样本集合；其中，当迭代归纳学习时，在第k+1次迭代时，标记样本集更新为：候选转导样本集合为 $S_{k+1}=\{\tilde{x}\∈S_k|h_1^{\left(k\right)}\≤f\left(\tilde{x}\right)\≤h_2^{\left(k\right)}and\tilde{x}\∉{\mathrm{\Ψ}}_k^{+}\∪{\mathrm{\Ψ}}_k^{-}\},$ 式中，和分别表示第k次迭代的转导样本集，分别表示第k次迭代的分类边界；

S1064、进行FPTSVM的转导学习阶段：将转导样本的惩罚参数设置为：式中，是初始的惩罚参数，C是初始的有标记样本的惩罚参数，k表示迭代的次数，0≤k≤G，G是最大的迭代次数；其中，在FPTSVM的转导学习阶段中，在第k次迭代前，如果任何一个转导样本的标记和前一次的标记不同，则该转导样本被重新作为未标记样本，并从当前的转导样本集中删除；在转导样本的选择过程中，如果两个转导样本的决策值相同，且相应的输入矢量相同，则删除其中一个样本；

S1065、不断迭代进行归纳学习和转导学习，直到达到指定的学习迭代次数或前后两次迭代后的非标记样本没有变化，得到最终的分类超平面。

S107、利用最终的分类超平面确定影像的变化区域和非变化区域。

本发明的效果可通过以下实验结果与分析进一步说明：

1、实验数据：

本发明的实验数据为LandsatTM数据，位于巴西的亚马逊森林地区的2幅遥感影像，获取时间分别为2000年7月和2006年7月，选用前4个波段，实验区大小为320像元×320像元，图2(a)和(b)分别为两个时相的真彩色遥感影像。变化参考图如图2(c)所示，共有16,826个变化像素。

2、实验方法：

方法1：基于SVM的归纳学习方法(ISVM)，即仅仅利用少量标记样本的分类方法。

方法2：意大利的BruzzoneL.等在文章“AnoveltransductiveSVMforsemisu-pervisedclassificationofremote-sensingimages”(IEEETransactionsonGeoscienceandRemoteSensing,2006,44(11):3363-3373.)中所提的方法。

方法3：本发明方法。

3、实验内容与结果分析：

用不同方法对如图2(a)和图2(b)所示的两幅不同时相的多光谱图像进行变化检测，其中图3(a)为现有方法1得到的变化检测结果图，图3(b)为现有方法2得到的变化检测结果图，图3(c)为本发明方法得到的变化检测结果图。从图3可以看出，现有方法1和现有方法2结果图中存在较多的伪变化信息，本发明方法伪变化信息很少，比较接近参考图2(c)。

表1各变化检测方法性能评价

对实验变化检测结果进行定量分析，检测性能用错检数FP、漏检数FN、总错误数OE和Kappa系数四个指标来衡量。FP、FN和OE越接近于0、Kappa系数越接近于1，表明变化检测方法的性能越好。结果如表1所示。

从表1中的四个评价指标可以看出，所有评价指标均优于其他方法。所提算法在总体检测性能上是最好的。

为了更进一步表明所提算法的快速性，以各算法运行10次的平均运行时间进行比对，如表2所示。从表2可以看出，本发明方法要比方法2快约4倍。从表1和表2可以看出，从检测效果和运行时间上来看，本发明方法是有效的。

表2各算法10次运行的平均时间(单位：秒)

	方法1	方法2	本发明方法
				运行10次的平均时间	0.01	11.50	2.90

因此，实施本发明，具有如下有益效果：

1)在基于FPTSVM的变化检测中，联合X_SA和原始多光谱影像的波段信息作为输入，利用这两个特征，优于仅仅使用原始多光谱影像的波段信息的SVM方法；

2)变化检测中，FPTSVM学习过程中，仅从SVM两个分界面内的未标记样本中选择侯选样本，选择的策略是基于区域的学习策略，即从离分解面一定距离内的区域选择样本，从而进一步加快学习速度。可以加快变化检测的速度。当达到一定的学习次数后，输出变化检测图像。该检测方法通过FPTSVM可以提高多时相遥感影像差异影像的分类精度及速度。

需要说明的是，在本文中，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者装置不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者装置所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括该要素的过程、方法、物品或者装置中还存在另外的相同要素。

专业人员还可以进一步意识到，结合本文中所公开的实施例描述的各示例的单元及算法步骤，能够以电子硬件、计算机软件或者二者的结合来实现，为了清楚地说明硬件和软件的可互换性，在上述说明中已经按照功能一般性地描述了各示例的组成及步骤。这些功能究竟以硬件还是软件方式来执行，取决于技术方案的特定应用和设计约束条件。专业技术人员可以对每个特定的应用来使用不同方法来实现所描述的功能，但是这种实现不应认为超出本发明的范围。软件模块可以置于随机存储器(RAM)、内存、只读存储器(ROM)、电可编程ROM、电可擦除可编程ROM、寄存器、硬盘、可移动磁盘、CD-ROM、或技术领域内所公知的任意其它形式的存储介质中。

对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的，本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下，在其它实施例中实现。因此，本发明将不会被限制于本文所示的这些实施例，而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。

Claims (6)

1.一种遥感影像半监督变化检测方法，其特征在于，包括：

获取两时相的原始光学遥感影像；

对所述原始光学遥感影像进行影像配准；

利用直方图调整方法对影像配准后的遥感影像进行辐射归一化校正；

根据辐射归一化校正后的遥感影像计算光谱角信息；

将所述辐射归一化校正后的遥感影像和所述光谱角信息组合，作为FPTSVM的输入；

采用FPTSVM方法进行学习，不断调整SVM的分类超平面，直到达到指定的学习迭代次数或前后两次迭代后的非标记样本没有变化；

利用最终的分类超平面确定影像的变化区域和非变化区域。

2.如权利要求1所述的遥感影像半监督变化检测方法，其特征在于，所述对所述原始光学遥感影像进行影像配准，具体包括：

采用ENVI14.8遥感软件对原始光学遥感影像进行几何粗校正；

利用自动匹配与三角剖分法对几何粗校正后的遥感影像进行几何精校正。

3.如权利要求1所述的遥感影像半监督变化检测方法，其特征在于，所述利用直方图调整方法对影像配准后的遥感影像进行辐射归一化校正，具体包括：

将影像配准后的遥感影像的像素灰度值分别调整为从而实现辐射归一化校正，

其中，影像配准后的遥感影像各个波段服从高斯分布，μ_f、σ_f分别为参考影像的像素均值和标准差，g_r、μ_r、σ_r分别为影像配准后的遥感影像的像素灰度值、均值和标准差。

4.如权利要求1所述的遥感影像半监督变化检测方法，其特征在于，所述根据辐射归一化校正后的遥感影像计算光谱角信息，具体包括：

根据辐射归一化校正后的遥感影像X₁和X₂计算光谱角信息X_SA，其中，

$X_{SA}\left(i,j\right)=\arccos \left(\underset{b=1}{\overset{B}{\Σ}}\left(X_{1b}\left(i,j\right)X_{2b}\left(i,j\right)\right)/\sqrt{\underset{b=1}{\overset{B}{\Σ}}{X}_{1b}^2\left(i,j\right)}\sqrt{\underset{b=1}{\overset{B}{\Σ}}{X}_{2b}^2\left(i,j\right)}\right),$

式中，B表示每一个时相的遥感影像的波段数目，(i,j)是影像的坐标，X_1b表示X₁的b波段影像，X_2b表示X₂的b波段影像。

5.如权利要求1所述的遥感影像半监督变化检测方法，其特征在于，所述将所述辐射归一化校正后的遥感影像和所述光谱角信息组合，作为FPTSVM的输入，具体包括：

将辐射归一化校正后的遥感影像X₁、X₂和光谱角信息组合X_SA组合为作为FPTSVM的输入。

6.如权利要求5所述的遥感影像半监督变化检测方法，其特征在于，所述采用FPTSVM方法进行学习，不断调整SVM的分类超平面，直到达到指定的学习迭代次数或前后两次迭代后的非标记样本没有变化，具体包括：

建立拉格朗日函数为

$\begin{array}{c}L(\mathrm{\α},{\mathrm{\α}}^{*},y^{*})=\underset{l=1}{\overset{n}{\Σ}}{\mathrm{\α}}_l+\underset{u=1}{\overset{m}{\Σ}}{\mathrm{\α}}_u^{*}-\frac{1}{2}(\underset{l=1}{\overset{n}{\Σ}}\underset{i=1}{\overset{n}{\Σ}}k\left(x_l,x_i\right)y_l{y}_i{\mathrm{\α}}_l{\mathrm{\α}}_i\\ +2\underset{l=1}{\overset{n}{\Σ}}\underset{u=1}{\overset{m}{\Σ}}k\left(x_l,x_u^{*}\right)y_l{y}_u^{*}{\mathrm{\α}}_l{\mathrm{\α}}_u^{*}+\underset{u=1}{\overset{m}{\Σ}}\underset{j=1}{\overset{m}{\Σ}}k\left(x_u^{*},x_j^{*}\right)y_u^{*}{y}_j^{*}{\mathrm{\α}}_u^{*}{\mathrm{\α}}_j^{*})\end{array}$

$Subjectto\left\{\begin{array}{c}\begin{array}{cc}0\≤{\mathrm{\α}}_l\≤C,& for1\≤l\≤n\\ 0\≤{\mathrm{\α}}_u^{*}\≤C^{*},& for1\≤l\≤m\end{array}\\ \underset{l=1}{\overset{n}{\Σ}}{y}_l{\mathrm{\α}}_l+\underset{u=1}{\overset{n}{\Σ}}{y}_u^{*}{\mathrm{\α}}_u^{*}=0\end{array}\right.$

式中，x_i为标记矢量，标记y_i∈{±1}，i＝1,2,…,n，n为标记矢量个数，-1表示没有变化类，+1表示变化类，为转导矢量，相应的伪标记为j＝1,2,…,m，m为转导矢量个数，k(·,·)表示核函数，C和C^*分别表示对有标记样本和转导样本的惩罚参数，α＝(α₁,α₂,…,α_n)和分别表示两个拉格朗日乘数；

建立分类决策函数为 $f\left(x\right)=\underset{l=1}{\overset{n}{\Σ}}{\mathrm{\α}}_l{y}_{l}k(x,x_l)+\underset{u=1}{\overset{m}{\Σ}}{\mathrm{\α}}_u^{*}{y}_{l}k(x,x_u^{*})+d;$ 式中d是一个常数，表示分类超平面的偏移量；

进行FPTSVM的归纳学习阶段：利用已知标记的样本，获取初始的超平面h₀和两个分类边界h₁、h₂，此时候选的转导样本集合式中，U₀是非标记样本集合；其中，当迭代归纳学习时，在第k+1次迭代时，标记样本集更新为：候选转导样本集合为 $S_{k+1}=\{\tilde{x}\∈S_k|h_1^{\left(k\right)}\≤f\left(\tilde{x}\right)\≤h_2^{\left(k\right)}and\tilde{x}\∉{\mathrm{\Ψ}}_k^{+}\∪{\mathrm{\Ψ}}_k^{-}\},$ 式中，和分别表示第k次迭代的转导样本集，分别表示第k次迭代的分类边界；

进行FPTSVM的转导学习阶段：将转导样本的惩罚参数设置为：式中，是初始的惩罚参数，C是初始的有标记样本的惩罚参数，k表示迭代的次数，0≤k≤G，G是最大的迭代次数；其中，在FPTSVM的转导学习阶段中，在第k次迭代前，如果任何一个转导样本的标记和前一次的标记不同，则该转导样本被重新作为未标记样本，并从当前的转导样本集中删除；在转导样本的选择过程中，如果两个转导样本的决策值相同，且相应的输入矢量相同，则删除其中一个样本；

不断迭代进行归纳学习和转导学习，直到达到指定的学习迭代次数或前后两次迭代后的非标记样本没有变化，得到最终的分类超平面。