一种基于压缩感知的高光谱混合像元分解的方法

摘要

一种基于压缩感知的高光谱混合像元分解的方法，涉及压缩感知领域与高光谱遥感领域。解决了现有采用传统的奈奎斯特采样定理对高光谱图像数据采集时，混合像元分解速度慢的问题。该方法首先输入观测矩阵Φ和压缩观测矩阵Y，利用压缩感知理论建立光谱混合模型Y＝ΦX

说明书

技术领域

本发明涉及压缩感知领域与高光谱遥感领域。

背景技术

高光谱成像技术是遥感领域发展起来的一种新型的对地观测技术，典型的硬件设备 为成像光谱仪。成像光谱仪通过分光技术，将电磁波信号分解为许多微小的、相邻的波段， 对应波段上的能量被不同的传感器接收。因此与传统的遥感成像技术相比，高光谱成像技 术具有图谱合一、光谱波段多、光谱分辨率高等特点，且在地表物质识别和分类、特征提 取等方面有较大的优势。高光谱图像的空间分辨率较低，因此不可避免的会产生混合像元， 混合像元的存在成为了高光谱遥感应用的进一步发展的障碍。如果对高光谱图像的混合像 元分解得到各端元所对应的丰度，即对高光谱遥感数据进行“盲”分解将有非常重要的意 义。

盲源分离技术是指在源信号和混合方式均未知的情况下，根据源信号的统计特性，仅 由观测到的混合信号恢复出源信号的过程。“盲”有两种解释：一种是除观测信号外，其 他的信息均未知，称为“全盲”；另一种是信号或者混合系统的某些先验信息已知，称为 “半盲”。针对混合像元分解的问题，盲源分离技术同时获得端元光谱和端元丰度，即实 现混合像元的盲分解。在处理混合像元分解问题之前，需要建立光谱混合模型。光谱混合 方式主要包括线性混合和非线性混合，线性光谱混合模型是目前应用最广泛的混合像元分 解模型，其表达式为：

$X=\underset{i=1}{\overset{P}{\Sigma }}{A}_i{s}_i=\mathrm{AS}---(1-1)$

式中，X＝[x₁,...x_l,...x_L]^T∈R^L×1为一个像元L个波段的混合像元光谱，是已知观测量； A∈R^L×P为端元光谱矩阵，其中每一列对应一个端元的光谱向量；P为端元数目，S∈R^P×1为该混合像元的丰度向量。丰度向量需要满足两个约束条件：非负性约束、和为1约束。 即

s_i≥0,i∈{1,2,...,P}    (2-1)

$\underset{i=1}{\overset{P}{\Sigma }}{s}_i=1,i\in \{\mathrm{1,2},...,P\}---(3-1)$

若高光谱图像具有N个像元，则公式(1-1)可以进一步扩展为公式(4-1)。

X＝AS    (4-1)

公式(4-1)中X∈R^L×N,其列向量表示每个像元在L个波段的混合光谱向量，S∈R^P×N构成了端元丰度矩阵。端元丰度矩阵满足以下约束条件：

s_i,j≥0,i∈{1,2,...,P},j∈{1,2,...,N}    (5-1)

$\underset{i=1}{\overset{P}{\Sigma }}{s}_{i,j}=1,j\in \{\mathrm{1,2},...,N\}---(6-1)$

基于盲源分离技术的高光谱图像混合像元分解问题即在端元矩阵和丰度矩阵均未知 的情况下，利用混合像元光谱的观测值来恢复端元光谱和端元丰度。

高光谱图像的波段数多的特点使高光谱图像的数据量巨大，若采用传统的奈奎斯特 采样定理对数据采集，将会对硬件的采样速度、传输速度以及存储能力产生巨大的压力。 压缩感知理论的提出解决了这一难题，该理论指出当信号满足稀疏性或者可压缩时，能以 远低于奈奎斯特采样率的速度对信号进行全局观测，然后采用压缩感知重构算法恢复源信 号。高光谱图像在小波域是稀疏的，因此可以将压缩感知理论应用到对高光谱图像的采样 中。具体的采样过程可以表示为下式：

Y＝ΦX^T＝Φ(AS)^T

在上述公式中，Φ∈R^M×N为观测矩阵,X∈R^L×N为混合像元光谱矩阵，Y∈R^M×L为压 缩观测矩阵。

发明内容

本发明是为了解决现有采用传统的奈奎斯特采样定理对高光谱图像数据采集时，混合 像元分解速度慢的问题，本发明提供了一种基于压缩感知的高光谱混合像元分解的方法。

一种基于压缩感知的高光谱混合像元分解的方法，该方法为，

步骤一，输入观测矩阵Φ和压缩观测矩阵Y，利用压缩感知理论，建立光谱混合模 型：

Y＝ΦX^T＝Φ(AS)^T    (1)

Φ∈R^M×N为M×N的观测矩阵，R为实数，

X∈R^L×N为L×N的混合像元光谱矩阵，

Y∈R^M×L为M×L的压缩观测矩阵，

S∈R^P×N为P×N的端元丰度矩阵，

A∈R^L×P为L×P的端元光谱矩阵，

步骤二，初始化，随机选取一个端元光谱矩阵A作为端元光谱矩阵A的估计值 且为L×P的矩阵，

令端元丰度矩阵S的估计值：

$\hat{S}={\left(\begin{array}{ccc}0& ...& 0\\ & ...& \\ 0& & 0\end{array}\right)}_{N\times P}---\left(2\right)$

其中，为N×P的矩阵，

步骤三，令迭代次数变量t的初始值为1；

步骤四，对端元丰度矩阵S的估计值和端元光谱矩阵A的估计值进行迭代处 理；

步骤五，在公式(1)两边同时乘以的伪逆，则公式(1)变形为

Y₁＝ΦS^T    (3)，

步骤六，设小波基为Θ，Σ为端元丰度的稀疏系数，则公式(3)变形为

Y₁＝Φ(ΘΣ)^T    (4)，

采用BP算法求解下式：

argmin||Σ||₁ s.t. Y₁＝Φ(ΘΣ)^T    (5)，

获得端元丰度矩阵S的估计值

$\hat{S}=---\mathrm{\Theta \Sigma }\left(6\right),$

步骤七，在概率密度函数为P(S_k)∝α exp(α|S_k|)的条件下，归一化端元光谱矩阵 A的估计值的列向量，更新端元光谱矩阵A的估计值，

$\hat{A}=A_1+\lambda \{-A_1(B{\hat{S}}^T+I)\}---\left(7\right),$

其中，A₁为上一次迭代得到的端元光谱矩阵的估计值，λ为迭代步长，I为单位矩 阵，B为向量集合，且B＝{B₁，B₂，B₃，......B_k}，为第 k个源信号的先验概率分布，k为整数，α为大于零的实数，

步骤八，如果相邻两次获得的端元光谱矩阵A的估计值中对应的每个元素的绝对 值之差小于0.1时，则停止迭代，执行步骤九，否则，令t＝t+1，且返回步骤四；

步骤九，输出端元丰度矩阵完成对高光谱混合像元的分解。

本发明带来的有益效果是，对端元丰度矩阵S的评价指标采用平均信噪比SNR_avg实 现，

${\mathrm{SNR}}_{\mathrm{avg}}=\frac{1}{P}\underset{k=1}{\overset{P}{\Sigma }}20\times {\log }_{10}\left(\frac{\left|\right|S_k{\left|\right|}_2}{{\left|\right|S_k-{\hat{S}}_k\left|\right|}_2}\right)---\left(8\right)$

实验中采用的高光谱图像的端元数目P＝3，波段数L＝40，高光谱的像元数N＝32× 32。在压缩感知采样时，压缩率为0.8，对混合像元分解的端元丰度矩阵的平均信噪比为 17.78dB，可以较好的实现分离，使得混合像元分解速度提高了40％以上。

具体实施方式

具体实施方式一：本实施方式所述的一种基于压缩感知的高光谱混合像元分解的方 法，该方法为，

步骤一，输入观测矩阵Φ和压缩观测矩阵Y，利用压缩感知理论，建立光谱混合模 型：

Y＝ΦX^T＝Φ(AS)^T    (1)

Φ∈R^M×N为M×N的观测矩阵，R为实数，

X∈R^L×N为L×N的混合像元光谱矩阵，

Y∈R^M×L为M×L的压缩观测矩阵，

S∈R^P×N为P×N的端元丰度矩阵，

A∈R^L×P为L×P的端元光谱矩阵，

步骤二，初始化，随机选取一个端元光谱矩阵A作为端元光谱矩阵A的估计值 且为L×P的矩阵，

令端元丰度矩阵S的估计值：

$\hat{S}={\left(\begin{array}{ccc}0& ...& 0\\ & ...& \\ 0& & 0\end{array}\right)}_{N\times P}---\left(2\right)$

其中，为N×P的矩阵，

步骤三，令迭代次数变量t的初始值为1；

步骤四，对端元丰度矩阵S的估计值和端元光谱矩阵A的估计值进行迭代处 理；

步骤五，在公式(1)两边同时乘以的伪逆，则公式(1)变形为

Y₁＝ΦS^T    (3)，

步骤六，设小波基为Θ，Σ为端元丰度的稀疏系数，则公式(3)变形为

Y₁＝Φ(ΘΣ)^T    (4)，

采用BP算法求解下式：

argmin||Σ||₁ s.t. Y₁＝Φ(ΘΣ)^T    (5)，

获得端元丰度矩阵S的估计值

$\hat{S}=---\mathrm{\Theta \Sigma }\left(6\right),$

步骤七，在概率密度函数为P(S_k)∝α exp(α|S_k|)的条件下，归一化端元光谱矩阵 A的估计值的列向量，更新端元光谱矩阵A的估计值

$\hat{A}=A_1+\lambda \{-A_1(B{\hat{S}}^T+I)\}---\left(7\right),$

其中，A₁为上一次迭代得到的端元光谱矩阵的估计值，λ为迭代步长，I为单位矩 阵，B为向量集合，且B＝{B₁，B₂，B₃，......B_k}，为第 k个源信号的先验概率分布，k为整数，α为大于零的实数，

步骤八，如果相邻两次获得的端元光谱矩阵A的估计值中对应的每个元素的绝对 值之差小于0.1时，则停止迭代，执行步骤九，否则，令t＝t+1，且返回步骤四；

步骤九，输出端元丰度矩阵完成对高光谱混合像元的分解。