光谱分辨率在10-2λ数量级范围内的光谱图像称为高光谱图像(Hyperspectral Image)。
高光谱遥感的发展得益于成像光谱技术的发展与成熟。成像光谱技术最大特点是将成像技术与光谱探测技术结合，在对目标的空间特征成像的同时，对每个空间像元经过色散形成几十个乃至几百个窄波段以进行连续的光谱覆盖。
这样形成的数据可以用“三维数据块”来形象地描述。其中x和y表示二维平面像素信息坐标轴，第三维（λ轴）是波长信息坐标轴。
高光谱图像集样本的图像信息与光谱信息于一身。图像信息可以反映样本的大小、形状、缺陷等外部品质特征，由于不同成分对光谱吸收也不同，在某个特定波长下图像对某个缺陷会有较显著的反映，而光谱信息能充分反映样品内部的物理结构、化学成分的差异。
高光谱图像数据
遥感影像具有四种分辨率：时间分辨率、空间分辨率、辐射分辨率和光谱分辨率。
其中，对于光谱分辨率而言，按照成像传感器波谱通道划分数目的多少，分为多光谱、高光谱和超光谱。
通常认为，多光谱波段数目在100以下，高光谱波段数目在100~10000之间，超光谱波段数目在10000以上。
高光谱图像就是好多个灰度图像叠加到一起，每一个灰度图代表了一个光谱波段。
举个例子来说，一般的二维图像表示，这个图像像素为255*255，也就是说这幅图像有255*255个像素，可以是灰度图或者彩色图，相关信息包含在每一个像素中。
而同样大小的高光谱图像，若包含200个光谱波段的信息，那就该表示为255*255*200，可以理解为每一个像素上有200维的光谱域信息，就类似于，200幅255*255的二维图像叠加在一起，200幅图像中相同位置像素的灰度值画成曲线表示出来便是这一像素点的光谱域信息了。
也就是说，高光谱图像不仅包含丰富的光谱域信息，同时也跟一般的二维图像一样，包含相同的空间域信息。
包含Indian Pines、Salinas、Pavia Centre and University、Cuprite、 Kennedy Space Center、Botswana、anomaly detection（7个地区），在每个地区都会有介绍类型数量、波段数，可直接下载。
这里，只下载了Botswana的数据，如下图所示：
高光谱图像数据的读取
1.Matlab读取.Mat格式的高光谱遥感数据：
load('Botswana.mat')InputMatImg=Botswana;b = size(InputMatImg);fprintf('输入图像宽度为 %d\n',b(1));fprintf('输入图像高度为 %d\n',b(2));fprintf('输入图像波段数为 %d\n',b(3));i=120;j=180;k=220;%自选三个波段InputImg_r= InputMatImg(:,:,i);%获取第i个波段InputImg_g= InputMatImg(:,:,j);%获取第j个波段InputImg_b= InputMatImg(:,:,k);%获取第k个波段InputImg_r= uint8(InputImg_r);%将i波段的灰度值转为0~255InputImg_g= uint8(InputImg_g);%将j波段的灰度值转为0~255InputImg_b= uint8(InputImg_b);%将k波段的灰度值转为0~255RGBImg=cat(3,InputImg_r,InputImg_g,InputImg_b);%将i、j、k三个波段进行合成figure;subplot(221);imshow(InputImg_r);title('红色波段');subplot(222);imshow(InputImg_g);title('绿色波段');subplot(223);imshow(InputImg_b);title('蓝色波段');subplot(224);imshow(RGBImg);title('合成波段');imwrite(InputImg_r,['MATBand',num2str(i),'.jpg']);imwrite(InputImg_g,['MATBand',num2str(j),'.jpg']);imwrite(InputImg_b,['MATBand',num2str(k),'.jpg']);imwrite(RGBImg,'compositeMATRGBimg.jpg');
2.python读取.Mat格式的高光谱遥感数据：
.mat是数据文件。mat数据文件是matlab的数据存储标准格式。mat数据文件是标准的二进制文件.
scipy是构建在numpy的基础之上的，它提供了许多的操作numpy的数组的函数。
https://scipy.io/
(点击下
方阅读原文可跳转
)
包提供了多种功能来解决不同格式的文件的输入和输出。
def testLoadmat(is_plot): """
作用： loadmat()函数用于加载.mat数据文件
参数: file_name: str .mat文件的文件名（当appendmat = True时，不用加.mat后缀），也可以打开类似文件的对象（file-like object）
mdict: dict, 可选，插入文件变量的字典
appendmat: bool , 可选，若为真则在文件名后添加后缀
byte_order: str / None, 可选 默认为none。表示从mat文件中猜测的字节顺序，可以是“native”，“=”，“little”，“”，“BIG”之一。
mat_dtype: bool，可选 若为真，则返回与加载到MATLAB中相同的dtype的数组，而不是保存数组的dtype。
squeeze_me: bool，可选，判断是否压缩单位矩阵的维数
chars_as_string: bool，可选，将char数组转string数组
matlab_compatible: bool，可选 返回被MATLAB读取的数组（当squeeze_me = false, chars_as_string = false, mat_dtype = true, struct_as_record = true）
struct_as_record: bool，可选 设置flag来判断加载MATLAB以numpy记录数组还是以原形式numpy数组（dtype为对象）。
verify_compressed_data_integrity: bool，可选，MATLAB文件的长度是否已确认。
variable_names: None / sequence 若为None（默认），则读取文件中的所有变量。否则variable_names将是一个string序列，表示需要从mat文件中读取的变量名。此时读取器将跳过不是这个变量名的序列，一定程度上减少读取时间。
返回: mat_dict: dict 以变量名为key，数组为values的字典
# 加载数据X = sio.loadmat('data/Botswana/Botswana.mat') print(X)此时输出的X是一个字典：{'__header__': b'MATLAB 5.0 MAT-file, Platform: GLNXA64, Created on: Thu Feb 20 15:12:42 2014', '__version__': '1.0', '__globals__': [], 'Botswana': array([[[3996, 3952, 3698, ..., 60, 53, 54], [8471, 8255, 8288, ..., 0, 0, 0]]], dtype=uint16)}若要获取其中的数据部分Botswana，则可直接输出X[‘Botswana’]，
Python代码如下：
import scipy.io as sio
# 加载高光谱数据和其对应的标签X = sio.loadmat('data/Botswana/Botswana.mat')['Botswana'] # 取字典中该项的值y = sio.loadmat('data/Botswana/Botswana_gt.mat')['Botswana_gt']
# 打印高光谱数据及其标签的形状print('Hyperspectral data shape: ', X.shape) # (1476, 256, 145)print('Label shape: ', y.shape) # (1476, 256)
高光谱图像数据的预处理及使用（python）
PCA（Principal Component Analysis）
是一种常见的数据分析方式，即主成分分析技术，又称主分量分析，常用于高维数据的降维，把多指标转化为少数几个综合指标,可用于提取数据的主要特征分量。
在统计学中，主成分分析PCA是一种简化数据集的技术。它是一个线性变换。这个变换把数据变换到一个新的坐标系统中，使得任何数据投影的第一大方差在第一个坐标(称为第一主成分)上，第二大方差在第二个坐标(第二主成分)上，依次类推。主成分分析经常用于减少数据集的维数，同时保持数据集的对方差贡献最大的特征。这是通过保留低阶主成分，忽略高阶主成分做到的。这样低阶成分往往能够保留住数据的最重要方面。
关于 PCA降维的数学原理可参考这篇文章，非常详细。
【机器学习】降维——PCA（非常详细） - 知乎 (zhihu.com)
(点击下
方阅读原文可跳转
)
代码如下：
def applyPCA(X, numComponents): newX = np.reshape(X, (-1, X.shape[2])) # 行数未知，列数为145 将高光谱数据展开成一维数据序列(377856, 145) # print(newX.shape) pca = PCA(n_components=numComponents, whiten=True) newX = pca.fit_transform(newX) newX = np.reshape(newX, (X.shape[0], X.shape[1], numComponents)) # 将降维后的数据按照指定的形状排列 return newX
# 对高光谱数据 X 应用 PCA 变换X_pca = applyPCA(X, numComponents=pca_components) print('Data shape after PCA: ', X_pca.shape)
sklearn.decomposition.PCA(n_components=None, copy=True, whiten=False)
意义：PCA算法中所要保留的主成分个数n，也即保留下来的特征个数n
类型：int 或者 string，缺省时默认为None，所有成分被保留。
类型：bool，True或者False，缺省时默认为True。
意义：表示是否在运行算法时，将原始训练数据复制一份。
若为True，则运行PCA算法后，原始训练数据的值不会有任何改变，因为是在原始数据的副本上进行运算；
若为False，则运行PCA算法后，原始训练数据的值会改，因为是在原始数据上进行降维计算。
PCA方法：fit_transform(X)
对部分数据先拟合fit，找到该part的整体指标，如均值、方差、最大值最小值等等，然后对该X进行转换transform，从而实现数据的标准化、归一化等等。
用X来训练PCA模型，同时返回降维后的数据。
newX=pca.fit_transform(X)，newX就是降维后的数据。
# 对单个像素周围提取 patch 时，边缘像素就无法取了，因此，给这部分像素进行 padding 操作def padWithZeros(X, margin=2): newX = np.zeros( (X.shape[0] + 2 * margin, X.shape[1] + 2 * margin, X.shape[2])) x_offset = margin y_offset = margin newX[x_offset:X.shape[0] + x_offset, y_offset:X.shape[1] + y_offset, :] = X return newX
# 在每个像素周围提取 patch ，然后创建成符合 keras 处理的格式def createImageCubes(X, y, windowSize=5, removeZeroLabels=True): # 给 X 做 padding margin = int((windowSize - 1) / 2) zeroPaddedX = padWithZeros(X, margin=margin) # split patches patchesData = np.zeros( (X.shape[0] * X.shape[1], windowSize, windowSize, X.shape[2])) patchesLabels = np.zeros((X.shape[0] * X.shape[1])) patchIndex = 0 for r in range(margin, zeroPaddedX.shape[0] - margin): for c in range(margin, zeroPaddedX.shape[1] - margin): patch = zeroPaddedX[r - margin:r + margin + 1, c - margin:c + margin + 1] patchesData[patchIndex, :, :, :] = patch patchesLabels[patchIndex] = y[r - margin, c - margin] patchIndex = patchIndex + 1 if removeZeroLabels: patchesData = patchesData[patchesLabels > 0, :, :, :] patchesLabels = patchesLabels[patchesLabels > 0] patchesLabels -= 1 return patchesData, patchesLabels
# 在每个像素周围提取 patch ，然后创建成符合 keras 处理的格式X_pca, y = createImageCubes(X_pca, y, windowSize=patch_size) print('Data cube X shape: ', X_pca.shape)print('Data cube y shape: ', y.shape
对单个像素周围提取 patch 时，边缘像素就无法取了，因此，给边缘像素进行 padding 操作。在原有数据的基础上给周围增加了7列0元素。
同时，对原有数据进行取样，取样大小为15*15.在一个特征点上，共取出了
1476*256=
377856个样本，最后将其中标签为0的去掉，共得到3248个样本，数据块大小为
(3248, 15, 15, 30)。
创建训练数据和测试数据集:
def splitTrainTestSet(X, y, testRatio, randomState=345): X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=testRatio, random_state=randomState, stratify=y) return X_train, X_test, y_train, y_test
Xtrain, Xtest, ytrain, ytest = splitTrainTestSet(X_pca, y, test_ratio)print('Xtrain shape: ', Xtrain.shape)print('Xtest shape: ', Xtest.shape)
train_test_split随机划分训练集和测试集：
X_train,X_test, y_train, y_test =train_test_split(
train_data,train_target,test_size=0.25, random_state=0,stratify=y)
# test_size：样本占比，如果是整数的话就是样本的数量
# random_state：是随机数的种子,随机数种子：其实就是该组随机数的编号，在需要重复试验的时候，保证得到一组一样的随机数。比如你每次都填1，其他参数一样的情况下你得到的随机数组是一样的。但填0或不填，每次都会不一样。
# stratify: 是为了保持split前类的分布，通常在这种类分布不平衡的情况下会用到stratify。比如有100个数据，80个属于A类，20个属于B类。如果train_test_split(... test_size=0.25, stratify = y), 那么split之后数据如下：
training: 75个数据，其中60个属于A类，15个属于B类。
testing: 25个数据，其中20个属于A类，5个属于B类。
用了stratify参数，training集和testing集的类的比例是 A：B= 4：1，等同于split前的比例（80：20）。
将stratify=X就是按照X中的比例分配 ，将stratify=y就是按照y中的比例分配 。
改变 Xtrain, Ytrain 的形状，以符合 keras 的要求：
Xtrain = Xtrain.reshape(-1, patch_size, patch_size, pca_components, 1)Xtest = Xtest.reshape(-1, patch_size, patch_size, pca_components, 1)print('before transpose: Xtrain shape: ', Xtrain.shape)print('before transpose: Xtest shape: ', Xtest.shape)
添加了一个维数，即通道数。
为了适应 pytorch 结构，数据要做转置：
Xtrain = Xtrain.transpose(0, 4, 3, 1, 2)
