CN104504391A - 一种基于稀疏特征和马尔科夫随机场的高光谱图像分类方法 - Google Patents

Abstract

本发明属于遥感信息处理技术领域，具体涉及一种基于稀疏特征和马尔科夫随机场的高光谱图像分类方法。本发明包括：读入高光谱图像数据；求解字典；求解稀疏特征；用概率支持向量机求解概率输出，并确定初始分类结果：通过马尔科夫随机场来求解样本概率；确定高光谱图像最终分类结果。本发明应用稀疏特征使得图像得到很好的表述，能够对图像中的邻域信息充分的利用，优化了分类图的视觉效果，提高了分类的精度等优点。

Description

一种基于稀疏特征和马尔科夫随机场的高光谱图像分类方法

技术领域

本发明属于遥感信息处理技术领域，具体涉及一种基于稀疏特征和马尔科夫随机场的高光谱图像分类方法。

背景技术

高光谱图像分类是一种重要的获取信息的手段，以高光谱图像为对象，利用统计模式识别技术，提取研究取识别模式的统计特征值，然后根据某种决策准则做出类别决策，实现对不同地物的判别，其目标是将图像中的每个像元划分给一个类别。不同的地物由于反射的电磁能量不同表现出波谱的差异性，这就是高光谱分类的理论依据。良好的高光谱图像分类有助于充分挖掘高光谱图像信息，高光谱数据的特点是数据量大、冗余多，维数较高，同时在波段之间存在着很强的相关性。传统高光谱处理中通常用到的信息为原始的光谱(OriginalSpectral,OS)特征和主成分分析(Principal Components Analysis,PCA)特征信息，如何有效地利用丰富的光谱信息，并同时确保处理精度获得了越来越广泛的关注。

传统的高光谱数据分类方法中存在以下几个问题：1、高光谱数据不能得到很好的表达致使分类精度不高。2、没有充分利用邻域信息。针对以上问题，本发明提出一种稀疏特征和马尔科夫随机场的高光谱图像分类方法。

发明内容

本发明的目的是提供一种具有更优的分类效果，分类视觉效果好，分类精度高的基于稀疏特征和和马尔科夫随机场的高光谱图像分类方法。

本发明的目的是这样实现的：

(1)读入高光谱图像数据：

读入三维的高光谱高维数据，进行维数转换从三维转换为二维数据，对所得的二维数据作归一化处理得到高光谱遥感图像数据集X，确定要处理的样本类别数为s；

(2)求解字典D：

高光谱遥感图像数据集 表示n维实数集，其中n为高光谱数据的行数，x_p表示X中的第p个列向量，h为其列总数，字典d_z表示D中的第z个列向量，W为其列总数，把稀疏表示表达为：

$\underset{D,\mathrm{\α}}{\min }\frac{1}{2}{\left|\right|X-\mathrm{D\α}\left|\right|}_F^2+\mathrm{\λ}{\left|\right|\mathrm{\α}\left|\right|}_{\mathrm{1,1}}s.t.{\left|\right|d_z\left|\right|}_2\≤1\∀z\∈1,\·\·\·,W,---\left(1\right)$

其中，为系数矩阵，λ＞0，参数λ平衡重构误差和稀疏性之间的折中关系；得到高光谱数据所对应的字典D；

(3)求解稀疏特征A：

将像元x∈X表示为字典D中原子的稀疏线性组合：

$\widehat{\mathrm{\α}}\left(x\right)=\arg \underset{\mathrm{\α}}{\min }{\left|\right|x-\mathrm{D\α}\left|\right|}_2^2+\mathrm{\λ}{\left|\right|\mathrm{\α}\left|\right|}_1;$

则像元x在字典D上的稀疏表示特征即为 得到高光谱数据的稀疏特征A；

(4)用概率支持向量机求解概率输出，并确定初始分类结果：

(4.1)设定训练集和测试集：

从稀疏特征A中将s个类别各选取一样本来构成训练集E，整个稀疏特征A设定为测试集；

(4.2)进行样本标号：

将具有监督信息的训练样本中属于第i(1≤i≤s)类样本对应的样本标号标记为+1，其余s-1的类别对应的样本标号标记为-1；

(4.3)概率支持向量机二分类：

A中的元素e的分类决策函数为f_i(e)为：

f_i(e)＝<w,e>+b

其中，截距b为 ${\left[\begin{array}{cc}0& 1_v^T\\ 1_v& E^{T}E\end{array}\right]}^{-1}\left[\begin{array}{c}0\\ y\end{array}\right]$ 中的第一个元素，斜率w为式

w＝Eβ

其中，y为对应于E中第i类位置处的值为+1、其余元素值为-1的s维列向量，1v为元素值全为1的s维列向量，参量β为 ${\left[\begin{array}{cc}0& 1_v^T\\ 1_v& E^{T}E\end{array}\right]}^{-1}\left[\begin{array}{c}0\\ y\end{array}\right]$ 中去掉第一个元素的列向量；

(4.4)通过二分类结果为f_i(e)，求解出分类概率输出po(i|e)，并对结果进行存储，

$\mathrm{po}\left(i\right|e)\&ap;\frac{1}{1+\exp ({\mathrm{af}}_i\left(e\right)+B)},$

其中，a和B分别为斜率和截距参数；

(4.5)重复步骤(4.1)-(4.3)，直至分类概率输出po(s|e)求解完毕；

(4.6)对概率输出进行归一化，得到最终概率输出P(i|e)：

$P\left(i\right|e)=\frac{\mathrm{po}\left(i\right|e)}{\underset{i=1}{\overset{s}{\mathrm{\&Sigma;}}}\mathrm{po}\left(i\right|e)}$

(4.7)确定初始分类结果：

由步骤(4.1)-(4.6)得s个概率输出，稀疏特征A中测试样本的归属类别取决于s个结果中对应的最大值，样本就偏向于这个类别并且最终归属为这一类，得到初始分类结果Y₀；

(5)通过马尔科夫随机场来求解样本概率：

(5.1)求解能量函数U(e)：

初始分类结果Y₀为马尔科夫随机场，求解其能量总函数U(e)：

$U\left(e\right)=\underset{\left\{m\right\}\&Element;C_1}{\mathrm{\&Sigma;}}{V}_1\left(e_m\right)+\underset{\{m,n\}\&Element;C_2}{\mathrm{\&Sigma;}}{V}_2(e_m,e_n)$

其中C₁为一阶邻域集团，C₂为二阶邻域集团，e为目标像素，e_m与e_n为邻域集团中的像素，m和n表示的为坐标位置，与分别为一阶邻域集团和二阶邻域集团的能量函数，其中一阶能量V₁(e_m＝i)＝-ρ，二阶能量 $V_2(e_m,e_n)=\left\{\begin{array}{cc}-{\mathrm{\&beta;}}_c& (e_m=e_n)\\ {\mathrm{\&beta;}}_c& (e_m\&NotEqual;e_n)\end{array}\right.,$ ρ和β_c分别为一阶邻域势参数和二阶邻域势参数；

(5.2)求解样本e出现的概率P(e)：

求解样本e出现的概率P(e)：

$P\left(e\right)=Z^{-1}\&times;\exp (-\frac{1}{T}U\left(e\right))$

其中是切分函数的归一化常量，Ω表示样本空间，T为温度常数这里设为1；

(6)确定高光谱图像最终分类结果；

通过最大后验概率准则确定样本的最终标号BH，最终输出最终分类结果Y：

$\mathrm{BH}=\arg \underset{e}{\max }P\left(i\right|e)P\left(e\right).$

本发明的有益效果在于其应用稀疏特征使得图像得到很好的表述，能够对图像中的邻域信息充分的利用，优化了分类图的视觉效果，提高了分类的精度等优点。

附图说明

图1为本发明的流程图；

图2为实验中的Indian Pines高光谱数据的真实图像；

图3为实验中在Indian Pines高光谱数据中所选取的原始各类的地物分布图；

图4为实验中在Indian Pines高光谱数据中所选取的样本名称和样本数量表；

图5为实验中采用的三种不同方法对数据集分类后与之对应的各类地物的分类分布图，其中5-a，5-b，5-c分别对应着数据集的OS+SVM分类分布图，PCA+SVM分类分布图和本发明方法分类分布图；

图6为实验中采用的三种不同方法对数据集分类的评价指标表。

具体实施方式

下面结合附图对本发明做出更详细的描述。

本发明提供一种基于稀疏特征和马尔科夫随机场的高光谱图像分类方法。具体包括以下步骤：1、读入高光谱图像数据。读入高光谱高维数据，进行维数转换，并对其作归一化处理得到X，其中含有样本类别数为s。2、求解字典D。对X进行字典学习，得到高光谱数据所对应的字典D。3、求解稀疏特征A。利用X和已经求得的字典D求解高光谱数据的稀疏特征A。4、用概率支持向量机(Support Vector Machines,SVM)求解概率输出，并确定初始分类结果。从A中将s个类别各选取一部分样本来构成训练集E，测试集设定为整个A。通过概率支持向量机来求解分类概率输出，并确定初始分类结果Y₀。5、通过马尔科夫随机场来求解样本概率。根据初始分类结果Y₀，通过马尔科夫随机场求解样本e出现的概率P(e)。6、确定高光谱图像最终分类结果。通过最大后验概率准则来确定最终分类结果Y。

本发明的详细步骤如下：

1、读入高光谱图像数据。

读入三维的高光谱高维数据，对其进行维数转换使其从三维转换为二维数据以方便后续的处理，并对所得的二维数据作归一化处理得到X，确定要处理的样本类别数为s。

2、求解字典D。

高光谱遥感图像数据集 表示n维实数集，其中n为高光谱数据的行数，x_p表示X中的第p个列向量，h为其列总数。字典(d_z表示D中的第z个列向量，W为其列总数)可以把稀疏表示表达为如下形式的优化问题:

$\underset{D,\mathrm{\&alpha;}}{\min }\frac{1}{2}{\left|\right|X-\mathrm{D\&alpha;}\left|\right|}_F^2+\mathrm{\&lambda;}{\left|\right|\mathrm{\&alpha;}\left|\right|}_{\mathrm{1,1}}s.t.{\left|\right|d_z\left|\right|}_2\&le;1\&ForAll;z\&Element;1,\&CenterDot;\&CenterDot;\&CenterDot;,W,---\left(1\right)$

其中，为系数矩阵，λ＞0，参数λ平衡重构误差和稀疏性之间的折中关系。

通过对式(1)优化求解，最终得到高光谱数据所对应的字典D。

3、求解稀疏特征A。

可以将像元x∈X表示为字典D中原子的稀疏线性组合：

$\widehat{\mathrm{\&alpha;}}\left(x\right)=\arg \underset{\mathrm{\&alpha;}}{\min }{\left|\right|x-\mathrm{D\&alpha;}\left|\right|}_2^2+\mathrm{\&lambda;}{\left|\right|\mathrm{\&alpha;}\left|\right|}_1---\left(2\right)$

则某像元x在字典D上的稀疏表示特征即为 通过对式(2)求解可以得到高光谱数据的稀疏特征A。

4、用概率SVM求解概率输出，并确定初始分类结果。

4.1、设定训练集和测试集。

从A中将s个类别各选取一部分样本来构成训练集E，整个A设定为测试集。

4.2、进行样本标号。

将具有监督信息的训练样本中属于第i(1≤i≤s)类样本对应的样本标号标记为+1，其余s-1的类别对应的样本标号标记为-1。

4.3、SVM二分类。

A中的元素e的分类决策函数为f_i(e)为式(3):

f_i(e)＝<w,e>+b    (3)

其中，截距b为 ${\left[\begin{array}{cc}0& 1_v^T\\ 1_v& E^{T}E\end{array}\right]}^{-1}\left[\begin{array}{c}0\\ y\end{array}\right]$ 中的第一个元素，斜率w为式(4):

w＝Eβ   (4)

其中，y为对应于E中第i类位置处的值为+1、其余元素值为-1的s维列向量，1v为元素值全为1的s维列向量，参量β为 ${\left[\begin{array}{cc}0& 1_v^T\\ 1_v& E^{T}E\end{array}\right]}^{-1}\left[\begin{array}{c}0\\ y\end{array}\right]$ 中去掉第一个元素的列向量。

4.4、通过5.2得到二分类结果为f_i(e)，通过式(5)求解出分类概率输出po(i|e)，并对结果进行存储。

$\mathrm{po}\left(i\right|e)\&ap;\frac{1}{1+\exp ({\mathrm{af}}_i\left(e\right)+B)}---\left(5\right)$

其中，a和B分别为斜率和截距参数。

4.5、重复5.1-5.3，直至po(s|e)求解完毕。

4.6、通过式子(6)对概率输出进行归一化，得到最终概率输出P(i|e)。

$P\left(i\right|e)=\frac{\mathrm{po}\left(i\right|e)}{\underset{i=1}{\overset{s}{\mathrm{\&Sigma;}}}\mathrm{po}\left(i\right|e)}---\left(6\right)$

4.7、确定初始分类结果。

由步骤4.1-4.6可得s个概率输出，A中测试样本的归属类别取决于s个结果中对应的哪个值最大，样本就偏向于这个类别并且最终归属为这一类，得到初始分类结果Y₀。

5、通过马尔科夫随机场来求解样本概率。

5.1、求解能量函数U(e)。

初始分类结果Y₀为马尔科夫随机场，通过式子(7)来求解其能量总函数U(e)。

$U\left(e\right)=\underset{\left\{m\right\}\&Element;C_1}{\mathrm{\&Sigma;}}{V}_1\left(e_m\right)+\underset{\{m,n\}\&Element;C_2}{\mathrm{\&Sigma;}}{V}_2(e_m,e_n)---\left(7\right)$

其中C₁为一阶邻域集团，C₂为二阶邻域集团，e为目标像素，e_m与e_n为其邻域集团中的像素，m和n表示的为坐标位置。与分别为一阶邻域集团和二阶邻域集团的能量函数，其中一阶能量V₁(e_m＝i)＝-ρ，二阶能量 $V_2(e_m,e_n)=\left\{\begin{array}{cc}-{\mathrm{\&beta;}}_c& (e_m=e_n)\\ {\mathrm{\&beta;}}_c& (e_m\&NotEqual;e_n)\end{array}\right.,$ ρ和β_c分别为一阶邻域势参数和二阶邻域势参数。

5.2、求解样本e出现的概率P(e)。

通过式(7)来求解样本e出现的概率P(e)。

$P\left(e\right)=Z^{-1}\&times;\exp (-\frac{1}{T}U\left(e\right))---\left(7\right)$

其中是切分函数的归一化常量，Ω表示样本空间，T为温度常数这里设为1。

6、确定高光谱图像最终分类结果。

通过最大后验概率准则即式(8)来确定样本的最终标号BH，最终输出最终分类结果Y:

$\mathrm{BH}=\arg \underset{e}{\max }P\left(i\right|e)P\left(e\right)---\left(8\right)$

通过步骤1至步骤6，最终输出本发明下的高光谱图像分类结果Y。图1给出了本发明方法的流程图。

为了说明本发明的有效性，特进行如下实验论证。

其中采用Indian Pines高光谱数据集来验证本发明方法的适用性。

美国印第安纳州Indian Pine实验区图像，其于1992年6月用AVIRIS传感器采集得到，空间分辨率为20m。原始图像共有220个波段，大小为144×144，共有16种地物分布，将原始的220个波段中受噪声影响较大的一些波段去除后选取200个波段作为仿真研究对象。考虑到地物个数、比例等原因从中16种地物中选取9类地物来进行实验，为方便对其分别标号为第1-9类。实验数据如图2所示，选取的原始各类的地物分布如图3所示，所选取的样本名称和样本数量如图4所示。从这9种实际地物分布样本中均匀抽取10％的数据作为训练样本。

在对高光谱图像进行分类时，本发明方法与经典方法即分别利用OS特征和PCA特征然后经SVM分类方法(两种方法分别标记为OS+SVM和PCA+SVM)进行对比。

运用三种方法分类后与之对应的各类地物的分类分布图如图5所示，可以很直观的看到相对于OS+SVM和PCA+SVM分类，本发明的方法分类效果图好。其中5-a，5-b，5-c分别对应着数据集的OS+SVM分类分布图，PCA+SVM分类分布图和本发明方法分类分布图；

三种分类方法的四个分类评价指标即总体分类精度，Kappa系数、平均漏分率和平均错分率如图6所示，其中由这些指标的定义可以知道总体分类精度和Kappa系数越高且平均漏分率和平均错分率越低时，图像的分类效果就越好。较之OS+SVM，针对本发明的方法中从总体分类精度和Kappa系数来看均要高出11％，从平均漏分率和平均错分率来均降低10％。较之PCA+SVM，针对本发明的方法中从总体分类精度和Kappa系数来看均要高出25％，从平均漏分率和平均错分率来均降低25％。

通过实验的对比分析可以进一步看出本发明方法优势所在。

Claims (1)

1.一种基于稀疏特征和马尔科夫随机场的高光谱图像分类方法，其特征在于：

(1)读入高光谱图像数据：

读入三维的高光谱高维数据，进行维数转换从三维转换为二维数据，对所得的二维数据作归一化处理得到高光谱遥感图像数据集X，确定要处理的样本类别数为s；

(2)求解字典D：

高光谱遥感图像数据集 表示n维实数集，其中n为高光谱数据的行数，x_p表示X中的第p个列向量，h为其列总数，字典d_z表示D中的第z个列向量，W为其列总数，把稀疏表示表达为：

$\underset{D,\mathrm{\&alpha;}}{\min }\frac{1}{2}{\left|\right|X-\mathrm{D\&alpha;}\left|\right|}_F^2+\mathrm{\&lambda;}{\left|\right|\mathrm{\&alpha;}\left|\right|}_{\mathrm{1,1}},s.t.{\left|\right|d_z\left|\right|}_2\&le;1,\&ForAll;z\&Element;1,\&CenterDot;\&CenterDot;\&CenterDot;,W,$

其中，为系数矩阵，λ＞0，参数λ平衡重构误差和稀疏性之间的折中关系；得到高光谱数据所对应的字典D；

(3)求解稀疏特征A：

将像元x∈X表示为字典D中原子的稀疏线性组合：

$\widehat{\mathrm{\&alpha;}}\left(x\right)=\arg \underset{\mathrm{\&alpha;}}{\min }{\left|\right|x-\mathrm{D\&alpha;}\left|\right|}_2^2+\mathrm{\&lambda;}{\left|\right|\mathrm{\&alpha;}\left|\right|}_1;$

则像元x在字典D上的稀疏表示特征即为 得到高光谱数据的稀疏特征A；

(4)用概率支持向量机求解概率输出，并确定初始分类结果：

(4.1)设定训练集和测试集：

从稀疏特征A中将s个类别各选取一样本来构成训练集E，整个稀疏特征A设定为测试集；

(4.2)进行样本标号：

将具有监督信息的训练样本中属于第i(1≤i≤s)类样本对应的样本标号标记为+1，其余s-1的类别对应的样本标号标记为-1；

(4.3)概率支持向量机二分类：

A中的元素e的分类决策函数为f_i(e)为：

f_i(e)＝<w,e>+b

其中，截距b为 ${\left[\begin{array}{cc}0& 1_v^T\\ 1_v& E^{T}E\end{array}\right]}^{-1}\left[\begin{array}{c}0\\ y\end{array}\right]$ 中的第一个元素，斜率w为式

w＝Eβ

其中，y为对应于E中第i类位置处的值为+1、其余元素值为-1的s维列向量，1_v为元素值全为1的s维列向量，参量β为 ${\left[\begin{array}{cc}0& 1_v^T\\ 1_v& E^{T}E\end{array}\right]}^{-1}\left[\begin{array}{c}0\\ y\end{array}\right]$ 中去掉第一个元素的列向量；

(4.4)通过二分类结果为f_i(e)，求解出分类概率输出po(i|e)，并对结果进行存储，

$\mathrm{po}\left(i\right|e)\&ap;\frac{1}{1+\exp ({\mathrm{af}}_i\left(e\right)+B)},$

其中，a和B分别为斜率和截距参数；

(4.5)重复步骤(4.1)-(4.3)，直至分类概率输出po(s|e)求解完毕；

(4.6)对概率输出进行归一化，得到最终概率输出P(i|e)：

$P\left(i\right|e)=\frac{\mathrm{po}\left(i\right|e)}{\underset{i=1}{\overset{s}{\mathrm{\&Sigma;}}}\mathrm{po}\left(i\right|e)}$

(4.7)确定初始分类结果：

由步骤(4.1)-(4.6)得s个概率输出，稀疏特征A中测试样本的归属类别取决于s个结果中对应的最大值，样本就偏向于这个类别并且最终归属为这一类，得到初始分类结果Y₀；

(5)通过马尔科夫随机场来求解样本概率：

(5.1)求解能量函数U(e)：

初始分类结果Y₀为马尔科夫随机场，求解其能量总函数U(e)：

$U\left(e\right)=\underset{\left\{m\right\}\&Element;C_1}{\mathrm{\&Sigma;}}{V}_1\left(e_m\right)+\underset{\{m,n\}\&Element;C_2}{\mathrm{\&Sigma;}}{V}_2(e_m,e_n)$

其中C₁为一阶邻域集团，C₂为二阶邻域集团，e为目标像素，e_m与e_n为邻域集团中的像素，m和n表示的为坐标位置，与分别为一阶邻域集团和二阶邻域集团的能量函数，其中一阶能量V₁(e_m＝i)＝-ρ，二阶能量 $V_2(e_m,e_n)=\left\{\begin{array}{cc}-{\mathrm{\&beta;}}_c& \left(e_m{=e}_n\right)\\ {\mathrm{\&beta;}}_c& (e_m\&NotEqual;e_n)\end{array}\right.,$ ρ和β_c分别为一阶邻域势参数和二阶邻域势参数；

(5.2)求解样本e出现的概率P(e)：

求解样本e出现的概率P(e)：

$P\left(e\right)=Z^{-1}\&times;\exp (-\frac{1}{T}U\left(e\right))$

其中是切分函数的归一化常量，Ω表示样本空间，T为温度常数这里设为1；

(6)确定高光谱图像最终分类结果；

通过最大后验概率准则确定样本的最终标号BH，最终输出最终分类结果Y：

$\mathrm{BH}=\arg \underset{e}{\max }P\left(i\right|e)P\left(e\right).$