CN109671019B - 一种基于多目标优化算法和稀疏表达的遥感影像亚像元制图方法 - Google Patents

Abstract

本发明涉及一种基于多目标优化算法和稀疏表达的遥感影像亚像元制图方法。本发明结合多目标优化理论，将亚像元制图问题转化为对数据保真项、稀疏约束项、空间先验项的多目标优化问题，在稀疏约束项上直接使用非凸的L0范数建模，并构建基于亚像元集群块的亚像元空间分布模式字典，实现亚像元级空间先验建模。然后利用多目标优化算法强大的寻优能力来同时优化三个目标项，求解使该模型达到最优的一组解。本发明能够解决现有的稀疏亚像元制图方法中L1范数稀疏性建模不精确问题，以及最优权重参数选取困难问题。本发明可以有效提高亚像元制图的适用性和精度。

Description

一种基于多目标优化算法和稀疏表达的遥感影像亚像元制图 方法

技术领域

本发明属于遥感影像亚像元制图领域，特别涉及一种基于多目标优化算法和稀疏表达的遥感影像亚像元制图方法。

背景技术

由于传感器瞬时视场角的影响以及地物分布的复杂多样性，遥感影像中普遍存在着混合像元。混合像元是指在一个像元中包含多种地物类别，并且各种地物具有不同的光谱响应特征。混合像元问题不仅影响遥感影像分类精度，而且是遥感技术定量发展的重要障碍，限制了遥感影像的应用范围。为了解决混合像元问题，学者提出了混合像元光谱分解技术，该方法可以估计混合像元中各种地物(端元)和地物所占比例(丰度)，但无法确定各个组成地物在混合像元内部的空间分布，因而造成遥感影像空间细节信息的丢失。为了进一步解译混合像元内部各种地物的空间分布情况，必须借助亚像元制图技术。

亚像元制图是由Atkinson于1997年首先提出，即将混合像元分割成更小单元的亚像元，并赋予这些亚像元特定的地物类别。亚像元制图是以混合像元光谱分解得到的丰度图为基础，在亚像元级上确定一个像元中不同的地物在空间上的可能分布，从而获得亚像元尺度的地物分类图。亚像元制图技术作为一种增强影像分类分辨率的有效手段，可以将影像的分类精度提升到亚像元级。在实际应用中，亚像元定位技术已成功运用于土地覆盖制图，湖泊边界提取，变化检测，目标探测等。

亚像元制图方法的理论基础为空间相关性假设理论，即认为距离较近的亚像元比距离较远亚像元更优可能属于同一种地物类别。由于亚像元制图问题属于病态的逆问题，在缺乏地物空间分布先验知识的情况下，空间相关性理论使得亚像元制图问题能够得到有效解决。已有的亚像元制图方法主要包括基于空间吸引力模型、像元置换算法、线性优化算法、差分进化算法、神经网络模型的制图方法、地统计学的亚像元制图方法、马尔科夫随机场的亚像元制图方法、基于最大后验概率模型的亚像元制图方法以及多智能体等亚像元制图方法等。虽然已有的亚像元制图方法能够获取地物的亚像元空间分布信息，然而由于制图过程中仅仅依赖空间相关性假设以及光谱分解获取的丰度图像作为制图依据，纳入考虑的仅为混合像元与混合像元之间、亚像元与混合像元之间的单一空间关系，所以上述方法获取的亚像元空间分布模式缺乏多样性，无法为欠定的亚像元制图问题提供足够的空间先验信息。

基于稀疏表达理论的亚像元制图方法通过在亚像元制图块字典中寻找最优的亚像元制图块组合，以实现丰度影像的稀疏表达。基于稀疏表达理论的亚像元制图方法从亚像元制图块的角度寻找最佳的匹配，考虑了亚像元图块自身的空间分布模式，较传统仅考虑亚像元级别的空间相关性的制图方法，能够从更高层面上考虑其空间分布模式信息。而且，稀疏表达理论通过稀疏约束能够使欠定的亚像元制图问题获得稳定解。

但是，稀疏亚像元制图方法仍然存在以下难点问题：(1)L1范数稀疏度建模不精确。在原始的稀疏亚像元制图模型中，稀疏约束项是使用L0范数进行精确建模，然而由于L0范数的求解属于非凸的NP-hard问题，传统梯度下降法无法求解。所以现有稀疏亚像元制图方法是采用L1范数对稀疏项进行近似建模。虽然在满足等距限制定理条件下，L0范数可以等价转换为L1范数进行求解。但是，由于遥感影像的复杂性，该条件往往难以完全满足，导致亚像元空间分布模式的稀疏性建模不够准确；(2)最优的权重参数选取困难。由于稀疏亚像元制图模型中存在数据保真项、稀疏约束项、空间约束项等多个优化项，所以现有稀疏亚像元制图方法引入了两个权重参数将多个优化项合并为一个目标函数，其中的权重参数的选取极大地影响着亚像元制图精度。虽然已有学者提出权重参数自适应选择策略，但是，其仍然存在自适应能力低、对初值敏感等问题。因此，稀疏亚像元制图方法仍存在最优权重参数选择困难问题。

综合分析当前的亚像元定位研究之后发现，目前亟需发明一种能够直接优化L0范数、无需选择权重参数的稀疏亚像元制图方法。多目标优化算法具有强大的优化能力，不仅能够求解非凸的L0范数，避免L1范数建模带来的误差，而且能够同时优化多个目标项，避免引入权重参数。因此非常适合与稀疏表达理论相结合，解决遥感影像亚像元制图的上述问题。

发明内容

本发明的目的在于提出一种基于多目标优化算法和稀疏表达的遥感影像亚像元制图方法。

本发明所提供的这种基于多目标优化算法和稀疏表达的遥感影像亚像元制图方法，对于一幅低分辨率的遥感影像，经过光谱分解得到低分辨率丰度图，然后学习过完备的空间分布模式先验字典，构建多目标稀疏亚像元制图模型，并采用多目标优化算法求解使制图模型最优的稀疏系数和亚像元级丰度图，最后按照类别决策策略投票产生遥感影像亚像元级分类图。

在本发明中我们所提出的基于多目标优化算法和稀疏表达的遥感影像亚像元制图方法具有以下三个显著特点。一是其学习了过完备的亚像元空间分布模式字典。过完备字典中包含了许多与输入的遥感影像空间分布类似的亚像元制图块，可以重构出满足稀疏性约束的、与真实地表分布接近的亚像元空间分布模式；第二个显著特点是多目标稀疏亚像元制图模型采用L0范数来对解的稀疏度进行准确建模，从而可以避免已有方法中的L1范数引起的潜在误差。例如当过完备字典不满足等距限制定理条件时，L1范数只能得到L0范数的近似解；第三个显著特点是采用多目标优化算法来求解制图模型的最优解，从而可以直接优化L0范数，并且无需人工选择权重参数，实现自动化的、高精度的亚像元制图。总而言之，我们所提出的基于多目标优化算法和稀疏表达的遥感影像亚像元制图方法不仅能够继承稀疏亚像元制图的优点，而且解决了稀疏亚像元制图的问题，具有以下优点：自动化程度高，方法内没有敏感参数，无需人工调参；稀疏项建模准确，方法能够直接优化原始稀疏亚像元制图中的L0范数，避免了凸放松带来的潜在误差；制图精度高，多目标优化方法具有强大的全局寻优能力，得到的解的质量有保证。

本发明提供基于多目标优化算法和稀疏表达的遥感影像亚像元制图方法，实现步骤如下：

步骤1，获取低分辨率的遥感影像，经过光谱分解得到遥感影像中各类地物的低分辨率丰度图；

步骤2，学习过完备的空间分布模式先验字典；

步骤3，构建多目标稀疏亚像元制图模型，并采用多目标优化算法求解使制图模型最优的稀疏系数和亚像元级丰度图；

步骤4，根据各类地物对应的亚像元级丰度图，按照类别决策策略投票产生遥感影像亚像元级分类图。

进一步的，步骤1的具体实现包括如下子步骤，

步骤1.1，使用手持Nuance近红外光谱仪获取一幅大小为[m,n]的低分辨率遥感影像；

步骤1.2，采集遥感影像中C类纯净地物的样本(Region Of Interest，ROI)，将样本的光谱曲线作为端元光谱；

步骤1.3，输入低分辨率遥感影像和步骤2得到的ROI样本文件，利用全约束最小二乘算法(Fully Constrained Least Squares，FCLS)自动求解并输出C类地物的低分辨率丰度影像Y＝{y¹,y²,...,y^C}。

进一步的，步骤2的具体实现包括如下子步骤，

步骤2.1，采用三次样条插值方法对低分辨率丰度影像Y＝{y¹,y²,...,y^C}进行尺度为S的上采样，得到大小为[M,N]的初始亚像元级丰度图Z＝{z¹,z²,...,z^C}，其中C为地物类别的个数，M＝m×S,N＝n×S；

步骤2.2，将初始亚像元级丰度图Z＝{z¹,z²,...,z^C}分割成C×m×n个互不重叠的、大小为[S,S]的影像块，作为空间分布模式的学习样本集；

步骤2.3，输入样本集，采用K-SVD算法学习影像块的空间分布模式，得到大小为[S²,b],b＞S²的过完备空间分布模式字典，其中，b是过完备字典中的原子个数，S代表亚像元制图尺度。

进一步的，步骤3中构建的多目标稀疏亚像元制图模型包括三个优化项：保真项、稀疏项、空间先验项，亚像元制图的结果是使该模型最小化的解，

min F(α)＝minimize{f₁(α),f₂(α),f₃(α)}

其中，y^c是地物c的低分辨率丰度图，α是地物c的高分辨率丰度影像，||α||₀是α的L0范数，A是过完备空间分布模式字典，D是大小为[1,S²]的降采样矩阵，其中每个元素的值都是1/S²，U(Aα)是对亚像元级丰度图Aα施加的空间约束。

进一步的，步骤3采用多目标优化算法求解使制图模型最优的稀疏系数和亚像元级丰度图，具体实现方式如下，

步骤3.1，采用基于稀疏系数矩阵的个体编码方式：在[minα^c,minα^c]范围内随机初始化低分辨率丰度图y^c的稀疏系数矩阵α^c作为一个个体，其中minα^c＝-1,maxα^c＝1；依次初始化NP个个体，形成种群

种群中的第k个个体可表达为

其中m,n分别代表低分辨率丰度图的行列数，b是过完备字典中的原子个数，S代表亚像元制图尺度，即一个像元划分为S×S个亚像元；

步骤3.2，均匀初始化多目标优化算法的NP个权重参数λ^j,j＝{1,2,...,NP}，每个权重向量与一个个体相对应；计算任意两个权重向量的相似度，统计每个权重向量λ^k的与其最相近的T个权重向量

得到个体

的T个邻域个体

步骤3.3，全局搜索：针对种群中的第k个个体

随机从其邻域权重向量

中选择三个权重向量λ^p1,λ^p2,λ^p3；然后，基于上述三个权重向量所对应的个体即

采用DE/rand/1策略产生子代个体

其中F和CR分别代表尺度因子和交叉概率；

其中i＝1,2,…,m×n×b，i_rand是1～NP的一个随机正整数；

步骤3.4，局部搜索：对子代个体

进行局部高斯变异，产生新的子代个体

其中

是均值为

标准差为δ的高斯分布，p_g是高斯变异的概率；

步骤3.5，稀疏变异：对子代个体

进行稀疏变异，产生更稀疏的系数矩阵

其中prob是控制稀疏变异的概率，itr是当前演化代数，maxitr是最大演化次数，

prob＝1-d^itr/maxitr

步骤3.6，采用基于子适应度的种群更新策略比较个体

和其子代个体

的优劣，将质量更好的个体保留到新种群POP^new中，对种群POP的所有个体重复步骤3.3-3.6；

步骤3.7，用新种群POP^new替换原始种群POP，进入新一轮演化(迭代步骤3.3-3.7)，直到演化代数达到maxitr代，输出一组最优个体

即可得到一组亚像元级丰度图

手动选择一幅目视效果好的

作为地物c的亚像元级丰度图；

步骤3.8，对所有地物的低分辨率丰度图Y＝{y¹,y²,...,y^C}，利用多目标优化算法求解出对应的稀疏系数矩阵和亚像元级丰度图，得到C个类别的高分辨率分布影像X＝{x¹,x²,...,x^C}。

进一步的，步骤3.6中采用基于子适应度的种群更新策略比较个体

和其子代个体

的优劣，将质量更好的个体保留到新种群POP^new中，具体实现方式如下，

步骤3.6.1，个体评价：由个体的编码策略可知，种群中的第k个个体

可分为m×n个部分，其中的第i个部分

代表低分辨率丰度影像y^c的第i个像元的b个稀疏系数，根据多目标稀疏亚像元制图模型计算每个部分所对应的目标函数值(保真项、稀疏项和空间先验项)作为子适应度；其中

代表种群中第k个个体的第i个部分，

代表第j个子问题的权重参数值，z^*代表参考值，Q代表拉普拉斯矩阵，

步骤3.6.2，基于子适应度的个体比较策略及种群更新策略：利用子代个体

更新与其对应的个体

的邻域个体

即分别比较

和

对应部分的子适应度SubFit(·)，然后保留其中具有更小子适应度的部分，得到新的个体

保留到新种群POP^new中，

步骤3.6.3，更新参考值z^*：参考值z^*是子适应度

当前所能得到的最优值，在初始阶段，z^*被指定为一个较大的值，利用迭代过程中计算的子适应度更新参考值z^*，即针对每个

如果

那么

进一步的，步骤4中，根据步骤3.8得到的亚像元级丰度图，按照“胜者为王”投票策略得到亚像元级的地物分类图，即亚像元z_i,i＝{1,2,…,M×N}的类别可以确定为

本发明方法具有以下显著效果：

1.本方法结合稀疏表达理论来解决遥感影像亚像元制图，基于亚像元制图块能够重构出更多样复杂的亚像元空间细节；

2.多目标优化算法能够同时优化亚像元制图模型的多个目标项，包括非凸的L0稀疏项，并且在演化过程能够使多个目标项之间信息互补。方法自动化程度高，无需人工调参，就能求解得到高精度的亚像元制图结果，能够提高亚像元制图的适用性和精度；

3.多目标优化算法中既有全局优化策略，也有局部优化策略，还有稀疏变异策略，所以本方法寻优能力强，收敛速度快，在15代以内就能找到最优的亚像元制图结果；

附图说明

图1是本发明的实施例1中输入的低分辨率遥感影像(a)和光谱解混得到的低分辨率丰度图(b)。

图2是本发明实施例1的步骤3.1中采用的基于稀疏系数矩阵的个体编码方式。

图3是本发明实施例1的步骤3.6采用的基于子适应度的种群更新策略示意图。

图4是本发明实施例1的亚像元制图结果——亚像元级分类图。

具体实施方式

下面通过实施例，进一步阐明本发明的突出特点和显著进步，仅在于说明本发明而决不限制本发明。

实施例1

(一)使用手持Nuance近红外光谱仪获取一幅24个波段、大小为[50,50]的低分辨率遥感影像，如图1(a)所示，其中包含3种地物(枯草、新鲜草、背景地物)。

1.1.利用ENVI 5.1的ROI tool采集枯草、新鲜草、背景地物的样本，将样本的光谱曲线作为端元光谱。

1.2.添加ENVI 5.1扩展工具包——全约束最小二乘算法(Fully ConstrainedLeast Squares，FCLS)Spectral Unmixing，输入低分辨率遥感影像和ROI样本文件，自动求解并输出各类地物的低分辨率丰度影像Y＝{y¹,y²,y³}，如图1(b)所示。

(二)学习过完备的空间分布模式字典。

2.1.采用三次样条插值方法对低分辨率丰度影像Y＝{y¹,y²,y³}进行尺度为3的上采样，得到大小为[150,150]的初始亚像元级丰度图Z＝{z¹,z²,z³}。

2.2.将初始亚像元级丰度图Z＝{z¹,z²,z³}分割成3×50×50个互不重叠的、大小为[3,3]的影像块，作为空间分布模式的学习样本集。

2.3.输入样本集，采用K-SVD算法学习影像块的空间分布模式，得到大小为[9,100]的过完备字典，包含100种空间分布模式。

(三)采用多目标优化算法求解使制图模型最优的稀疏系数和亚像元级丰度图。

3.1.采用基于稀疏系数矩阵的个体编码方式：在[-1,1]范围内随机初始化低分辨率丰度图y^c的稀疏系数矩阵α^c作为一个个体。依次初始化50个个体，形成种群

种群中的第k个个体可表达为

其原理图如图2所示。

3.2.均匀初始化多目标优化算法的50个权重参数λ^j,j＝{1,2,…,50}，每个权重向量与一个个体相对应。计算任意两个权重向量的相似度，统计每个权重向量λ^k的与其最相近的8个权重向量

得到个体

的8个邻域个体

3.3.全局搜索：针对种群中的第k个个体

随机从其邻域权重向量

中选择三个权重向量λ^p1,λ^p2,λ^p3；然后，基于上述三个权重向量所对应的个体即

采用DE/rand/1策略产生子代个体

其中F＝0.8和CR＝0.1分别代表尺度因子和交叉概率。

其中i＝1,2,…,50×50×100，i_rand是1～50的一个随机正整数。

3.4.局部搜索：对子代个体

进行局部高斯变异，产生新的子代个体

其中

是均值为

标准差为0.1的高斯分布。

3.5.稀疏变异：对子代个体

进行稀疏变异，产生更稀疏的系数矩阵

其中itr是当前演化代数。

prob＝1-0.001^itr/15

3.6.采用基于子适应度的种群更新策略比较个体

和其子代个体

的优劣，将质量更好的个体保留到新种群POP^new中。对种群POP的所有个体重复步骤3.3-3.6。

3.6.1.个体评价：由个体的编码策略可知，种群中的第k个个体

可分为50×50个部分，其中的第i个部分

代表低分辨率丰度影像y^c的第i个像元的100个稀疏系数。根据多目标稀疏亚像元制图模型计算每个部分所对应的目标函数值(保真项、稀疏项和空间先验项)作为子适应度。其中

代表种群中第k个个体的第i个部分，

代表第j个子问题的权重参数值，z^*代表参考值，Q代表拉普拉斯矩阵。

3.6.2.基于子适应度的个体比较策略及种群更新策略：利用子代个体

更新与其对应的个体

的邻域个体

即分别比较

和

对应部分的子适应度SubFit(·)，然后保留其中具有更小子适应度的部分，其原理图如图3所示。得到新的个体

保留到新种群POP^new中。

3.6.3.更新参考值z^*：参考值z^*是子适应度

当前所能得到的最优值。在初始阶段，z^*被指定为一个较大的值。利用迭代过程中计算的子适应度更新参考值z^*，即针对每个

如果

那么

3.7.用新种群POP^new替换原始种群POP，进入新一轮演化(迭代步骤3.3-3.7)。直到演化代数达到15代，输出一组最优个体

即可得到一组亚像元级丰度图

手动选择一幅目视效果好的

作为地物c的亚像元级丰度图。

3.8.对所有地物的低分辨率丰度图Y＝{y¹,y²,y³}都用多目标优化算法求解出对应的稀疏系数矩阵和亚像元级丰度图，求解的亚像元级丰度图即地物类别对应的高分辨率分布影像X＝{x¹,x²,x³}。

3.9.在得到亚像元级丰度图X＝{x¹,x²,x³}之后，根据“胜者为王”投票策略得到亚像元级的地物分类图，即亚像元z_i,i＝{1,2,...,150×150}的类别可以确定为

最终的制图效果如图4所示。

本文中所描述的具体实施例仅仅是对本发明精神作举例说明。本发明所属技术领域的技术人员可以对所描述的具体实施例做各种各样的修改或补充或采用类似的方式替代，但并不会偏离本发明的精神或者超越所附权利要求书所定义的范围。

Claims (1)

1.一种基于多目标优化算法和稀疏表达的遥感影像亚像元制图方法，其特征在于，包括以下步骤：

步骤1，获取低分辨率的遥感影像，经过光谱分解得到遥感影像中各类地物的低分辨率丰度图；

步骤1的具体实现包括如下子步骤，

步骤1.1，使用手持Nuance近红外光谱仪获取一幅大小为[m,n]的低分辨率遥感影像；

步骤1.2，采集遥感影像中C类纯净地物的样本，将样本的光谱曲线作为端元光谱；

步骤1.3，输入低分辨率遥感影像和步骤2得到的ROI样本文件，利用全约束最小二乘算法自动求解并输出C类地物的低分辨率丰度影像Y＝{y¹,y²,...,y^C}；

步骤2，学习过完备的空间分布模式先验字典；

步骤2的具体实现包括如下子步骤，

步骤2.1，采用三次样条插值方法对低分辨率丰度影像Y＝{y¹,y²,...,y^C}进行尺度为S的上采样，得到大小为[M,N]的初始亚像元级丰度图Z＝{z¹,z²,...,z^C}，其中C为地物类别的个数，M＝m×S,N＝n×S；

步骤2.2，将初始亚像元级丰度图Z＝{z¹,z²,...,z^C}分割成C×m×n个互不重叠的、大小为[S,S]的影像块，作为空间分布模式的学习样本集；

步骤2.3，输入样本集，采用K-SVD算法学习影像块的空间分布模式，得到大小为[S²,b],b>S²的过完备空间分布模式字典，其中，b是过完备字典中的原子个数，S代表亚像元制图尺度；

步骤3，构建多目标稀疏亚像元制图模型，并采用多目标优化算法求解使制图模型最优的稀疏系数和亚像元级丰度图；

步骤3中构建的多目标稀疏亚像元制图模型包括三个优化项：保真项、稀疏项、空间先验项，亚像元制图的结果是使该模型最小化的解，

min F(α)＝minimize{f₁(α),f₂(α),f₃(α)}

其中，y^c是地物c的低分辨率丰度图，α是地物c的高分辨率丰度影像，||α||₀是α的L0范数，A是过完备空间分布模式字典，D是大小为[1,S2]的降采样矩阵，其中每个元素的值都是1/S2，U(Aα)是对亚像元级丰度图Aα施加的空间约束；

步骤4，根据各类地物对应的亚像元级丰度图，按照类别决策策略投票产生遥感影像亚像元级分类图。

Images