CN102521603A - 基于条件随机场的超光谱图像分类方法 - Google Patents

Abstract

本发明公开了一种基于条件随机场的超光谱图像分类方法，其实现步骤为：(1)初始化；(2)输入分别与待分类超光谱图像地物类别分别对应的训练样本集；(3)训练一元势团信息量的参数；(4)训练二元势团信息量的参数；(5)输入待分类的测试数据；(6)确定一元势团信息量的后验概率；(7)确定二元势团信息量的后验概率；(8)信息更新；(9)确定类标。本发明采用条件随机场模型避免了对观测数据的明确建模；既利用了空间邻域信息，又避免了空间邻域信息的复杂建模，具有计算复杂度低，分类准确率高，算法适应性广泛的优点。可应用于遥感图像的地物分类、目标检测及识别、军事侦察等领域。

Description

基于条件随机场的超光谱图像分类方法

技术领域

本发明属于图像处理技术领域，更进一步涉及目标识别领域中基于条件随机场的超光谱图像分类方法。该方法可应用于目标识别图像信息的获取、超光谱图像目标识别，能够对图像的不同区域进行准确分类。

背景技术

遥感图像处理领域的一项重要任务，就是地物分类。而地物分类方面的一个焦点就是超光谱图像的分类。超光谱图像可通过数百个连续且细分的光谱波段对地表区域同时成像，进而获得三维图像数据。高光谱图像信息通过较窄的波段区间，较多的波段数量而提供，包含丰富的光谱信息，且表现出非常强的谱相关性和空间相关性，这为特征提取、目标识别和分类提供了可能，使得地物分类具有更高的精度和鲁棒性。

目前，有很多超光谱数据分类的方法被提出。其中一个主要的研究方向就是将统计的方法应用于超光谱图像的分类。部分现有的理论中存在一些显著问题，如需要对观测到的超光谱图像像素值做明确的建模、仅考虑了像素的谱矢量信息而忽略其空间相关性。

北京航空航天大学在其专利申请“基于免疫进化策略的高光谱图像分类方法”(专利申请号：201010139249.2，公开号：CN101826160A)中提出了一种基于免疫进化策略的高光谱图像分类方法。该方法虽然能够自适应地选取不同场景下不同地物分类所需的最佳波段组合，但仍存在的不足是，忽略了高光谱图像的空间邻域信息，分类准确度不高。

发明内容

本发明的目的在于克服上述已有技术的不足，提出了一种基于基于条件随机场的超光谱图像分类方法。本发明简化了超光谱图像的分类过程，充分利用了超光谱图像的谱相关性及空间相关性，且避免了对观测数据的建模。

本发明实现上述目的的思路是：采用条件随机场模型，结合谱矢量信息(一元势团信息)和空间邻域信息(二元势团信息)，直接对分类的后验概率进行建模，通过计算后验概率来确定分类结果。

本发明实现的具体步骤如下：

(1)初始化，对原始数据所有像素值进行初始化，使其满足均值为0、方差为1；

(2)输入分别与待分类超光谱图像地物类别分别对应的训练样本集；

(3)采用多项式对数归一(MLR)法训练一元势团信息量的参数

3a)确定一元势团信息量的黑塞矩阵(Hessian matrix)边界；

3b)确定一元势团信息量的梯度矢量(Gradient vector)；

3c)确定一元势团信息量的参数；

(4)采用多项式对数归一(MLR)法训练二元势团信息量的参数

4a)确定二元势团信息量的黑塞矩阵(Hessian matrix)边界；

4b)确定二元势团信息量的梯度矢量(Gradient vector)；

4c)确定二元势团信息量的参数；

(5)输入待分类的测试数据；

(6)确定一元势团信息量的后验概率；

(7)确定二元势团信息量的后验概率；

(8)信息更新，求一元势团信息量后验概率矢量和二元势团信息量后验概率矢量的平均矢量，一元势团信息量的后验概率，表示为L维概率矢量，L为超光谱图像的类别数，二元势团信息量的后验概率，取两相邻像元同属于每一类的概率，表示为L维概率矢量；

(9)确定类标

9a)从后验概率平均矢量中选取最大值；

9b)将最大值对应的类标作为最终类标。

本发明与现有技术相比具有以下优点：

第一，由于本发明采用条件随机场模型，直接将后验概率建模为吉布斯(Gibbs)分布，克服了现有技术对观测到的超光谱数据的明确建模，使得本发明可用于各种具有不同统计特性的超光谱图像的分类，具有普遍适用性。

第二，由于本发明采用条件随机场模型，结合了空间邻域信息，克服了现有技术忽略了超光谱图像的空间邻域信息，使得本发明所分类的图像正确率较高。

第三，由于本发明采用条件随机场模型，无论是在标记样本还是观测数据中，空间邻域信息可通过条件随机场模型直接获取，克服了现有技术需对观测数据空间邻域信息做复杂的建模，使得本发明的计算复杂度低。

附图说明

图1是本发明的流程图；

具体实施方式

参照图1，本发明的具体实施步骤如下。

步骤1.初始化，按如下步骤操作：

第一步，按照下式求原始数据对应矩阵RHSI所有像素的平均值：

$k=\underset{i}{\mathrm{\Σ}}\underset{j}{\mathrm{\Σ}}\underset{h}{\mathrm{\Σ}}\left(\mathrm{RHSI}\right)/N$

其中，k为所有像素的平均值，∑表示叠加，i表示矩阵行数，j表示矩阵列数，h表示波段数，

为所有行叠加，

为所有列叠加，

为所有波段叠加，RHSI为原始数据对应的矩阵，N＝i＊j＊h表示像素点总个数；

第二步，对像素点总个数N开平方；

第三步，对每一个像素点，减去像素平均值k，再除以像素点总个数开平方所得的值，得到新的数据矩阵，此矩阵的值应当满足均值为0、方差为1。

步骤2.输入分别与所有需求的地物类别对应的待分类超光谱图像训练样本集。

步骤3.采用多项式对数归一(MLR)法训练一元势团信息量的参数

首先按照下式预确定一元势团信息量的黑塞矩阵(Hessian matrix)边界：

$B\≡-(1/2)[I-{11}^T/L]\⊗\mathrm{\Σ}{\tilde{y}}_i{\tilde{y}}_i^T$

其中，B为一元势团信息量的黑塞矩阵边界，符号≡表示定义为，I为L-1维单位矩阵，L为超光谱图像的类别数，1＝[1，1，...，1]^T是L-1维矢量，上标T表示转置，

表示Kronecker矩阵乘积，∑表示叠加，

为像元i对应的谱矢量；

其次按照下预式确定一元势团信息量的梯度矢量(Gradient vector)：

$g\left(w\right)=\mathrm{\Σ}({\tilde{x}}_i^{\′}-p_i\left(w\right))\⊗{\tilde{y}}_i$

其中，g(w)为一元势团信息量的梯度矢量，w为一元势团信息量的参数，∑表示叠加，

中第l维元素为1，l为像元i的类标，其余元素为0，共L-1维，L为超光谱图像的类别数，

表示像元i的类标为l的概率，

为像元i对应的谱矢量；

最后按照下式预确定一元势团信息量的参数：

w^(t+1)＝w^(t)-B^-1g(w^(t))

其中，w为一元势团信息量的参数，上标^(·)表示迭代次数，t即迭代次数，本发明实施例中，迭代次数为300次。B^-1表示一元势团信息量的黑塞矩阵边界B的逆矩阵，g(w)表示一元势团信息量的梯度矢量。

步骤4.采用多项式对数归一(MLR)法训练二元势团信息量的参数

首先按照下式确定二元势团信息量的黑塞矩阵(Hessian matrix)边界：

$B\≡-(1/2)[I-{11}^T/L]\⊗\mathrm{\Σ}{\mathrm{\μ}}_{\mathrm{ij}}{\mathrm{\μ}}_{\mathrm{ij}}^T$

其中，B为二元势团信息量的黑塞矩阵的边界，符号≡表示定义为，I为L维单位矩阵，L为超光谱图像的类别数，1＝[1，1，...，1]^T是L维矢量，上标^T表示转置，

表示Kronecker矩阵乘积，∑表示叠加，μ_ij为像元组(i，j)的谱特征矢量，本发明实施例中，μ_ij(y)为连接像元i与像元j的谱矢量所得，像元j为像元i的四邻域像元；

其次按照下式确定二元势团信息量的梯度矢量(Gradient vector)：

$g\left(v\right)=\mathrm{\Σ}({\tilde{x}}_{\mathrm{ij}}^{\′}-p_{\mathrm{ij}}\left(v\right))\⊗{\mathrm{\μ}}_{\mathrm{ij}}$

其中，g(v)为二元势团信息量的梯度矢量，v为二元势团信息量的参数，∑表示叠加，

中第l维元素为1，l为像元i与像元j共同的类标，其余元素为0，共L维，L为超光谱图像的类别数，

其中

为像元组(i，j)类标为(l_i，l_j)的概率，μ_ij为像元组(i，j)的谱特征矢量，本发明实施例中，μ_ij(y)为连接像元i与像元j的谱矢量所得，像元j为像元i的四邻域像元；

最后按照下式确定二元势团信息量的参数：

v^(t+1)＝v^(t)-B^-1g(v^(t))

其中，v为二元势团信息量的参数，上标^(·)表示迭代次数，t即迭代次数，本发明实施例中，迭代次数为300次，B^-1表示二元势团信息量的黑塞矩阵的边界B的逆矩阵，g(v)为二元势团信息量的梯度矢量。

步骤5.输入待分类的测试数据。

步骤6.确定一元势团信息量的后验概率，按照下式求得：

$P(x_i=l|y,w)=\left\{\begin{array}{cc}\frac{\exp \left(w_l^T{y}_i\right)}{1+{\mathrm{\&Sigma;}}_{k=1}^{L-1}\exp \left(w_k^T{y}_i\right)}& l<L\\ \frac{1}{1+{\mathrm{\&Sigma;}}_{k=1}^{L-1}\exp \left(w_k^T{y}_i\right)}& l=L\end{array}\right.$

其中，P(x_i＝l|y，w)表示像元i类标为l的后验概率，y为待分类谱矢量，w为一元势团信息量的参数，exp为指数函数，w_l为第l类对应的参数，l∈{1，2，...，L}，L为超光谱图像所需分为的类别数，上标^T表示转置，y_i为像元i对应的谱矢量，∑表示叠加，w_k为第k类对应的参数。

步骤7.确定二元势团信息量的后验概率，按照下式求得：

$P((x_i,x_j)\&equiv;k|{\mathrm{\&mu;}}_{\mathrm{ij}}\left(y\right),v)=\left\{\begin{array}{cc}\frac{\exp \left(v_{\mathrm{kk}}^T{\mathrm{\&mu;}}_{\mathrm{ij}}\left(y\right)\right)}{1+{\mathrm{\&Sigma;}}_{n=1}^L\exp \left(v_{\mathrm{nn}}^T{\mathrm{\&mu;}}_{\mathrm{ij}}\left(y\right)\right)}& k\&le;L\\ \frac{1}{1+{\mathrm{\&Sigma;}}_{n=1}^L\exp \left(v_{\mathrm{nn}}^T{\mathrm{\&mu;}}_{\mathrm{ij}}\left(y\right)\right)}& k=L+1\end{array}\right.$

其中，P((x_i，x_j)≡k|μ_ij(y)，v)表示像元组(i，j)类标为k的后验概率，k∈{1，2，...，L+1}，L为超光谱图像的类别数，L+1表示像元组(i，j)中两个像元类别不同，μ_ij(y)为像元组(i，j)的谱特征矢量，本发明实施例中，μ_ij(y)为连接像元i与像元j的谱矢量所得，像元j为像元i的四邻域像元，v为二元势团信息量的参数，exp为指数函数，v_kk为第k类对应参数，上标^T表示转置，∑表示叠加，v_nn为第n类对应参数。

步骤8.信息更新，一元势团信息量的后验概率，表示为L维概率矢量，L为超光谱图像的类别数，二元势团信息量的后验概率，取两相邻像元同属于每一类的概率，表示为L维概率矢量，求一元势团信息量后验概率矢量和二元势团信息量后验概率矢量的平均矢量。

步骤9.从后验概率平均矢量中选取最大值，将最大值对应的类标作为最终类标。

本发明的效果可通过以下仿真进一步说明。

1、仿真条件

本发明的仿真试验是在主频2.27GHZ的intel Core i5 CPU、内存3.06GB的硬件环境和MATLAB 7.0的软件环境下进行的。

2、仿真内容

本发明仿真试验的数据来源为KSC(Kennedy Space Center)数据，它是由NASAAVIRIS在1996年3月对KSC区域在大约20km高空采集到的，空间分辨率约为18m，AVIRIS的带宽范围为400-2500nm，间隔10nm，原始数据有220个波段。对KSC数据，去除低信噪比的波段之后共有176个波段，训练数据由KSC的红外成像和Landsat的TM成像得出。

下表为本发明仿真实验中选用的数据类别、数目及其所占比例。

	类别	数目(比例)
			1	灌木	761(14.6％)
2	柳木沼泽	243(4.66％)
			3	棕榈	256(4.92％)
4	棕榈/橡树	252(4.84％)
			5	松木	161(3.07％)
6	橡树/阔叶林	229(4.38％)
			7	硬木沼泽	105(2.00％)
8	禾草沼泽	431(8.27％)
			9	虎尾草沼泽	520(9.99％)
10	香蒲沼泽	404(7.76％)
			11	盐碱沼泽	419(8.04％)
12	泥泞地	503(9.66％)
			13	水域	927(17.8％)
	总数	5211

本发明仿真内容为，从上表的每一类数据中随机选取一部分做训练样本，构成训练样本集，剩余部分构成测试数据集，对测试数据集进行分类。

3、仿真效果分析

下表为本发明仿真中对测试数据进行分类，并通过与已知信息的比对得到的分类正确率如下：

正确率	15％	30％	45％	60％	75％
						1	0.80341	0.87351	0.94511	0.9377	0.98953
2	0.89372	0.9037	0.95556	0.9798	0.98361
						3	0.66055	0.7305	0.86525	0.86538	0.89063
4	0.21333	0.58741	0.46853	0.79412	0.79365
						5	0.67153	0.76404	0.86517	0.8125	0.70732
6	0.63077	0.49606	0.55906	0.58511	0.67241
						7	0.43333	0.50877	0.63158	0.7193	0.92593
8	0.66485	0.80172	0.86638	0.98225	1
						9	0.83484	0.87762	0.92308	0.99519	0.98462
10	0.86047	0.9375	0.87946	0.89634	0.9703
						11	0.85154	0.85652	0.86087	0.97605	1
12	0.93925	0.90647	0.92446	0.99015	0.98413
						13	0.97589	0.98235	0.96471	0.94879	1
Total	0.79523	0.84709	0.87321	0.91979	0.94721

上表中第一行表示随机选取的训练样本占数据总量(训练数据和测试数据)的百分比，第一列为类标，最后一行为总的正确率。从上表中可以看出，在训练样本百分比较小时，数据较多的类别能达到很好的分类效果，如第13类，部分数据较少的类别由于训练样本光谱特征的覆盖率较高，也能达到较高的分类正确率，如第2类；随着训练样本占数据总量百分比的增加，无论是单独类中的数据，还是总体数据，分类正确率均有提高，能达到较高的标准。总体来说，本发明所用方法对超光谱图像的分类起到了较好的效果。

Claims (10)

1.一种基于条件随机场的超光谱图像分类方法，包括如下步骤：

(1)初始化，对原始数据所有像素值进行初始化，使其满足均值为0、方差为1；

(2)输入分别与待分类超光谱图像地物类别分别对应的训练样本集；

(3)采用多项式对数归一(MLR)法训练一元势团信息量的参数

3a)确定一元势团信息量的黑塞矩阵(Hessian matrix)边界；

3b)确定一元势团信息量的梯度矢量(Gradient vector)；

3c)确定一元势团信息量的参数；

(4)采用多项式对数归一(MLR)法训练二元势团信息量的参数

4a)确定二元势团信息量的黑塞矩阵(Hessian matrix)边界；

4b)确定二元势团信息量的梯度矢量(Gradient vector)；

4c)确定二元势团信息量的参数；

(5)输入待分类的测试数据；

(6)确定一元势团信息量的后验概率；

(7)确定二元势团信息量的后验概率；

(8)信息更新，求一元势团信息量后验概率矢量和二元势团信息量后验概率矢量的平均矢量，一元势团信息量的后验概率，表示为L维概率矢量，L为超光谱图像的类别数，二元势团信息量的后验概率，取两相邻像元同属于每一类的概率，表示为L维概率矢量；

(9)确定类标

9a)从后验概率平均矢量中选取最大值；

9b)将最大值对应的类标作为最终类标。

2.根据权利要求1所述的基于条件随机场的超光谱图像分类方法，其特征在于：步骤(1)所述的初始化的步骤是：

第一步，按照下式求原始数据对应矩阵RHSI所有像素的平均值：

$k=\underset{i}{\mathrm{\&Sigma;}}\underset{j}{\mathrm{\&Sigma;}}\underset{h}{\mathrm{\&Sigma;}}\left(\mathrm{RHSI}\right)/N$

其中，k为所有像素的平均值，∑表示叠加，i表示矩阵行数，j表示矩阵列数，h表示波段数，为所有行叠加，

为所有列叠加，

为所有波段叠加，RHSI为原始数据对应的矩阵，N＝i＊j＊h表示像素点总个数；

第二步，对像素点总个数N开平方；

第三步，对每一个像素点，减去像素平均值k，再除以像素点总个数开平方所得的值，得到新的数据矩阵，此矩阵的值应当满足均值为0、方差为1。

3.根据权利要求1所述的基于条件随机场的超光谱图像分类方法，其特征在于：步骤3a)所述的一元势团信息量的黑塞矩阵(Hessian matrix)边界按照下式求得：

$B\&equiv;-(1/2)[I-{11}^T/L]\&CircleTimes;\mathrm{\&Sigma;}{\tilde{y}}_i{\tilde{y}}_i^T$

其中，B为一元势团信息量的黑塞矩阵边界，符号≡表示定义为，I为L-1维单位矩阵，L为超光谱图像的类别数，1＝[1，1，...，1]^T是L-1维矢量，上标^T表示转置，

表示Kronecker矩阵乘积，∑表示叠加，

为像元i对应的谱矢量。

4.根据权利要求1所述的基于条件随机场的超光谱图像分类方法，其特征在于：步骤3b)所述的一元势团信息量的梯度矢量(Gradient vector)按照下式求得：

$g\left(w\right)=\mathrm{\&Sigma;}({\tilde{x}}_i^{\&prime;}-p_i\left(w\right))\&CircleTimes;{\tilde{y}}_i$

其中，g(w)为一元势团信息量的梯度矢量，w为一元势团信息量的参数，∑表示叠加，

中第l维元素为1，l为像元i的类标，其余元素为0，共L-1维，L为超光谱图像的类别数，

表示像元i的类标为l的概率，

为像元i对应的谱矢量。

5.根据权利要求1所述的基于条件随机场的超光谱图像分类方法，其特征在于：步骤3c)所述的一元势团信息量的参数按照下式求得：

w^(t+1)＝w^(t)-B^-1g(w^(t))

其中，w为一元势团信息量的参数，上标^(·)表示迭代次数，t即迭代次数，B^-1表示一元势团信息量的黑塞矩阵边界B的逆矩阵，g(w)表示一元势团信息量的梯度矢量。

6.根据权利要求1所述的基于条件随机场的超光谱图像分类方法，其特征在于：步骤4a)所述的二元势团信息量的黑塞矩阵(Hessian matrix)边界按照下式求得：

$B\&equiv;-(1/2)[I-{11}^T/L]\&CircleTimes;\mathrm{\&Sigma;}{\mathrm{\&mu;}}_{\mathrm{ij}}{\mathrm{\&mu;}}_{\mathrm{ij}}^T$

其中，B为二元势团信息量的黑塞矩阵的边界，符号≡表示定义为，I为L维单位矩阵，L为超光谱图像的类别数，1＝[1，1，...，1]^T是L维矢量，上标^T表示转置，

表示Kronecker矩阵乘积，∑表示叠加，μ_ij为像元组(i，j)的特征矢量。

7.根据权利要求1所述的基于条件随机场的超光谱图像分类方法，其特征在于：步骤4b)所述的二元势团信息量的梯度矢量(Gradient vector)按照下式求得：

$g\left(v\right)=\mathrm{\&Sigma;}({\tilde{x}}_{\mathrm{ij}}^{\&prime;}-p_{\mathrm{ij}}\left(v\right))\&CircleTimes;{\mathrm{\&mu;}}_{\mathrm{ij}}$

其中，g(v)为二元势团信息量的梯度矢量，v为二元势团信息量的参数，∑表示叠加，

中第l维元素为1，l为像元i与像元j共同的类标，其余元素为0，共L维，L为超光谱图像的类别数，

其中

为像元组(i，j)类标为(l_i，l_j)的概率，μ_ij为像元组(i，j)的特征矢量。

8.根据权利要求1所述的基于条件随机场的超光谱图像分类方法，其特征在于：步骤4c)所述的二元势团信息量的参数按照下式求得：

v^(t+1)＝v^(t)-B^-1g(v^(t))

其中，v为二元势团信息量的参数，上标^(·)表示迭代次数，t即迭代次数，B^-1表示二元势团信息量的黑塞矩阵的边界B的逆矩阵，g(v)为二元势团信息量的梯度矢量。

9.根据权利要求1所述的基于条件随机场的超光谱图像分类方法，其特征在于：步骤(6)所述的一元势团后验概率按照下式求得：

$P(x_i=l|y,w)=\left\{\begin{array}{cc}\frac{\exp \left(w_l^T{y}_i\right)}{1+{\mathrm{\&Sigma;}}_{k=1}^{L-1}\exp \left(w_k^T{y}_i\right)}& l<L\\ \frac{1}{1+{\mathrm{\&Sigma;}}_{k=1}^{L-1}\exp \left(w_k^T{y}_i\right)}& l=L\end{array}\right.$

其中，P(x_i＝l|y，w)表示像元i类标为l的后验概率，y为待分类谱矢量，w为一元势团信息量的参数，exp为指数函数，w_l为第l类对应的参数，l∈{1，2，...，L}，L为超光谱图像所需分为的类别数，上标T表示转置，y_i为像元i对应的谱矢量，∑表示叠加，w_k为第k类对应的参数。

10.根据权利要求1所述的基于条件随机场的超光谱图像分类方法，其特征在于：步骤(7)所述的二元势团后验概率按照下式求得：

$P((x_i,x_j)\&equiv;k|{\mathrm{\&mu;}}_{\mathrm{ij}}\left(y\right),v)=\left\{\begin{array}{cc}\frac{\exp \left(v_{\mathrm{kk}}^T{\mathrm{\&mu;}}_{\mathrm{ij}}\left(y\right)\right)}{1+{\mathrm{\&Sigma;}}_{n=1}^L\exp \left(v_{\mathrm{nn}}^T{\mathrm{\&mu;}}_{\mathrm{ij}}\left(y\right)\right)}& k\&le;L\\ \frac{1}{1+{\mathrm{\&Sigma;}}_{n=1}^L\exp \left(v_{\mathrm{nn}}^T{\mathrm{\&mu;}}_{\mathrm{ij}}\left(y\right)\right)}& k=L+1\end{array}\right.$

其中，P((x_i，x_j)≡k|μ_ij(y)，v)表示像元组(i，j)类标为k的后验概率，k∈{1，2，...，L+1}，L为超光谱图像的类别数，L+1表示像元组(i，j)中两个像元类别不同，μ_ij(y)为连接像元i与像元j的谱矢量所得的特征矢量，像元j为像元i的四邻域像元，v为二元势团信息量的参数，exp为指数函数，v_kk为第k类对应参数，上标^T表示转置，∑表示叠加，v_nn为第n类对应参数。

Images