一类基于图像块距离的高光谱图像流形降维方法

1.一类基于图像块距离的高光谱图像流形降维方法，其特征在于，采用改进的LLE算法或改进的ISOMAP算法；其中使用一种新的图像块距离度量，该距离度量将观测像素看作高维流形上的一点，同时兼顾到观测点周围的空间结构：
对于大小为 且波段数为 的高光谱数据集 来说，假设以观
测像素 和 为中心的 空间邻域 和 内的像素点组成的集合分别为
， ，其中 ，集合 和 的大小均为 ；
则观测像素 和 之间的图像块距离（IPD）定义为：
（1）
其中，图像块大小为 ， ，
，距离函数 为两观测像素之间的光谱角距离、欧式距
离或者Kullback–Leibler散度距离， 相应的定义分别如下：
（2）
（3）
（4）
其中， ， ， ， ， 、 、 和 分别表示向量 、
、和 的第 个元素；
上述距离度量表征了以观测像素为中心的图像块之间的差别；
设已知水平宽度、垂直宽度和波段数分别为 、 和 的高光谱图像数据矩阵，邻接点个数 ，数据集内在维数 ；
所述改进的LLE算法的步骤如下：
步骤一：计算图像块距离矩阵
（a）图像镜像扩展：对应于大小为 的图像块，原始数据扩充为
，使得处于边缘和四角的像素也可以使用图像块距离；
（b）对于原始数据集上的任意两个观测像素点 和 ，根据式（1）计算图像块距离，得到图像块距离矩阵 ，其中 ；
步骤二：在全样本点中寻找每个样本点的 个邻近点；
步骤三：根据原始的LLE算法，通过求解最优化问题得到约束的权重矩阵，得到低维嵌入坐标；
所述的改进ISOMAP算法的步骤如下：
步骤一：计算图像块距离矩阵
（a）图像镜像扩展：对应于大小为 的图像块，原始数据扩充为
，使得处于边缘和四角的像素也可以使用图像块距离；
（b）对于原始数据集上的任意两个观测像素点 和 ，根据式（1）计算图像块距离，得到图像块距离矩阵 ；
步骤二：在全样本点中寻找每个样本点的 个邻近点，并基于 构建邻域图，计算最短路径矩阵和执行MDS算法，得到最优嵌入结果。

一类基于图像块距离的高光谱图像流形降维方法\n技术领域\n[0001] 本发明属于遥感图像处理技术领域，具体涉及一种高光谱遥感图像非线性降维方法。\n背景技术\n[0002] 遥感是本世纪六十年代发展起来的新兴综合技术，与空间、电子光学、计算机、地理学等科学技术紧密相关，是研究地球资源环境的最有力的技术手段之一。高光谱遥感是将成像技术与光谱技术相结合的多维信息获取技术。高光谱成像仪在电磁波谱的数十至数百个非常窄且连续的光谱段上同时探测目标的二维几何空间与一维光谱信息。高光谱图像中，每一个观测像素都可以提取出一条完整连续的光谱曲线，为地物信息的提取和分析提供了极其丰富的信息，有助于更加精细的地物分类和目标识别。然而，波段数的增多必然导致了信息的冗余和数据处理复杂性的增加。同时，较高的光谱分辨率在增强地物细微差别分辨能力的同时，也带来了维数灾难（Hughes现象），这种现象严重影响了高光谱图像的处理效果。高光谱图像特征降维成为解决这种现象的常用方法，它对于高光谱图像分类等应用具有特殊的意义，在高光谱图像处理中具有十分重要的作用[1]。\n[0003] 高光谱数据降维技术是以图像特征提取为目的，利用低维数据来有效地表达高维数据特征的数据处理技术。它在有效地保留了图像信息的同时也大大减少了信息的冗余，更有利于信息的快速提取。常见的高光谱图像降维算法可以分为线性降维和非线性降维两大类[2]~[7]。主成分分析（Principal Component Analysis，PCA）[2]是一种最常用的线性降维方法。它的主要目标是通过线性变换寻找一组最优的单位正交向量基，并用它们的线性组合来重构原样本，以使重建后的样本和原样本的误差最小。其它代表性的线性降维算法还有独立成分分析（Independent Component Analysis，ICA）[3]，线性判别分析（Linear Discriminant Analysis，LDA）[4]等算法。\n[0004] 下面介绍与本发明相关的一些概念：\n[0005] 流形学习算法\n[0006] 流形学习算法是一种常见的非线性降维方法，它是基于这样的假设：高维数据在特征空间中对应的点分布在“低维流形”上。因此，流形学习算法实现降维的目的是寻找原始数据在“低维流形”上的嵌入坐标。代表性的流形学习算法有局部线性嵌入（Locally Linear Embedding，LLE）[5]算法，等距映射（Isometric Feature Mapping，ISOMAP) [6]算法和拉普拉斯特征映射（Laplacian Eigenmap，LE）[7]算法等。作为一种局部性保持算法，LLE认为数据流形具有局部线性，即一个数据点可以通过其邻域完全重建，于是可以通过在降维空间中尽可能保持其局部线性特征来实现降维。ISOMAP算法则是一种通过保持流形上两点间的测地线距离来保持数据集的全局几何特性。它保证了降维结果的稳健性和全局最优性，但是其运算复杂度较高。邻接点个数 （或者邻域距离 ）和低维数据的维度（固有维度） 是LLE和ISOMAP算法的两个主要参数。两种算法的步骤如下：\n[0007] 算法：局部线性嵌入算法（Locally Linear Embedding, LLE）\n[0008] 输入：邻域大小 ，内在维度 ，高维原始数据集 ，其中为数\n据点个数， 为数据维度。\n[0009] 输出：低维映射空间数据集\n[0010] 步骤一： 构建原始数据集的欧式距离矩阵 ，并寻找每个样本点的 个邻近点。\n[0011] 步骤二： 通过最小化如下目标函数\n[0012] (5)\n[0013] 得到约束的权值矩阵 ，其中 表示 的第 个邻近点。\n[0014] 步骤三： 通过最小化目标函数\n[0015] （6）\n[0016] 计算最优嵌入结果 , 同样 表示 的第 个邻近\n点。\n[0017] 算法：等距映射算法（Isometric Feature Mapping , ISOMAP）\n[0018] 输入：高维原始数据集 ，邻域大小 或者邻域距离\n，内在维度 ；\n[0019] 输出：低维映射空间数据集 ；\n[0020] 步骤一： 构建邻域图。根据原始数据集的欧式距离矩阵 ，寻找每个样本点的 个邻近点（ 邻域方法）或者距离小于 的点（ 邻域方法），从而得到原始数据的邻域图。\n[0021] 步骤二：计算最短路径。根据邻接图中两点及其邻域之间的欧氏距离计算两点之间的最短路径，从而得到最短路径矩阵 。该最短路径即为估计得到的测地线距离。\n[0022] 步骤三：利用多维尺度分析（MDS）方法得到保存测地线距离的最优低维度映射结果 。\n发明内容\n[0023] 本发明的目的在于提出一类新的可以有效地去除原始数据集中的信息冗余、保持数据集的局部空间结构的基于图像块距离的高光谱图像流形降维方法。\n[0024] 本发明提出的基于图像块距离的高光谱图像流形降维方法，采用一种新的距离度量——图像块距离，并将其应用于流形学习算法，从而得到一组新的非线性降维方法。本发明将高光谱图像数据看成包含低维嵌入结构的高维流形，结合一种新的距离度量和流形学习算法实现高光谱数据的降维。这是一种非线性降维算法，可以有效地去除原始数据集中的信息冗余（包括光谱维和空间维冗余）。与传统的基于线性模型的线性降维算法（如PCA算法）相比，该发明显著地考虑到了数据集的空间结构，同时也可以很好的处理高光谱数据集中存在的非线性现象。与数据驱动的流形学习算法相比，图像块距离的引入使得改进的算法更好的发现嵌入在高维流形中的低维结构，同时在降维过程中也保持了数据集的局部空间结构。\n[0025] 本发明包括四种改进的算法：基于图像块角距离和欧式距离的LLE算法（Image Patch Angle Distance and Euclidean Distance-based LLE, IPAD-LLE and IPED-LLE）、基于图像块角距离和欧式距离的ISOMAP算法（Image Patch Angle Distance and Euclidean Distance-based ISOMAP, IPAD-ISOMAP and IPED-ISOMAP）。与线性降维算法（如PCA算法）以及其他流形学习算法（如原始的LLE和ISOMAP算法）相比，在应用于高光谱图像分类时，使用本发明可以得到相对较高的分类精度。\n[0026] 本发明提出的基于图像块距离的高光谱图像流形降维方法，具体内容如下：\n[0027] 1、图像块距离度量\n[0028] 在高光谱图像中，观测像素矢量包含特定位置的所有波段的像素值，不仅是频谱相关而且是空间相关[8]。例如，如果一个给定的观测像素属于某一特定的类，那么它的空间上邻近的点属于该类的概率就特别高。基于此，在本节中我们提出一种新的距离度量并将其应用于流形学习算法中。作为一种距离，它既将观测像素看作高维流形上的一点，又兼顾到观测点周围的空间结构。\n[0029] 假设以观测像素 为中心的大小为 的方形空间区域为 ，其中 为\n大于1的奇数。我们称这样的邻域为 的 空间邻域。考虑到观测像素的空间关系时，我们利用观测像素的空间邻域计算观测像素之间距离。因为这种新的距离量度表征了图像块之间的差别，我们称之为图像块距离（Image Patch Distance）。而且，基于这种距离测度的流形算法也被称为基于图像块距离的流形学习算法。\n[0030] 对于大小为 且波段数为 的高光谱数据集 来说，任意\n两观测像素 和 之间的“距离”定义为 ，该距离可以为光谱角距离、欧式距离或者Kullback–Leibler散度距离等，其中观测像素间的光谱角距离定义为：\n[0031] （7）\n[0032] 光谱角距离既可以表征两个观测像素的距离，又可以避免强度等因素产生的误差。\n[0033] 假设以观测像素 和 为中心的 空间邻域 和 内的像\n素点组成的集合分别为 ， ，\n其中 而且集合 和 的大小均为 。则观测像素 和 之间的图像块\n距离定义为：\n[0034] （8）\n[0035] 其中 ， 。该距离表征了以观测像素为中心的大小为\n的图像块之间的差别。\n[0036] 需要说明的是，为了在图像的边缘和四角也可以使用到新的图像块距离，在初始化阶段首先要做的是镜像扩展图像，即将原始高光谱数据 以边缘为轴\n镜像扩展扩充为 ，然后在此基础上得到原始数据\n集对应的图像块距离矩阵 ，其中 为观测像素个数。最后，根据预置的邻\n域个数选择邻域并应用于LLE或者ISOMAP算法实现非线性降维。\n[0037] 此外，虽然在式（7）中我们采用的是光谱角距离，但是距离函数 可以方便的换成其他距离定义，如欧氏距离、Kullback–Leibler 散度距离等，它们的定义分别如下：\n[0038] （9）\n[0039] （10）\n[0040] 其中 ， ， ， ， 、 、 和 分别表示向量 、 、\n和 的第 个元素。\n[0041] 算法复杂度\n[0042] 本文提出的图像块距离矩阵的计算有两部分组成。第一部分来自观测点距离矩阵的计算，例如，如果使用的是欧式距离（定义如（9））则为欧式距离矩阵的计算。对于经过扩展的高光谱数据集 来说，如果采用光谱角距离，则计算光谱角\n距离矩阵所需flops为 ，而如果采用欧氏距离，计算欧式距离矩阵所需flops为 ，其中 。它们都和 相关。第二部分则来自\n于式（2）中涉及的最大最小计算，以及大量的内存索引操作。由此可知，本节提出的图像块距离矩阵的计算与数据量以及波段数有关，这与欧式距离矩阵的复杂度是一致的，而且从理论上来说，计算图像块距离矩阵的复杂度并不明显高于欧式距离矩阵的计算复杂度，也就是引入图像块距离矩阵后，运算复杂度并没有明显的提升。\n[0043] 根据上述图像块距离度量，本发明采用改进的LLE算法和改进的ISOMAP算法，具体分别介绍如下：\n[0044] 1 、改进的LLE算法\n[0045] 已知水平宽度、垂直宽度和波段数分别为 、 和 的高光谱图像数据矩阵，邻接点个数 ，数据集内在维数 ，所述改进LLE算法的步骤如下：\n[0046] 步骤一：计算图像块距离矩阵\n[0047] （a）图像镜像扩展：对应于大小为 的图像块，原始数据扩充为\n，使得处于边缘和四角的像素也可以使用图像块距离；\n[0048] （b）对于原始数据集上的任意两个观测像素点 和 ，根据式（8）计算图像块距离 ，得到图像块距离矩阵 。\n[0049] 步骤二：在全样本点中寻找每个样本点的 个邻近点。\n[0050] 步骤三：类似于原始的LLE算法，通过求解最优化问题得到约束的权重矩阵，进而得到低维嵌入坐标。\n[0051] 改进的ISOMAP算法\n[0052] 已知水平宽度、垂直宽度和波段数分别为 、 和 的高光谱图像数据矩阵，邻接点个数 ，数据集内在维数 ，所发明的改进ISOMAP算法的步骤如\n下：\n[0053] 步骤一：计算图像块距离矩阵\n[0054] （a）图像镜像扩展：对应于大小为 的图像块，原始数据扩充为\n，使得处于边缘和四角的像素也可以使用图像块距\n离；\n[0055] （b）对于原始数据集上的任意两个观测像素点 和 ，根据式（8）计算图像块距离 ，得到图像块距离矩阵 。\n[0056] 步骤二：在全样本点中寻找每个样本点的 个邻近点。\n[0057] 步骤三：类似于ISOMAP算法，基于 构建邻域图、计算最短路径矩阵和执行MDS算法得到最优嵌入结果。\n[0058] 本发明的优点\n[0059] 本发明包含一类新的基于图像块距离的高光谱图像流形降维方法。改进的流形算法结合高光谱数据的物理特性，更好地发现嵌入在高光谱数据集中的低维流形结构，在应用到高光谱分类时，可以更好地去除高光谱数据集中冗余的空间信息和光谱维度信息。实际高光谱数据实验表明，本文所提出改进的流形学习算法在应用到高光谱图像的分类时能够取得良好的效果，并且性能大大优于其它已有的高光谱图像降维方法。因此该算法具有较大的实际意义。\n附图说明\n[0060] 图1 (a) Indiana Pine 数据集中波段70的灰度图, (b) 包含16类第五的真实地物分布图。\n[0061] 图2不同降维维度下七种算法对应的分类结果：（a） KNN,（b） SVM。\n[0062] 图3不同训练集大小下七种算法对应的分类结果：（a） KNN,（b） SVM。\n[0063] 图4 七种算法对应的分类图 (d=20)：(a) PCA, (b) LLE, (c) ISOMAP, (d) IPAD-LLE, (e) IPAD-ISOMAP, (f) IPED-LLE, (g) IPED-ISOMAP。\n具体实施方式\n[0064] 下面，分别用仿真数据和实际遥感图像数据为例说明本发明的具体的实施方式:\n[0065] 我们将四种改进算法，即IPAD-LLE、IPED-LLE、IPAD-ISOMAP 和 IPED-ISOMAP，与PCA和原始的LLE、ISOMAP算法进行比较。它们都是应用于高光谱数据降维的常见的且性能较优的算法。为了比较不同的降维算法的性能，我们在降维基础之上，利用分类算法对降维结果执行分类操作，通过分析分类的精度评价这七种算法。采用的分类算法为最近邻分类法（K-Nearst Neighbourhood，KNN）[9]和支撑向量机（Support Vector Machine，SVM）[10]。评价分类结果的指标为总体分类精度（Overall Accuracy, OA）。\n[0066] 所使用的Indiana Pine 数据集为1992年6月拍摄于Indiana州西北方向的一片农田的高光谱遥感数据 。该数据集大小为145×145（共21025个观测像素），共包含220个波段，波谱范围从400nm~2450nm，光谱分辨率为10nm，空间分辨率为17m。覆盖该区域的主要是各种农作物（包括玉米、小麦、大豆、干草堆）和天然植被（树林、草地等）以及一些人工用地。该数据集被广泛地用于遥感图像的研究中。Purdue大学给出一份关于该区域的实地调查报告供参考[11]。在进行处理之前，将信噪比太低的波段以及水吸收波段（包括波段\n104~108、150~163和220）去除，剩下200个波段被用于进一步处理。图1（a）给出该数据的波段70的灰度图而图1（b）为该地区对应的16类地物真实分布。\n[0067] 图2给出了不同维度下的七种算法的性能。实验中流形学习算法所采用的邻接点个数为8。改进的算法中，图像块的大小为 。使用的SVM代码来自与LIBSVM库[12]，其中参数为： ， ，而KNN代码则来自于MATLAB的Bioinformatics 工具箱，其中用于分类的临近点的个数为1。执行分类时，我们随机选取1/4数据作为训练样本，剩下的作为测试样本。\n[0068] 由图2可知，使用新的邻域选择算法之后，改进的LLE和ISOMAP算法的性能在不同的降维维度下都得到了大幅度的提高且均高于PCA和原始的LLE、ISOMAP算法。这是由于一方面该实际数据集中存在较为明显的非线性，而改进的流形学习算法可以很好地找到嵌入其中的线性结构，因而得到较优的结果。另一方面，原始的LLE和ISOMAP算法，由于仅仅考虑到数据集的内在几何结构而未顾及实际图像的空间信息，因而在应用于实际数据时效果较差，甚至还不如PCA算法。\n[0069] 为了测试不同训练集大小下，各个算法的降维结果对应的分类精度，我们在给出了总体精度关于训练集大小的曲线图（如图3）。其他条件如上。由图3可知，改进的算法的性能要优于PCA和原始的LLE、ISOMAP算法。此外结合图2和图3，我么还可以看出IPAD-ISOMAP算法的性能总是最优的。\n[0070] 为了更精确地说明分类的效果，表1给出了降维维度为20、图像块大小为 、训练样本所占比例为1/4时，16类地物的具体分类结果，采用的分类算法为KNN。此外，我们给出相应的分类图（图4）。由表1和图4可知，与PCA、原始LLE和ISOMAP算法相比，改进的算法得到的结果中，各个类的分类精度都得到了较大的提高。\n[0071] 表 1 七种算法对应的KNN分类结果比较（d=20）\n[0072] \n[0073] 参考文献\n[0074] [1]Chang C.I. Hyperspectral imaging: techniques for spectral detection and classification [M]. New York: Plenum, 2003.\n[0075] [2]Jolliffe L.T. Principal Component Analysis [M]. Springer, 2nd ed. \n2002.\n[0076] [3]Hyvarinen, Karhunen J, and Oja E. Independent Component Analysis [M]. New York: Wiley, 2001.\n[0077] [4]Belhumeur P.N, Hespanha J.P, and Kriegman D.J. Eigenfaces vs. Fisherfaces: Recognition Using Class Specific Linear Projection [J]. IEEE Trans. on Pattern Analysis and Machine Intelligence, 19(7):711-720, July 1997.[0078] [5]Sam T.R and Lawrence K.S. Nonlinear Dimensionality Reduction by Locally Linear Embedding [J].Science, 2323-2326, Dec 22 2000.\n[0079] [6]Tenenbaum J. B, de Silva.V, and Langford J.C. A global geometric framework for nonlinear dimensionality reduction [J]. Science, \n290(5500):2319–2323, 2000.\n[0080] [7]Belkin M, Niyogi P. Laplacian eigenmaps for dimensionality reduction and data representation [J].Neural Computation, 15(6), 2003.\n[0081] [8]Mohan A, Sapiro G, and Bosch E. Spatially Coherent Nonlinear Dimensionality Reduction and Segmentation of Hyperspectral Images [J]. IEEE Trans. on Geosci. Remote Sens., 4(2):206-210, Apr. 2007.\n[0082] [9]Cover T, and Hart P. Nearest neighbor pattern classification [J]. IEEE Trans. on Inf. Theory, 13(1): 21-27, Jan. 1967.\n[0083] [10]Melgani F, and Bruzzone L. Classification of hyperspectral remote sensing images with support vector machines [J]. IEEE Trans. on Geosci. Remote Sens., 42(8):1778-1790, Aug. 2004.\n[0084] [11]Landgrebe D. Multispectral data analysis: A signal theory perspective. West Lafayette: School of Electrical & Computer Engineering, Purdue University, pp.56-89, 1998.\n[0085] [12]http://www.csie.ntu.edu.tw/cjlin/libsvm/。