基于3-D小波变换和稀疏张量的高光谱图像特征抽取方法

摘要

本发明公开了基于3-D小波变换和稀疏张量的高光谱图像特征抽取方法，包括以下步骤：步骤(1)：采用数据归一化方法平衡数据本身对判别特征抽取的影响；步骤(2)：采用3-D离散小波变换从归一化后的数据中抽取光谱域和空间域特征；步骤(3)：通过将小波变换特征表示成二阶特征张量形式，保持了特征之间良好的结构相关性；步骤(4)：通过稀疏张量判别方法实现特征稀疏化；步骤(5)：将稀疏后的特征重新表示成以为向量形式。本发明能有效提高整个分类系统的分类精度。

说明书

技术领域

本发明属于高光谱图像数据处理与应用领域，尤其涉及一种基于3-D小波变换和稀疏张 量的高光谱图像特征抽取方法。

背景技术

高光谱图像将目标的空间信息和光谱信息融为一体，对目标空间成像的同时，对每个空 间像元采集几十乃至几百连续的波段的光谱数据。高光谱图像在识别与精确分类方面具有突 出的优势，已被广泛成功应用于医学诊断，农业检测，矿物探测，环境监测等领域中。

高光谱数据存在数据量大，冗余度高和维数灾难等问题，要实现高光谱图像分类问题， 首先需要进行判别特征的提取。现有的特征提取方法，如DWT(Discrete Wavelet Transform)、 EMPs(Extended morphological profiles)、EAPs(Extended attribute profiles)等方法，对高光 谱图像的所有波段或者几个主成分进行特征变换，然后把得到的特征整合成一个长向量，不 仅造成特征向量维度过高，还丢失了大量的结构信息。

发明内容

本发明的目的就是为了解决上述问题，提供一种基于3-D小波变换和稀疏张量的高光谱 图象特征抽取方法，本发明能有效提高整个分类系统的分类精度。

为了实现上述目的，本发明采用如下技术方案：

基于3-D小波变换和稀疏张量判别分析的高光谱特征抽取方法，包括以下步骤：

首先，采用数据归一化方法平衡数据本身对判别特征抽取的影响，然后，采用3-D离散 小波变换从归一化后的数据中抽取光谱域和空间域特征；再次，通过将小波变换特征表示成 二阶特征张量形式，保持了特征之间良好的结构相关性，最后，通过稀疏张量判别方法实现 特征稀疏化，并将稀疏后的特征重新表示成以为向量形式。

基于3-D小波变换和稀疏张量判别分析的高光谱特征抽取方法，包括以下步骤：

步骤(1)：对给定的高光谱图像数据立方体C进行数据归一化处理，得到归一化后的高 光谱图像数据立方体C_n；

步骤(2)：对步骤(1)归一化后的高光谱图像数据立方体C_n进行3-D离散小波变换， 得到不同尺度下的小波变换系数立方体C_k(k＝1,2,...15)；

步骤(3)：张量表示。基于步骤(2)中所有的小波变换系数立方体C_k(k＝1,2,...15)， 对任意位置(i,j)的像元，分别从C_k(k＝1,2,...15)中抽取对应位置的小波变换系数向量 k＝1,2,...,15，然后在以(i,j)为中心的3×3邻域内对C_k(i,j,·)取均值，为位置 (i,j)的像元构建二阶特征张量P是波段数；

步骤(4)：采用稀疏张量判别分析法，对步骤(3)中的二阶特征张量T_i,j进行稀疏化， 并将二阶特征张量投影到低维特征张量L₁≤P,L₂≤15；

步骤(5)：将步骤(4)的低维特征向量重新表示成向量形式。

所述步骤(1)中数据归一化的方法为：

给定高光谱图像数据立方体X是高光谱图像的宽度、Y表示高光谱图像高度， P是高光谱图像的波段数，表示空间坐标为(i,j)的样本(光谱向量)，归一化 方法如下：

$\left(\begin{array}{c}C(i,j,k)=\frac{C(i,j,k)}{{\sigma }_k},i=1,...,X,j=1,...,Y,k=1,...,P\\ {\mu }_k=\frac{1}{X*Y}\underset{i=1}{\overset{X}{\Sigma }}\underset{j=1}{\overset{Y}{\Sigma }}C(i,j,k)\\ {\sigma }_k^2=\frac{1}{X*Y}\underset{i=1}{\overset{X}{\Sigma }}\underset{j=1}{\overset{Y}{\Sigma }}{[C(i,j,k)-{\mu }_k]}^2\end{array}\right)---\left(1\right)$

C(i,j,k)是样本C(i,j,·)的第k个特征，μ_k是所有样本第k个特征的均值，是所有样本第 k个特征的方差，C(i,j,k)是归一化后的新特征值，其中i＝1,2,...,X表示样本在图像中水平方 向的位置，j＝1,2,...,Y表示样本在图像中垂直方向上的位置，k＝1,2,...,P表示第k个波段。

所述步骤(2)的3-D离散小波变换为：

步骤(2-1)：将任意样本对应的光谱向量视为离散信号，离散信号的快速小波变换，公 式如下：

$\left(\begin{array}{c}{c}_{j,k}=\underset{m}{\Sigma }{h}_0(m-k)c_{j-1,m}\\ d_{j,k}=\underset{m}{\Sigma }{h}_1(m-k)c_{j-1,m}\end{array}\right)---\left(2\right)$

其中，c_j,k是尺度系数c_j的第k个系数，d_j,k是小波系数d_j的第k个系数，j表示小波变 换的阶数，h₀(·)为低通滤波器，h₁(·)为高通滤波器，m是平移因子，k∈[1,l]，l是c_j的长度；

步骤(2-2)：选择Haar二进小波作为快速小波变换的母波，则步骤(2-1)中的高通滤波 器为$h_0\left(k\right)=(\frac{1}{\sqrt{2}},\frac{1}{\sqrt{2}}),$低通滤波器为$h_1\left(k\right)=(-\frac{1}{\sqrt{2}},\frac{1}{\sqrt{2}}).$

步骤(2-3)：设计两层三阶小波变换，令初始尺度系数c₀＝C(i,j,·)，(i＝1,...,X j＝1,...,Y)， 依次对步骤(1)归一化后的高光谱图像数据立方体C_n中的每个像素点对应的光谱向量进行 离散小波变换，得到15个小波变换系数立方体C_k(k＝1,...,15)。

所述步骤(3)张量表示的方法为：

步骤(3-1)：基于步骤(2)中小波变换系数立方体C_k(k＝1,...,15)，位置(i,j)像素点对应 15个同样长度的小波变换系数向量C₁(i,j,·),C₂(i,j,·),...,C₁₅(i,j,·)，以(i,j)像素点为中心在此 像素点的3×3邻域内取平均：

${\bar{C}}_k(i,j,·)=\frac{1}{9}\underset{a=i-1}{\overset{i+1}{\Sigma }}\underset{b=j-1}{\overset{j+1}{\Sigma }}\left|C_k(a,b,·)\right|k=1,...,15---\left(3\right)$

k＝1,...,15是3×3邻域的均值，(i,j),i＝1,...,X,j＝1,...,Y，a表示3×3邻域水平方向 坐标值，b是3×3邻域垂直方向坐标值，C_k(a,b,·)是对应3×3模板中位置(a,b)的小波变换系 数向量。

步骤(3-2)：为位置为(i,j)的像素点构建一个二阶特征张量

T_i,j＝[C₁(i,j,·),C₂(i,j,·),...,C₁₅(i,j,·)]        (4)

所述步骤(4)稀疏张量判别分析方法为：

步骤(4-1)：基于步骤(3)中构建的二阶特征张量，随机抽取训练样本，得到 χ_i表示一个训练样本，为二阶特征张量，N是样本总数；

步骤(4-2)：计算特征张量的mode-t展开矩阵

t＝1时，mode-1展开矩阵：

${\chi }_i^{\left(1\right)}=X_i^{\left(1\right)}---\left(5\right)$

t＝2时，mode-2展开矩阵：

${\chi }_i^{\left(2\right)}=X_i^{\left(2\right)}=X_i^{\left(1\right)T}---\left(6\right)$

中上标表示张量χ_i的mode-t展开，是张量χ_i的mode-t展开矩阵，t＝1、2。 是张量χ_i的mode-1展开矩阵,x＝1,...,P,y＝1,...,15。是张量χ_i的mode-2展开矩阵，等于的转置。

步骤(4-3)：分别计算所有样本的特征张量的mode-t类内和类间散列矩阵：

$\left(\begin{array}{c}{S}_w^{\left(t\right)}=\underset{j=1}{\overset{N_c}{\Sigma }}\underset{i=1}{\overset{N_{c_j}}{\Sigma }}(X_i^{\left(t\right)}-{\bar{X}}_j^{\left(t\right)}){(X_i^{\left(t\right)}-{\bar{X}}_j^{\left(t\right)})}^T\\ S_B^{\left(t\right)}=\underset{j=1}{\overset{N_c}{\Sigma }}{N}_{c_j}({\bar{X}}_j^{\left(t\right)}-{\bar{X}}^{\left(t\right)}){({\bar{X}}_j^{\left(t\right)}-{\bar{X}}^{\left(t\right)})}^T\end{array}\right)---\left(7\right)$

表示所有样本的特征张量mode-t类内散列矩阵，表示所有样本的特征张量mode-t 类间散列矩阵，表示第j类内样本的特征张量mode-t展开矩阵的均值，表示所有样 本的特征张量mode-t展开矩阵的均值，参数t＝1、2，N_c表示类别总数，N_cj表示第j类中的 样本数。

步骤(4-4)：初始化稀疏映射矩阵l₁≤P,l₂≤15为任意列正交矩阵， 按式(8)迭代优化，得到最优的稀疏映射矩阵t＝1,2：

t＝1,2,μ是一个预设的常量，u_tj是U_t的第j列，||·||和|·|表示L₂范数和L₁范数，α_t是预设的 系数，β_tj是L₁范数的系数。

步骤(4-5)：张量特征的稀疏化映射，映射关系如下：

$T_{i,j}^s=T_{i,j}\times {}_{1}U_1^{*}\times {}_{2}U_2^{*}---\left(9\right)$

T_i,j是步骤(3)中构建的二阶特征张量，t＝1,2为步骤(4-4)优化后的稀疏映射矩 阵，是稀疏投影的结果，×_t,表示张量与矩阵的mode-t内积，t＝1,2。

所述步骤(4-5)中mode-t内积的方法为：

任意二阶特征张量以及步骤(4-4)的稀疏映射矩阵L₁≤P, L₂≤15，定义中间张量

张量mode-1内积：

$\left(\begin{array}{c}B=T_{i,j}\times {}_{1}U_1^{*}\\ B(x,y)=\underset{k=1}{\overset{P}{\Sigma }}{T}_{i,j}(k,y)*U_1^{*}(k,x),x=1,...,L_1,y=1,...,15\end{array}\right)---\left(10\right)$

张量mode-2内积：

$\left(\begin{array}{c}{T}_{i,j}^s=B\times {{{}_{2}U}^{*}}_2\\ T_{i,j}^s(x,y)=\underset{k=1}{\overset{15}{\Sigma }}B(x,k)U_2^{*}(k,y),x=1,...,L_1,y=1,...,L_2\end{array}\right)---\left(11\right)$

所述步骤(5)低维特征向量重新表示成向量形式方法为：

将按行或者列展开成一维向量形式。

本发明的有益效果：

1.本方法通过二层三阶离散小波变换，抽取良好的光谱细节变化特征；

2.本方法通过将不同层不同阶的小波变换系数表示为特征张量形式，保持了良好的特 征之间的结构相关性；

3.本方法采用稀疏张量判别分析方法，对特征张量进行稀疏映射，抽取充足的判别特 征。

4.本方法具有良好的稳定性，计算复杂度低，判别特征充分，分类精度高等优点。

附图说明

图1 本发明流程图；

图2 离散信号快速小波变换示意图；

图3 3-D离散小波变换示意图。

具体实施方式

下面结合附图与实施例对本发明作进一步说明。

如图1所示，一种基于3-D小波变换和稀疏张量判别分析的高光谱特征抽取方法的过程 是：

(1)数据归一化。给定高光谱图像数据立方体X、Y表示高光谱图像的空间 维度，P是波段数，表示空间坐标为(i,j)的样本(光谱向量)，采用如下归一 化方法：

$\left(\begin{array}{c}C(i,j,k)=\frac{C(i,j,k)}{{\sigma }_k},i=1,...,X,j=1,...,Y,k=1,...,P\\ {\mu }_k=\frac{1}{X*Y}\underset{i=1}{\overset{X}{\Sigma }}\underset{j=1}{\overset{Y}{\Sigma }}C(i,j,k)\\ {\sigma }_k^2=\frac{1}{X*Y}\underset{i=1}{\overset{X}{\Sigma }}\underset{j=1}{\overset{Y}{\Sigma }}{[C(i,j,k)-{\mu }_k]}^2\end{array}\right)$

C(i,j,k)是样本C(i,j,·)的第k个特征，是所有样本第k个特征的均值和方差， C(i,j,k)是归一化后的新特征值。

(2)3-D离散小波变换。将任意样本对应的光谱向量视为离散信号，考虑离散信号的快 速小波变换，如图2所示为离散信号的快速小波变换的过程，L和H分别代表低通和高通， 第j-1阶尺度系数c_j-1分别通过低通、高通滤波器得到第j阶尺度系数c_j和小波系数d_j。

本发明选择Haar二进小波作为快速小波变换的母波，则低通、高通滤波器函数如下：

$h_0\left(k\right)=(\frac{1}{\sqrt{2}},\frac{1}{\sqrt{2}})$

$h_1\left(k\right)=(-\frac{1}{\sqrt{2}},\frac{1}{\sqrt{2}})$

本发明设计两层三阶3-D离散小波变换，如图3所示，每层三阶3-D小波变换都是一个 完全生成树结构，前一阶的小波变换输出的尺度系数和小波系数作为下一阶小波变换的输入， 两层相同的结构级联，令初始尺度系数c₀＝C(i,j,·)，(i＝1,...,X j＝1,...,Y)，，依次对C_n中 的每个像素点对应的光谱向量进行小波变换，得到15个小波变换系数立方体C_k(k＝1,...,15)。

(3)二阶特征张量表示。如图3所示任意光谱向量c₀，经过3-D离散小波变换以后，得 到一组i＝1，...15，本发明采用二阶特征张量表示3-D离散小波变换后的特征，二阶 特征张量即为矩阵，则对应c₀可以得到二阶特征张量：P是波段数。

(4)稀疏张量判别分析。基于(3)中构建的二阶特征张量，按类别随机抽取30％的样 本作为训练样本，首先进行张量mode-展开，t＝1,2，然后分别求取mode-1和mode-2下的张 量类内散列矩阵和类间散列矩阵，预设参数μ∈[0,1]，以任意列正交矩阵初始化稀疏投影矩 阵U_t,t＝1,2，运用Elastic net方法求解公式(8)，得到最优的稀疏投影矩阵t＝1,2。

(5)特征张量稀疏投影映射，其基本发方法是求解步骤(3)中构建的二阶特征张量与 步骤(4)中求解的最优的稀疏投影矩阵t＝1,2的mode-t内积，得到二阶稀疏特征张量。

(6)一维向量表示。将步骤(5)中得到的二阶稀疏特征张量按行或者列展开成以为向 量。

上述虽然结合附图对本发明的具体实施方式进行了描述，但并非对本发明保护范围的限 制，所属领域技术人员应该明白，在本发明的技术方案的基础上，本领域技术人员不需要付 出创造性劳动即可做出的各种修改或变形仍在本发明的保护范围以内。