| 高光谱图像子空间的波段选择 |  [PDF全文] |
收稿日期: 2017-11-27
2. College of Information and Communication Engineering, Dalian Minzu University, Dalian 116600, China
与传统的遥感数据相比，高光谱遥感数据具有近乎连续的光谱特征和纳米级的光谱分辨率。光谱分辨率的提高一方面可以增强地物分辨能力，另一方面也带来了数据冗余，导致数据分析时产生较高的计算复杂度以及维数灾难现象(Pal和Foody，2010)。因而高光谱遥感数据分析前，降维预处理常常是十分必要的。
降维方法一般可分为两类：特征提取和波段选择。特征提取是用映射的方法将原始数据变换为较少的新特征，包括主成分分析PCA(Principal Component Analysis)(Agarwal 等，2007)、独立成分分析ICA(Independent Component Analysis)(Wang和Chang，2006)、局部线性嵌入LLE(Locally Linear Embedding)(Li 等，2012)等方法。与特征提取不同，波段选择可以在保留原始数据的物理意义及光谱特性的同时降低数据维度，因而被广泛应用。其中，依目的波段搜索方法相较其他无监督的依单波段指标排序(Chang和Wang，2006)、波段聚类(Sun 等，2015)等方法，前者在有、无监督情况下均适用，因此本文主要研究此类方法。
现有的波段搜索方法分为最优搜索算法与次优搜索算法。最优搜索算法是利用穷举法在光谱空间中遍历所有波段组合，虽然可以找到最优解，但实施困难。因而实际应用中，一般使用次优搜索算法(Pudil 等，1994)。次优搜索算法以准则函数为评价依据，通过特定的方法从原始特征集中选择一组性能较好，但不一定是最好的特征组合。已有的次优搜索算法包括顺序前向选择法、顺序后向选择法、顺序前向浮动选择法和顺序后向浮动选择法(Serpico和Bruzzone，2001)、最速上升法、遗传算法(赵冬和赵光恒，2009)、拟态物理学算法(王立国和魏芳洁，2013)、人工蜂群算法ABC(Artificial Bee Colony algorithm)(王立国 等，2015)等。顺序前向选择法和顺序后向选择法的运算速度较快，但这两种方法所选波段组合冗余较多，性能较差，并且一旦某一波段被选(或剔除)，就无法更改。顺序前向浮动选择法和顺序后向浮动选择法是对顺序前向选择法和顺序后向选择法的改进，能够动态地改变选入或者剔除的波段数，但计算复杂。最速上升法较顺序前向浮动选择法可获得性能更好的波段组合，但计算复杂度提高，并且结果对初始值敏感，随着维度的增加其性能会呈现下降趋势，鲁棒性较差。遗传算法、拟态物理学算法、人工蜂群算法均属于智能优化方法，这类方法虽可以针对目标函数获得具有倾向性的波段组合，但一般要设置较多参数，导致模型复杂性提高。
鉴于以上方法的优势与不足，本文提出一种只需输入待选波段子集大小，而无需设置其他参数的波段搜索策略，同时为进一步降低算法的计算复杂度，将子空间划分引入该方法中，通过逐子空间搜索替换的方式获取波段子集来满足监督与无监督两种情况下的降维需求。
2、算法描述本文方法预先对波段子空间进行划分，得到相关性较弱的几个波段子空间，而后利用提出的搜索方法根据相应搜索准则在各波段子空间内搜索波段子集。
(2.1) 子空间划分 2.1.1. 划分依据更加精细的光谱划分是高光谱数据的特点，这使得数据在光谱域具有较强相关性，这里的谱间相关性是指：空间上同一位置的地物，相邻波段的波段图像具有相似性。产生谱间相似的原因主要为：同一地物在相邻波段的光反射率是非常相近的，因此产生了一定的相关性。这种相关性可以用相关系数矩阵来描述(谷延锋和张晔，2003)，以AVIRIS采集的印第安纳农林数据为例，计算其相关系数矩阵和相关系数向量，并将得到的矩阵和向量进行可视化(图1)：
|
| 图 1 Indian数据谱间相关性的可视化 Figure 1 Visualization of spectral correlation of Indian Pines |
图1(a)是以灰度图像的形式呈现，由灰度图像的取值特点可知，越明亮的区域相关系数越大，而明亮区域主要集中于主对角线，因此可以说明相邻波段间的相关性更强，而从图1(b)可以直观看到相关性较强的各个波段范围。依据高光谱数据的这种波段聚集特性，可以将这些具有较大相关性的波段进行划分，以便加速后续处理。
2.1.2. 自动子空间划分方法自动子空间划分的方法主要是根据波段相关系数矩阵图像分块特性及近邻可传递相关性来进行高光谱数据空间划分，具体步骤如下：
步骤1 将二维波段图像转换为一维的波段向量；
步骤2 计算所有波段的相关系数得到高光谱数据的相关系数矩阵R，其定义为R=[r1,1，r1,2，…，r1, j；r2,1，r2,2，…，r2, j；…；r j,1，r j,2，…，r j,j]；
步骤3 从相关矩阵中提取出近邻可传递相关矢量rNTR，其定义为rNTR=[r1,2，r2,3，…r i,i+1，…r j−2,j−1，r j−1,1]T，对近邻可传递相关矢量进行处理得到c−1个局部相关的极小值；
步骤4 根据得到的c−1个极小值将高光谱数据空间划分c个适合的数据子空间。
经过划分后可以得到不同维度的子空间，每个子空间内的波段数据具有相近的光谱特性。
(2.2) 搜索准则波段搜索准则与波段的评价手段相关，针对不同的评价方式可以选择不同的搜索准则或目标函数，对于高光谱遥感数据，波段子集的评价主要包括信息量、各波段间相关性以及类别可分性。其中以最佳指数、波段子集平均相关性和整体分类精度为目前常用评价指标，也可根据实际应用进行调整。
为验证搜索算法的有效性，将以上述各评价指标分别作为目标搜索准则，下面介绍各准则的计算方法。
(1)平均相关性。波段子集平均相关性的计算方式如式(1)(王立国和魏芳洁，2013)所示。
| $\bar R = \frac{1}{{{M^2}}}\sum\limits_{i = 1}^{N - 1} {\sum\limits_{j = i + 1}^N {{R_{ij}}} } $ | (1) |
式中，N为波段总数，M为样本总数，Rij为波段i和波段j的相关系数，其计算方式如式(2)(谷延锋和张晔，2003)所示。
| ${R_{ij}} = \frac{{\sum\limits_{k = 1}^M {\left({{x_{ik}} - {{\bar x}_i}} \right)} \left({{x_{jk}} - {{\bar x}_j}} \right)}}{{\sqrt {\sum\limits_{k = 1}^M {{{\left({{x_{ik}} - {{\bar x}_i}} \right)}^2}\sum\limits_{k = 1}^M {{{\left({{x_{jk}} - {{\bar x}_j}} \right)}^2}} } } }}$ | (2) |
式中，
(2)最佳指数。最佳指数(Chavez 等，1982)采用了波段的标准差与波段间相关性的比值进行计算，其计算方式如式(3)所示。
| ${O_{{\rm{OIF}}}} = \frac{{\sum\limits_{i = 1}^N {{S_i}} }}{{\sum\limits_{i = 1}^N {\sum\limits_{j = 1}^N {\left| {{R_{ij}}} \right|} } }}$ | (3) |
式中，Si是第i波段的标准差，Rij同上。波段的标准差是衡量波段信息量的一种手段，标准差越大则信息量越大。另如上文所述，波段间相关性越小冗余越小，综合这两个值来看，波段子集的最佳指数越大越好。
(3)总体分类精度。高光谱数据的总体分类精度计算方式如式(4)(童庆禧 等，2006)所示。
| ${O_{{\rm{OA}}}} = \frac{{\sum\limits_{i = 1}^C {{m_{ii}}} }}{M}$ | (4) |
式中，mii为
2.2节中介绍了波段选择所追求的几种目标，波段搜索最终就是实现在所选波段子集下某目标函数值的最大或者最小化。为此，算法首先在各波段子空间提取一个波段(本文选择每个子空间中方差最大的波段)，构成大小为K的波段子集，并相应地计算由该初始子集下的目标函数值。然后，固定K−1个波段子空间，只变化余下的一个子空间中的波段，用该子空间中的波段依次替换当前选择的波段，如果某个替换能够得到更大(或更小)的目标函数值，那么该次替换就作为有效替换得以保留，否则作为无效替换而被淘汰。当得到有效替换后，跳出该子空间，进入下一子空间重复这样的择优替换过程，直到没有任何替换能够使得目标函数值增加(或减小)为止。此时，所选择的波段组合则为最终输出的波段子集。
值得注意的是，当所需波段数目大于根据2.1节中计算得到的子空间个数时，需要对子空间进行二次划分。这里，将2.1节介绍的划分称为一级划分，一级划分得到的子空间称为一级子空间；对一级子空间进行再划分的过程称为二级划分，得到的各子空间称为二级子空间。二级划分方法分为两步：首先，按照各一级子空间大小的比例分配每个一级子空间所需选择的波段数目；然后，对每个一级子空间按待选波段数目进行均分，得到二级子空间。
(2.4) 算法总体流程现将2.1—2.3节所介绍内容总结为下面的步骤：
步骤1 根据所需的波段数目K对数据在光谱空间进行波段子空间划分；
步骤2 在各子空间内各选取一个方差最大的波段进行波段子集初始化；
步骤3 计算初始波段子集下的目标函数值G；
步骤4 固定K−1个子空间的已选波段，在剩余的一个子空间中按顺序对波段子集进行更新，计算新的波段子集下目标函数值Gg；
步骤5 比较G与Gg的大小，若G<Gg则返回步骤4，若G>Gg则保留该次替换，同时按标号顺序更换进行迭代替换的波段子空间；
步骤6 重复步骤5，直至没有任何一次替换使得目标函数改变，输出最终的波段子集，算法结束。
3、实验与分析为验证本文方法的可行性与有效性，进行了以下实验:与基于ABC算法的波段选择方法进行有监督与无监督两种情况下的对比；与基于反向传播AP(Affinity Propagation)算法的波段选择，自适应波段选择ABS(Adaptive Band Selection)算法进行无监督情况对比。基于AP算法的波段选择是基于聚类方法的波段选择，而ABS算法则是同时考量了信息量与相关性的波段选择方法。
实验环境为AMD双核处理器，主频2.47 Hz，有效内存8 GB，开发环境为Matlab R2008a。实验数据为去除噪声波段的200波段的AVIRIS印第安纳农林数据和103波段的ROSIS帕维亚大学数据：
(1)印第安纳农林数据。该数据的波长范围为0.4—2.5 µm，空间分辨率为17 m，共有144×144个像素点。数据中剔除背景共包含16类地物，主要农作物是生长期的玉米和大豆，结合地面实际测量数据，其中7种地物样本量过少，对于该数据不具有代表性，因此选取另9种样本数目较多的代表性地物作为实验用地物。
(2)帕维亚大学数据。该数据的波长范围为0.43—0.86 µm，空间分辨率为1.3 m，共有610×340个像素点，共包含9类地物，实验中9种地物均作为实验用地物。
两组数据所对应的假彩色图像与真值数据如图2、图3所示，9种地物类型及数目如表1所示。
|
| 图 2 Indian数据 Figure 2 Land-covers for Indian Pines |
|
| 图 3 PaviaU数据 Figure 3 Land-covers for University of Pavia |
|
表 1 地物类别及数目表 Table 1 Land-covers for Indian Pines and University of Pavia |
印第安纳农林数据按照2.1节进行子空间划分后得到5个子空间，分别为：(1—36)、(37—79)、(80—103)、(104—144)、(145—200)，于每个子空间内各选择一个波段，获得相关性较低的5波段组合。
对于多光谱数据通过穷举法可以获取最优波段子集，但对于高光谱数据，这是一个惊人的计算过程。以印第安纳农林数据为例，5个子空间的大小分别为：36、43、24、41、56，若遍历所有的组合方式，则要进行85300922次运算，这在实际应用中很难实现。而本文方法对目标函数进行搜索时可通过少量计算得到较优的波段子集，下面进行具体的实验分析。
(1)搜索次数比较。从本文方法搜索3种目标函数时所计算的组合次数中数据可以看出(表2)，本文方法通过计算远少于穷举组合的次数的波段组合方式即可以收敛到一个较优的波段组合。
|
|
表 2 搜索次数比较 Table 2 Comparison of searching times |
(2)指标评价比较。将本文方法与ABC算法，AP算法与ABS算法所得波段子集的各目标函数值与计算时间列于表3、表4中(表中数据均是进行20次实验的统计平均结果)。表中“总体分类精度(目标函数)”一项为本文方法与ABC算法在有监督情况下将分类精度作为目标函数进行波段选择时的总体分类精度。将总体分类精度作为目标函数时，监督数据的选择方式为各地物均匀抽取200个像素点，再将200个像素点均分为2份，每份100个像素分别作为此过程的训练与测试数据。AP算法与ABS算法为无监督波段选择方法，因此在这一项不做比较。表中“总体分类精度(最终评价)”一项为各方法所获取的波段子集在整个数据集上进行地物分类的结果，这里训练样本同样均匀抽取200像素点，剩余像素点为测试数据。而表4中只对本文方法与ABC方法针对不同目标函数进行搜索时间的对比，而不对只根据单一准则进行搜索的ABS算法与AP算法进行比较。
|
|
表 3 各方法的指标评价 Table 3 Indexes evaluation of different method |
|
|
表 4 波段搜索时间对比 Table 4 Comparison of computational time for searching band |
从表3中可以看出针对不同需求进行搜索的本文方法与ABC算法相较AP算法与ABS算法可以得到更适应需求的波段组合。这是因为AP算法与ABS算法均是按照单一的准则进行子集选取的，在面向不同的需求时，无法灵活应对。
通过表4进一步分析本文方法与ABC算法的搜索效率可以看出本文方法在搜索时间上的花费更少，效率更高。并且本文方法所搜索得到的波段组合具有唯一解，而基于ABC算法所得到的波段组合则具有随机性，这是由ABC算法的随机初始化特性所决定的。
(3.2) 帕维亚大学数据实验帕维亚大学数据经过子空间划分得到3个子空间，分别为：(1—73)、(74—84)、(85—103)，于每个子空间内各选择一个波段，获得相关性较低的3波段组合。这3个子空间的大小分别为：73、11、19，若穷举所有波段组合方式，则有15257种组合。
采用同样的对比实验进行本文方法有效性的验证，将实验结果分别列于表5、表6、表7中。
|
|
表 5 搜索次数比较 Table 5 Comparison of searching times |
|
|
表 6 各方法的指标评价 Table 6 Indexes evaluation of different method |
|
|
表 7 波段搜索时间对比 Table 7 Comparison of computational time for searching band |
(1)搜索次数比较。将本文方法搜索4种目标函数时所计算的组合次数列于表5中，从表中数据可以看出，本文方法通过计算远少于穷举组合的次数的波段组合方式即可以收敛到一个较优的波段组合。
(2)指标评价比较。与印第安纳农林数据实验的方式相同，将本文方法与ABC算法，AP算法与ABS算法进行比较分析，每种方法各进行20次实验，将统计平均的结果列于表5、表6中。从表中可以看出本文方法与ABC算法在有针对性的应用时相较另外两种方法能获得性能更优的波段子集。同时，在搜索效率上本文方法与ABC算法相比更高。
通过对两个实际高光谱数据集的实验结果可以看出，本文方法在两组数据上可以得到一致的性能评价，证明了本文方法的有效性。
4、结 论本文针对高光谱数据光谱空间特征冗余问题，提出了一种子空间下的波段选择方法，该方法可以按照实际需要快速地选择出一个性能优良的波段子集，算法设置参数极少，并且通过子空间划分既可以降低波段组合的冗余度，同时又能加速搜索过程。针对两种实际高光谱数据集进行实验，结果证明了该方法的可行性与有效性。
未来的研究工作可以从以下两个方向进行展开：(1)本文在所需波段数目较多时，采用的是按比例对一级子空间进行平均再分的方式，这种方式简单易行，但是没有进一步挖掘一级子空间的光谱特性来进行二次划分，因此对初步子空间划分后的各波段区间进行再分析是一个重要的研究方向；(2)文中各子空间中波段替换方式是按顺序进行的，然而面对不同应用时波段子集的初始化应采用不同的方法以加速收敛，因此寻求更灵活的替换机制同样是一个待研究的问题。
志 谢 此实验的数据为普度大学与帕维亚大学所提供的公开数据集，在此衷心地表达感谢。
| [1] | Agarwal A, El-Ghazawi T, El-Askary H and Le-Moigne J. 2007. Efficient hierarchical-PCA dimension reduction for hyperspectral image//Proceedings of 2007 IEEE International Symposium on Signal Processing and Information Technology. Giza, Egypt: IEEE: 353-356 |
| [2] | Chang C I and Wang S. Constrained band selection for hyperspectral imagery[J]. IEEE Transactions on Geoscience and Remote Sensing, 2006, 44 (6) : 1575 –1585. DOI: 10.1109/TGRS.2006.864389 |
| [3] | Chavez P S, Berlin G L and Sowers L B. Statistical method for selecting Landsat MSS ratios[J]. Journal of Applied Photographic Engineering, 1982, 8 (1) : 23 –30. |
| [4] | 谷延锋, 张晔. 基于自动子空间划分的高光谱数据特征提取[J]. 遥感技术与应用, 2003, 18 (6) : 384 –387. Gu Y F and Zhang Y. Feature extraction based on automatic subspace partition for hyperspectral images[J]. Remote Sensing Technology and Application, 2003, 18 (6) : 384 –387. DOI: 10.11873/j.issn.1004-0323.2003.6.384 |
| [5] | Li W, Prasad S, Fowler J E and Bruce L M. Locality-preserving dimensionality reduction and classification for hyperspectral image analysis[J]. IEEE Transactions on Geoscience and Remote Sensing, 2012, 50 (4) : 1185 –1198. DOI: 10.1109/TGRS.2011.2165957 |
| [6] | Pal M and Foody G M. Feature selection for classification of hyperspectral data by SVM[J]. IEEE Transactions on Geoscience and Remote Sensing, 2010, 48 (5) : 2297 –2307. DOI: 10.1109/TGRS.2009.2039484 |
| [7] | Pudil P, Novovičová J and Kittler J. Floating search methods in feature selection[J]. Pattern Recognition Letters, 1994, 15 (11) : 1119 –1125. DOI: 10.1016/0167-8655(94)90127-9 |
| [8] | Serpico S B and Bruzzone L. A new search algorithm for feature selection in hyperspectral remote sensing image[J]. IEEE Transactions on Geoscience and Remote Sensing, 2001, 39 (7) : 1360 –1367. DOI: 10.1109/36.934069 |
| [9] | Sun K, Geng X R and Ji L Y. Exemplar component analysis: a fast band selection method for hyperspectral imagery[J]. IEEE Geoscience and Remote Sensing Letters, 2015, 12 (5) : 998 –1002. DOI: 10.1109/LGRS.2014.2372071 |
| [10] | 童庆禧, 张兵, 郑兰芬. 2006. 高光谱遥感——原理、技术与应用. 北京: 高等教育出版社 Tong Q X, Zhang B and Zheng L F. 2006. Hyperspectral Remote Sensing. Beijing: Higher Education Press |
| [11] | Wang J and Chang C I. Independent component analysis-based dimensionality reduction with applications in hyperspectral image analysis[J]. IEEE Transactions on Geoscience and Remote Sensing, 2006, 44 (6) : 1586 –1600. DOI: 10.1109/TGRS.2005.863297 |
| [12] | 王立国, 魏芳洁. 结合APO算法的高光谱图像波段选择[J]. 哈尔滨工业大学学报, 2013, 45 (9) : 100 –106. Wang L G and Wei F J. Artificial physics optimization algorithm combined band selection for hyperspectral imagery[J]. Journal of Harbin Institute of Technology, 2013, 45 (9) : 100 –106. DOI: 10.11918/j.issn.0367-6234.2013.09.018 |
| [13] | 王立国, 赵亮, 刘丹凤. 基于人工蜂群算法高光谱图像波段选择[J]. 哈尔滨工业大学学报, 2015, 47 (11) : 82 –88. Wang L G, Zhao L and Liu D F. Artificial bee colony algorithm-based band selection for hyperspectral imagery[J]. Journal of Harbin Institute of Technology, 2015, 47 (11) : 82 –88. DOI: 10.11918/j.issn.0367-6234.2015.11.014 |
| [14] | 赵冬, 赵光恒. 基于改进遗传算法的高光谱图像波段选择[J]. 中国科学院研究生院学报, 2009, 26 (6) : 795 –802. Zhao D and Zhao G H. Band selection of hyperspectral image based on improved genetic algorithm[J]. Journal of the Graduate School of the Chinese Academy of Sciences, 2009, 26 (6) : 795 –802. |