CN103440505A - 空间邻域信息加权的高光谱遥感图像分类方法 - Google Patents

Abstract

本发明针对现有高光谱遥感图像分类技术，在图像分类结果存在“杂点”的问题上，提出了一种空间邻域信息加权的高光谱遥感图像分类方案。首先对预处理后的图像数据进行部分主成分分析，获取变换后的特征矩阵，然后对特征矩阵利用空间邻域信息加权的FCM算法实施快速聚类，得到高光谱图像的分类结果。本发明提供的高光谱遥感图像分类方案，充分结合了主成分分析特征降维的功效和高光谱图像丰富的空间邻域信息，在保证算法效率的同时，提高了高光谱图像的分类效果，与传统方法法相比，不仅能够降低计算量，而且可以改善分类效果，克服“同物异谱”和噪声引起的“杂点”现象，获得很好的分类效果。

Description

空间邻域信息加权的高光谱遥感图像分类方法

技术邻域

本发明属于高光谱遥感图像处理领域，具体涉及一种空间邻域信息加权的高光谱遥感图像分类方法。

背景技术

高光谱遥感技术是一种新兴遥感技术，它借助成像光谱仪将成像技术和细分光谱技术革命性地结合在一起，与传统多光谱遥感相比，高光谱遥感具有光谱分辨率高、图谱合一、在某一光谱段范围内连续成像的优势。这些特点和优势使得高光谱图像在军事目标的探测、植被的精细分类、地质岩矿物的识别、海洋检测、环境检测、以及城市规划等众多领域得到广泛应用。高光谱图像数据是一个光谱图像的立方体，通常有以下三种表示方法：

1)图像空间：它是一种最直观的高光谱数据信息的表现方式，反映了地物的分布和变化以及不同地物之间的空间关系。图像空间最重要的用途就是将图像中每一个像素与其地面位置对应起来，为高光谱影像处理与分析提供空间知识，但是波段之间的相互关系很难从图像中反映出来。

2)光谱空间：利用目标的光谱响应与波长之间的变化关系来描述高光谱数据所蕴含的信息。图像中每一个像元矢量对应一条近似连续的光谱曲线，它是电磁波能量对波长的函数，反映了其所代表的目标的辐射光谱信息。曲线中每个点的数值是相应成像波段上传感器对目标光谱辐射与反射能量的响应值。

3)特征空间：高光谱影像中的每一个像元对应着多个成像波段的反射值，这些反射值可以用一个多维矢量表示。在特征空间中，不同的目标分布在不同区域，并且有不同的分布特性，这就有利于定量地描述目标的光谱辐射特性及其在特征空间内的变化规律。

高光谱遥感技术中一个经典、具有探索性的研究领域就是遥感图像的分类。由于目前尚无通用的分类理论，而现有的大多数算法只针对具体问题，所以人们至今仍在不断研究和发掘新的分类理论与分类算法。高光谱图像分类方法，可以归结为两类：一类是基于光谱空间的分类方法，利用反映地物光谱特征的光谱曲线来识别地物；另一类是基于特征空间的分类方法，利用地物在特征空间的统计特性来建立分类模型。这两类方法通常只注重地物在光谱空间或特征空间上的差异，而忽视了像元“同物异谱”的现象和像元在图像空间维的结构，这种处理方式造成的后果就是图像分类后会有很多“杂点”，这是由噪声和地物分布在数据获取尺度上的不均匀引起的。一般情况下，地物在空间上都是连续分布的，或者说如果某点处为A类地物，那么在此点的邻域上是A类地物的可能性最大，这是一个非常重要的信息，利用好这个信息是消除“杂点”现象、提高分类精度的关键所在。

高光谱图像具有数据量大、冗余度高和波段间相关性强的特点，直接对原始数据处理并不经济，所以通常在进行分类之前会对遥感数据做特征降维处理，常用的高光谱数据降维方法是主成分分析(Principal Component Analysis)算法。主成分分析是考察多个变量间相关性的一种多元统计方法，研究如何通过少数几个主成分来揭示多个变量间的内部结构，即从原始变量中导出少数几个主成分，使它们尽可能多地保留原始变量的信息，且彼此间互不相关。这对高光谱数据有非常重要的意义：一方面，能够使图像远离噪声，提高数据质量；另一方面，去除图像中的无价值波段，可以降低计算量，提高处理效率。

发明内容

本发明所要解决的技术问题是，针对高光谱遥感图像数据量大、冗余高和现有分类技术在分类结果中存在“杂点”的问题上，提出了一种空间邻域信息加权的高光谱遥感图像分类方法。

本发明解决上述技术问题的技术方案是：一种空间邻域信息加权的高光谱遥感图像分类方法，包括：预处理阶段和分类阶段，预处理阶段，采用PPCA变换对高光谱图像波段进行降维，降低图像波段信息的冗余；分类阶段，采用图像空间邻域信息加权的模糊C均值算法，消除遥感图像分类结果中存在的“杂点”。具体技术方案为：

读取遥感图像，对图像特征矢量进行PPCA变换实现对特征矢量矩阵降维，从波段数据中获得部分线性不相关的矢量作为像元样本；从图像矩阵中提取像元样本的邻域像元，获取特征隶属度，根据特征隶属度计算像元样本所有邻域像元对应第某类特征隶属度的中值，获得特征矩阵中所有样本对应每个类别的空间邻域信息，根据空间邻域信息更新模糊划分矩阵和聚类中心，由模糊划分矩阵得到各像元分类结果，构造伪彩色图。

其中，对特征矢量矩阵降维包括：选取样本并转换成数据矩阵，计算样本数据的相关系数矩阵，计算相关系数矩阵的特征根，计算特征根的累计贡献率，根据累计贡献率选取特征维数，获取特征矢量。

具体为：提取图像中心像元选取样本，并转换成2维数据矩阵B，获得数据矩阵S的相关系数矩阵R；根据特征方程|R-λI|＝0，计算R的特征根λ_i，并使其从大到小排列：λ₁≥λ₂≥…≥λ_ω，得到对应的特征向量a＝{a₁，a₂，…a_ω}；根据公式：

计算累计贡献率，找到使累计贡献率第一次到达累计贡献率上限值K的特征维数L，提取前L个主成分对应的特征向量β＝{a₁，a₂，…a_L}，根据公式获取图像的特征矩阵X＝β^T·B，其中的特征矢量作为样本。

根据特征隶属度计算像元样本所有邻域像元对应第某类特征隶属度的中值具体为：从图像矩阵中提取样本x_j的邻域像元；获得x_j隶属于第i类的特征隶属度和空间邻域信息

获得特征隶属度

和空间邻域信息

加权后的模糊相似度 $u_{\mathrm{ij}}=(1-\mathrm{\α})\×v_{\mathrm{ij}}^1+\mathrm{\α}\×v_{\mathrm{ij}}^2,$ 以及聚类中心 $z_i={\mathrm{\Σ}}_{j=1}^p{\left(u_{\mathrm{ij}}\right)}^m{x}_j/{\mathrm{\Σ}}_{j=1}^p{\left(u_{\mathrm{ij}}\right)}^m,(i=\mathrm{1,2}...N),$ 更新模糊划分矩阵U＝{u_ij}和聚类中心Z＝{z_i}，其中，根据公式： $v_{\mathrm{ij}}^2=\frac{1}{2}(\max \left(\underset{m=1}{\overset{C}{{\mathrm{SNU}}_i}}\left(x_j\right)\right)+\min \left(\underset{m=1}{\overset{C}{{\mathrm{SNU}}_i}}\left(x_j\right)\right))$ 计算x_j对应第i类的空间邻域信息，

和

为样本x_j所有邻域像元对应第i类特征隶属度的最大值和最小值。

本发明对高光谱图像波段进行降维，降低图像波段信息的冗余，有效地提高了算法的效率；采用图像空间邻域信息加权的模糊C均值算法，消除了遥感图像分类结果中存在的“杂点”。

附图说明

图1九宫格选取中心法过程的示意图；

图2是本发明方案流程图。

具体实施方式

基于特征空间的分类方法是近年的热门研究方向，它可分为非监督分类和监督分类方法。模糊C均值聚类(FCM)是由Bezdek提出的一种目标函数法，是经典的非监督分类算法，也是模糊聚类算法中最著名、运用最广的算法。FCM算法是一种柔性的模糊划分，它能在一定程度上克服普通C均值算法硬划分的缺点。基本思想就是使得被划分到同一簇的对象之间相似度最大，而不同簇之间的相似度最小。

本文中主要涉及到的变量意义如下：IM表示3维高光谱图像数据矩阵；B表示与IM对应的2维数据矩阵；S表示部分主成分分析的样本矩阵；p代表部分主成分分析样本的个数；ω代表波段数；R表示S的相关系数矩阵；I为单位矩阵；λ表示R的特征根；K表示累计贡献率的上限值；L表示累计贡献率第一次到达K值的特征维数；X表示部分主成分分析变换后的特征矩阵；N表示聚类数目；C表示邻域窗口大小；m表示模糊指数；ε表示目标函数精度；t表示当前迭代的次数；T_m表示最大迭代次数；V¹表示特征隶属度；V²表示空间邻域信息；U表示模糊划分矩阵；Z表示聚类中心。

以下使用具体事例和附图对本发明作进一步说明。如图2为本发明方法流程示意图，本发明提出空间邻域信息加权的高光谱遥感图像分类方法，主要包括：预处理阶段，采用部分主成分分析(Partly Principal Component Analysis，简称PPCA变换)方法对高光谱图像波段进行降维，降低图像波段信息的冗余，有效地提高了算法的效率；分类阶段，采用图像空间邻域信息加权的模糊C均值算法(Fuzzy C-Means algorithm based on Space Neighbor，简称SNFCM算法)，消除遥感图像分类结果中存在的“杂点”。具体为：

预处理阶段：读取遥感图像，截取感兴趣区域；选取样本并转换成数据矩阵；为了排除数量级和量纲不同带来的影响，对样本数据进行标准化处理；计算标准化后样本数据的相关系数矩阵；计算相关系数矩阵的特征根，并使其从大到小排列，同时得到对应的特征向量；计算贡献率和累计贡献率，根据累计贡献率选取特征维数；计算原图像数据的主成分，获取遥感图像经PPCA变换后的特征矢量。

图像分类阶段：首先设置图像分类类别数、模糊指数、目标函数精度和最大迭代次数，并随机产生初始特征隶属度矩阵，计算初始聚类中心；然后根据特征隶属度矩阵计算空间邻域信息矩阵，加权空间信息矩阵和特征隶属度矩阵，得到隶属度函数；更新隶属度矩阵和聚类中心，直到前后两次目标函数差值小于目标函数精度或者迭代次数达到最大设置值；最后统计分类结果，给每一个图像像元作类别标示。

以下具体描述本发明的实现方式：

获取遥感图像图像特征矢量。读取高光谱图像，截取感兴趣区域图像，形成3维图像数据矩阵，根据图像尺寸和波段数，将3维矩阵通过变维变换成2维矩阵(每个波段数作为第1维数据，波段的图像数据作为第2维数据)，以便后续操作。

本实施例方案采用美国喷气实验室提供的高光谱图像数据源，这些数据源被绝大多数研究者使用，研究结果具有可比性。以常见的AVIRIS高光谱图像作为示例，它含有224个波段，每个波段的图像尺寸是512×614。读取高光谱图像，在空间维(每个波段作为一个空间维)上，截取大小为256×256的图像块，形成数据矩阵IM(256,256,224)；将矩阵IM的每个空间维逐列重排，构成一个行向量，最后一共形成224个行向量，这样3维矩阵IM就转化为2维矩阵B(224,256×256)。

1.对图像特征矢量进行PPCA变换。对图像特征矢量进行PPCA变换就是对特征矢量矩阵作降维处理，用少数的几个线性不相关的矢量来代替原有波段数据参与后续运算。本发明采用部分主成分分析(PPCA变换)方法，只利用高光谱图像的部分像元进行主成分分析，更好地减少了计算量，具体可采用如下方法实现：

1)在空间维中对图像矩阵IM采用九宫格取中心法提取样本，并转换成2维矩阵。如附图1所示，每间隔一个像元提取一个像元做样本(如图中黑色块代表)，提取顺序可自左向右、自上向下。

2)提取到的样本矩阵为：

$S=\left(\begin{array}{cccc}{s}_{11}& s_{12}& ...& s_{1p}\\ s_{21}& s_{22}& ...& s_{2p}\\ .& .& & .\\ .& .& ...& .\\ .& .& & .\\ s_{\mathrm{\ω}1}& s_{\mathrm{\ω}2}& ...& s_{\mathrm{\ωp}}\end{array}\right)\stackrel{\mathrm{\Δ}}{=}(s_1{,s}_2,...s_p)---\left(1\right)$

其中，p代表样本个数，ω代表波段数。为了排除数量级和量纲不同带来的影响，首先对样本数据根据公式：

进行标准化处理，获得标准化数据表(s′_ij)_ω×p，式中，s_ij为第j个样本第i个波段的数据；

和σ_i分别为样本第i个波段的均值和标准差。为简单起见，将标准化数据表(s′_ij)_ω×p仍记为(s_ij)_ω×p，

3)计算样本矩阵的相关系数矩阵R=(r_ij)_ω×p，其中：(2)

4)根据特征方程|R-λI|=0，计算R的特征根λ_i，并使其从大到小排列：λ₁≥λ₂≥…≥λ_ω，同时得到对应的特征向量a₁，a₂，…a_ω。其中，

$a_1=\left(\begin{array}{c}{a}_{11}\\ a_{21}\\ .\\ .\\ .\\ a_{\mathrm{\ω}1}\end{array}\right),a_2=\left(\begin{array}{c}{a}_{12}\\ a_{22}\\ .\\ .\\ .\\ a_{\mathrm{\ω}2}\end{array}\right),...,a_{\mathrm{\ω}}=\left(\begin{array}{c}{a}_{1\mathrm{\ω}}\\ a_{2\mathrm{\ω}}\\ .\\ .\\ .\\ a_{\mathrm{\ω\ω}}\end{array}\right)---\left(3\right)$

5)设定累计贡献率的上限值K，计算特征根的贡献率和累计贡献率：贡献率 $e_i={\mathrm{\λ}}_i/\underset{k=1}{\overset{\mathrm{\ω}}{\mathrm{\Σ}}}{\mathrm{\λ}}_k,$ 累计贡献率 $E_m=\underset{k=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\λ}}_m/\underset{k=1}{\overset{\mathrm{\ω}}{\mathrm{\Σ}}}{\mathrm{\λ}}_k(m=\mathrm{1,2},...,\mathrm{\ω}),$ 找到使累计贡献率第一次到达K值的特征维数L(即m值，代表选取样本的主成分个数)。

6)利用前L个主成分计算原图像数据B的特征矩阵X={x_ij}，其中：

x_ij=a_1iB_1j+a_2iB_2j+…+a_ωiB_ωj(i=1，2，…，L，j=1，2，…，p)  (4)

2.利用SNFCM算法对特征矩阵进行聚类。

由于传统FCM遥感图像聚类算法只考虑图像光谱特性，忽视了空间邻域信息，因此，本方案提出适用于高光谱图像的空间邻域信息加权FCM算法(简称SNFCM算法)，从图像矩阵IM中提取样本x_j的邻域像元，获取特征隶属度，根据特征隶属度计算样本所有邻域像元对应第i类特征隶属度的中值

获得特征矩阵中所有样本对应每个类别的空间邻域信息

根据空间邻域信息更新模糊划分矩阵和聚类中心，由模糊划分矩阵U={u_ij}得到各像元分类结果，构造伪彩色图。

具体步骤如下：

1)设置图像分类类别数N、邻域窗口大小C(通常取9或25)、模糊指数m、目标函数精度ε和最大迭代次数T_m；

2)设置邻域信息加权系数α。随机产生初始特征隶属度V¹(N×p矩阵，且必须满足每个样本对应各类别隶属度的和为1)；

3)从图像矩阵IM中提取样本x_j(即特征矩阵X中的每一列的数据)的邻域像元，记第m个邻域像元为SN_m(x_j)(m=1，2，…，C；注：图像有些边缘像元没有C个邻域像元，应另作处理)。从V¹中提取编号为m的邻域像元对应第i(i=1，2，N)类的特征隶属度，记为SNU_im(x_j)，则C个邻域像元对应第i类的特征隶属度集合可记为：

4)求样本x_j所有邻域像元对应第i类特征隶属度的最大值

和最小值 $\min \left(\underset{m=1}{\overset{C}{{\mathrm{SNU}}_i}}\left(x_j\right)\right);$

5)求x_j对应第i类的空间邻域信息，即样本x_j所有邻域像元对应第i类特征隶属度的中值

即：

$v_{\mathrm{ij}}^2=\frac{1}{2}(\max \left(\underset{m=1}{\overset{C}{{\mathrm{SNU}}_i}}\left(x_j\right)\right)+\min \left(\underset{m=1}{\overset{C}{{\mathrm{SNU}}_i}}\left(x_j\right)\right))---\left(5\right)$

6)重复上述步骤3)至5)，求特征矩阵中所有样本对应每个类别的空间邻域信息

并归一化；

7)更新模糊划分矩阵U={u_ij}和聚类中心Z={z_i}，其中：

$u_{\mathrm{ij}}=(1-\mathrm{\α})\×v_{\mathrm{ij}}^1+\mathrm{\α}\×v_{\mathrm{ij}}^2---\left(6\right)$

$z_i={\mathrm{\Σ}}_{j=1}^p{\left(u_{\mathrm{ij}}\right)}^m{x}_j/{\mathrm{\Σ}}_{j=1}^p{\left(u_{\mathrm{ij}}\right)}^m,(i=\mathrm{1,2}...N)---\left(7\right)$

是矢量x_j隶属于第i类的特征隶属度和空间邻域信息

加权后的模糊相似度，满足u_ij∈[0，1]且

z_i为第i类聚类中心，聚类中心Z={z₁，z₂，z₃…z_N}是L×N矩阵。

8)计算目标函数值J(t)，若前后两次目标函数值|J(t)-J(t-1)|<ε或t>T_m就转到步骤9)，否则t=t+1，更新特征隶属度矩阵

并转到步骤3)；其中：

$J\left(t\right)={\mathrm{\&Sigma;}}_{j=1}^p{\mathrm{\&Sigma;}}_{i=1}^N{\left(u_{\mathrm{ij}}\right)}^m{d}^2(x_j,z_i)---\left(8\right)$

$v_{\mathrm{ij}}^1=\frac{1}{{\mathrm{\&Sigma;}}_{k=1}^N{(d(x_j,z_i)/d(x_j,z_k))}^{2/(m-1)}}---\left(9\right)$

m∈(1，∞)控制模糊度的权重指数，它的选取与数据集本身的特点有关；d(x_j，z_i)=||x_j-z_i||为相似性测度，一般选用欧氏距离。

9)由模糊划分矩阵U={u_ij}得到各像元分类结果，构造伪彩色图。

在MATLAB R2009b版本平台下，结合表格中的数据和附图对本发明方案实例进行详细说明。

使用规格为512×614×224的高光谱图像数据，每个像元以两个字节的带符号整型数据表示，并通过仿真实验将本发明方案与K-Means算法、FCM算法进行比较。本发明方案的具体实施步骤如下：

构造图像特征矢量：

1)读取高光谱图像数据源，在空间维中截取1至256行、257至512列、224个波段的图像作为感兴趣图像区域，并得到高光谱图像数据3维矩阵IM(包含256行，256列，224个波段)，简记为IM；

2)将图像矩阵IM每个空间维(每个波段)上的数据按以下方式逐列重排：第二列第一个数据接在IM第一列的最后一个数据后；第三列第一个数据接在IM第二列的最后一个数据后；…；依次类推，构成65536个数据的行向量。最后一共形成224个行向量(即IM中波段数作为矩阵B的第1维数据，而每个波段上的图像数据作为矩阵B的第2维数据)，这样3维矩阵IM就转化为2维矩阵B(224，256×256)。

对图像特征矢量进行部分主成分分析(PPCA变换)：

1)对图像矩阵IM自左向右、自上向下，采用九宫格取中心法提取样本矢量，得到样本矩阵S(大小为224×6554)。

2)计算第i(i=1，2，…，224)个波段所有样本的均值：

(s_ij为第j个样本第i个波段的数据，j=1，2，…，65536)；

3)计算第i个波段所有样本的标准差：

4)利用上述步骤2)和3)求得所有波段的均值和标准差，然后对样本数据进行标准化处理： $s_{\mathrm{ij}}^{\&prime;}=\frac{s_{\mathrm{ij}}-\stackrel{\&OverBar;}{s_i}}{{\mathrm{\&sigma;}}_i}(i=\mathrm{1,2},...,224;j=\mathrm{1,2},...,6554);$

5)为简单起见，将标准化数据表（s′_ij)_224×6554仍记为(s_ij)_224×6554，按公式(10)计算第i个和第j个波段之间的相关系数r_ij(i，j=1，2，…，224)，得到所有波段相互之间的相关系数矩阵R=(r_ij)_224×224；

$r_{\mathrm{ij}}=\frac{1}{6554}\underset{k=1}{\overset{6554}{\mathrm{\&Sigma;}}}\frac{(s_{\mathrm{ik}}-\stackrel{\&OverBar;}{s_i})(s_{\mathrm{jk}}-\stackrel{\&OverBar;}{s_j})}{{\mathrm{\&sigma;}}_i{\mathrm{\&sigma;}}_j}---\left(10\right)$

6)设I为单位矩阵。根据特征方程|R-λI|=0，计算R的特征根λ_i，并使其从大到小排列：λ₁≥λ₂≥…≥λ₂₂₄，同时得到对应的特征向量a₁，a₂，…a₂₂₄；

7)设定累计贡献率的上限值K=0.99，计算贡献率：

8)计算累计贡献率： $E_m=\underset{k=1}{\overset{m}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_m/\underset{k=1}{\overset{224}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_k(m=\mathrm{1,2},...,224),$ 当累计贡献率第一次超过K值时，m=5，即主成分个数L=5；

9)利用前5个主成分计算原图像矩阵IM经PPCA变换后的特征矩阵X={x_ij}，其中：

x_ij=a_1iB_1j+a_2iB_2j+…+a_224iB_224j(i=1，2，…，5，j=1，2，…，65536)   (11)

利用SNFCM算法对特征矩阵进行聚类：

1)设置图像分类类别数N=3、邻域窗口大小C=9、模糊指数m=2、目标函数精度ε=10^-6和最大迭代次数T_m=30；

2)令t＝1，设置邻域信息加权系数α=0.5。利用Matlab软件的随机函数产生初始特征隶属度 $V^1=\left\{v_{\mathrm{ij}}^1\right\}(i=\mathrm{1,2,3};j=\mathrm{1,2},...,65536),$ 行向量代表每个样本对应同一类别的特征隶属度，列向量代表每个样本对应每个类别的特征隶属度；

3)将V¹归一化(必须满足每个列向量的和为1)：

4)按自上向下、自左向右的顺序从图像矩阵IM中提取样本x_j(特征矩阵X中的每一列的数据)的领域像元，记第m个邻域像元为SN_m(x_j)(m=1，2，…，8)；

注：图像中有些边缘像元没有8个邻域像元，只取与它相邻的像元即可，比如：左上角的像元只有右方、下方、右下方3个邻域像元。

5)从V¹中提取编号为m的邻域像元对应第i(i＝1，2，3)类的特征隶属度，记为SNU_im(x_j)，则8个邻域像元对应第i类的特征隶属度集合可记为：

6)求样本x_j所有邻域像元对应第i类特征隶属度的最大值

和最小值 $\min \left(\underset{m=1}{\overset{8}{{\mathrm{SNU}}_i}}\left(x_j\right)\right);$

7)计算x_j对应第i类的空间邻域信息

它的定义为：x_j所有邻域像元对应第i类特征隶属度的中值，即：

$v_{\mathrm{ij}}^2=\frac{1}{2}(\max \left(\underset{C}{\mathrm{SN}}\left(v_{\mathrm{ij}}^1\right)\right)+\min \left(\underset{C}{\mathrm{SN}}\left(v_{\mathrm{ij}}^1\right)\right))---\left(12\right)$

8)重复上述步骤4)至7)，求X={x_ij}中所有样本对应每个类别的空间领域信息 $V^2=\left\{v_{\mathrm{ij}}^2\right\};$

9)将V²归一化(必须满足每个列向量的和为1)：

10)利用式(13)加权特征隶属度矩阵V¹和空间邻域信息V²，得到模糊划分矩阵U={u_ij}；

$u_{\mathrm{ij}}=0.5\&times;(v_{\mathrm{ij}}^1+v_{\mathrm{ij}}^2)---\left(13\right)$

(3×65536矩阵)是矢量x_j隶属于第i类的特征隶属度

和空间邻域信息加权后的模糊相似度，满足u_ij∈[0，1]且

11)根据式(14)计算聚类中心Z={z_i}(5×3矩阵)，z_i为第i类聚类中心；

$z_i={\mathrm{\&Sigma;}}_{j=1}^{65536}{\left(u_{\mathrm{ij}}\right)}^2{x}_j/{\mathrm{\&Sigma;}}_{j=1}^{65536}{\left(u_{\mathrm{ij}}\right)}^2,(i=\mathrm{1,2}...3)---\left(14\right)$

12)当t＞1时，按式(15)计算样本x_j与聚类中心z_i的距离，并与对应的模糊隶属度u_ij进行乘积加权，然后根据式(16)计算目标函数值J(t)，若t＞30或|J(t)-J(t-1)|＜10^-6就转到步骤14)，否则t=t+1；

$d(x_j,z_i)=\sqrt{\underset{h=1}{\overset{5}{\mathrm{\&Sigma;}}}{(x_{\mathrm{hj}}-z_{\mathrm{hi}})}^2}---\left(15\right)$

$J\left(t\right)={\mathrm{\&Sigma;}}_{j=1}^{65536}{\mathrm{\&Sigma;}}_{i=1}^3{\left(u_{\mathrm{ij}}\right)}^2{d}^2(x_j,z_i)---\left(16\right)$

13)根据式(17)更新特征隶属度矩阵

并转到步骤4)；

$v_{\mathrm{ij}}^1=\frac{1}{{\mathrm{\&Sigma;}}_{k=1}^3{(d(x_j,z_i)/d(x_j,z_k))}^2}---\left(17\right)$

14)聚类过程结束，利用U={u_ij}求取每个像元对应的类别号，记录图像分类结果，构造伪彩色图。

高光谱图像分类的主观评价是目前比较常用的评价高光谱图像分类质量的方法，它是从最终分类后的高光谱图像对图像的观察者的心理影响的角度，来量度分类图像优劣的一种方法。主观评价的标准是由专业人员根据事先规定的评价尺度或者自身的经验出发，对最终分类后的高光谱图像质量做出判断。表1给出了国际上通用的五级妨碍标准以及质量标准，这种方法简称为主观评价标准的五级评分制。

表1高光谱图像质量主观评价等级V值采用五级评分制

评分标准	质量标准	妨碍标准
			优	非常好	丝毫不能觉察图像的质量变坏
良	好	能觉察图像的质量变坏，但并不讨厌
			中	一般	清楚地觉察图像质量变坏，但可以容忍
差	差	明显觉察图像质量变坏，令人讨厌
			劣	非常差	极其明显觉察图像质量变坏，很讨厌

以某AVIRIS高光谱图像的部分区域作为仿真测试源数据，通过算法用时和分类效果图将本发明方案与K-Means算法、FCM算法进行比较，整个实验在InterCore双核2.20GHz、1.8GB内存，WindowsXP、MATLAB环境下实现。表1列出了三种方法在五种相同迭代次数下分别对同一高光谱图像分类的算法用时。其中，高光谱图像分类类别数为3；FCM算法的权重系数m=2；本发明方案选用前5个主成分作为特征矢量进行分类，权重系数m=2，空间邻域信息加权系数α=0.5。

表2三种算法的时耗比较

由表2的实验结果可以看出，本发明方案效率远远高于K-Means算法，平均用时约为传统K-Means算法的2.89％，几乎可以忽略不计；虽然本发明方案中空间邻域信息加权的步骤会导致增加计算量，但是由于部分主成分分析的作用减少了样本的特征个数，整体方案复杂度较传统FCM算法仍然有所降低，平均用时约为传统FCM算法的87.78％。

传统K-Means算法和FCM算法的分类效果相当，主观评分为“差”，两种算法分类后的伪彩色图中都存在着“杂点”，三种地物类别的边界比较模糊。经过本发明方案处理后的结果图，主观评分为“优”，明显可以看出，“杂点”现象得到了很好的抑制，各类别的边界变清晰了，而且随着加权系数的增大，抑制效果越来越好，但到一定程度的时候，分类效果就变化不大了。结合表1可以得到一个结论：本发明提出的部分主成分分析方法和空间邻域信息加权的FCM算法，使得整体方案在降低波段冗余、提高算法效率的同时，可以充分利用空间邻域信息来抑制传统算法分类结果中的“杂点”现象，具有很高的工程应用价值。随着高光谱图像应用的范围越来越广泛，高光谱图像处理技术一直备受关注，遥感分类技术更是其中研究的重点领域，利用上述优势，本发明方案可以很好应用于植被覆盖、水体检测等工程领域。

Claims (8)

1.一种空间邻域信息加权的高光谱遥感图像分类方法，其特征在于，该方法包括：读取遥感图像，对图像特征矢量进行PPCA变换实现对特征矢量矩阵降维，从波段数据中获得部分线性不相关的矢量作为像元样本；从图像矩阵中提取像元样本的邻域像元，获取特征隶属度，根据特征隶属度计算像元样本所有邻域像元对应第某类特征隶属度的中值，作为特征矩阵中所有样本对应每个类别的空间邻域信息，根据空间邻域信息更新模糊划分矩阵和聚类中心，由模糊划分矩阵得到各像元分类结果，构造伪彩色图。 

2.根据权利要求1所述的方法，其特征在于，对特征矢量矩阵降维包括：选取样本并转换成数据矩阵，计算样本数据的相关系数矩阵，计算相关系数矩阵的特征根，计算特征根的累计贡献率，根据累计贡献率选取特征维数，获取特征矢量。 

3.根据权利要求1所述的方法，其特征在于，根据特征隶属度计算像元样本所有邻域像元对应第某类特征隶属度的中值具体为：根据公式： 计算x_j对应第i类特征隶属度的中值，其中，和为样本x_j所有邻域像元对应第i类特征隶属度的最大值和最小值。 

4.根据权利要求1所述的方法，其特征在于，根据空间邻域信息更新模糊划分矩阵和聚类中心具体包括：从图像矩阵中提取样本x_j的邻域像元；获得x_j隶属于第i类的特征隶属度

和空间邻域信息获得特征隶属度

和空间邻域信息

加权后的模糊相似度

以及聚类中心 

更新模糊划分矩阵U={u_ij}和聚类中心Z={z_i｝。 

5.根据权利要求2所述的方法，其特征在于，计算相关系数矩阵的特征根具体为：根据数据矩阵S的相关系数矩阵R，调用特征方程|R-λI|=0计算R的特征根λ_i。 

6.根据权利要求2所述的方法，其特征在于，计算特征根的累计贡献率，根据累计贡献率选取特征维数具体为：将特征根从大到小排列，得到对应的特征向量a={a₁，a₂，…a_ω｝，根据公式：

计算累计贡献率。 

7.根据权利要求2或6所述的方法，其特征在于，根据累计贡献率选取特征维数，获取特征矢量具体为：找到使累计贡献率第一次到达累计贡献率上限值K的特征维数L，提取前L个主成分对应的特征向量β={a₁，a₂，…a_L}，调用公式X=β^T·B获取图像的特征矩阵，即根据公式： 

X_ij=a_1iB_1j+a_2iB_2j+…+a_ωiB_ωj(i＝1，2，…，L，j=1，2，…，p)计算原图像数据的主成分。 

8.根据权利要求4所述的方法，其特征在于，设定邻域窗口大小C，以当前样本像元为中心，与它四周相邻的C个像元即为邻域像元。 

Images