ENVI 混合像元分解

一、混合像元分解的过程

在影像已经完成预处理的前提下（如几何校正、大气校正、去噪等），混合像元分解的一般的过程：首先获取端元波谱（从图像上、波谱库中或者其他来源），然后选择一种分解模型在每个像素中获取每个端元波谱的相对丰度图，最后从丰度图上提取不同组成比例的像元。

混合像元分解基本概念：
1.端元（endmember）
组成混合像元的纯净地物被称作是端元。

2.混合像元产生的原因：1）传感器的采样频率小于或等于地表空间变异的频率，使得单个像元对应地面的面积较大，故包含了多种地物类型；2）相对于地表的变异程度，传感器的空间分辨率过于精细，使得有些像元不可避免地横跨地物边界。

3.确定端元类型的方法：1）从光谱库中选择；2）直接从影像上提取；3）结合光谱库和影像上提取的波谱，先从影像上选择相对单一的像元的光谱，然后对光谱进行分解，通过与光谱库中的光谱进行对比，确定端元的类型和波谱值。

4.最大噪声比变换（MNF）：是一种连续使用主成分分析（PCA）对高光谱遥感数据进行压缩和降维的方法。

5.像元纯度指数（PPI）：是表征多波段遥感图像中每个像元“纯度”的指标，该值越大，说明对应的像元越接近纯像元。原理：它将 N 维像元点投影到一个随机的单位向量上，如果该像元纯度大，则应该更接近单位向量的端点，否则位于单位向量的内部。通过这种投影方式迭代多次后，纯度大的像元靠近端点的概率值大，纯度小的像元靠近端点的概率值小，于是可得到一幅反映端元纯度大小的影像。

二、混合像元分解算法

混合像元分解算法包括数据降维、端元选取和反演三个步骤。

1.数据降维
尽管数据降维不是混合像元分解算法的一个必需步骤，但由于大多数算法都将其作为一个流程，我们也将其当作一个步骤。常用的降维算法有主成分分析（Principle Component Analysis，PCA）、最大噪声比变换（Maximum Noise Fraction，MNF）和奇异值分解（Singular Value Decomposition，SVD）。

(1) 主成分分析：遥感图像各波段之间经常是高度相关的，因此所有的波段参加分析是不必要的。PCA就是一种去除波段之间相关性的变换。PCA通过对原数据进行线性变换，获得新的一组变量，即主成分。其中前几个主成分包含了原数据主要方差，同时各个主成分之间是不相关的。

(2) 最大噪声比变换：最大噪声比变换（Maximum Noise Fraction，MNF）由Green等（1989）提出，该变换通过引入噪声协方差矩阵以实现对噪声比率的估计。首先，通过一定方式（比如对图像进行高通滤波）获取噪声的协方差矩阵，然后将噪声协方差矩阵对角化和标准化，即可获得对图像的变换矩阵，该变换实现了噪声的去相关和标准化，即变换后的图像包含的噪声在各个波段上方差都为1，并且互不相关。最后对变换后的图像再做主成分变换，从而实现了MNF变换，此时得到的图像的主成分的解释方差量对应于该主成分的信噪比大小。

(3) 奇异值分解：奇异值分解（Singular Value Decomposition, SVD）也是遥感图像处理中常用的变换，与PCA类似，SVD能够找出包含原始数据大部分方差的特征方向，不同的是，SVD特别适合于波段间高度相关的数据，而PCA在这种情况下很有可能会失败。

2.端元选取
选取合适的端元是成功的混合像元分解的关键。端元选取包括确定端元数量以及端元的光谱。理论上，只要端元数量m小于等于L+1（L表示波段数），线性方程组就可以求解。然而实际上由于端元波段间的相关性，选取过多的端元会导致分解结果更大的误差，尽管此时残差会减少。在能够描述一个场景内光谱的大部分方差的前提下，越少的端元数量是越好的选择。对于城市地区，最常用的端元选取方式是由Ridd等（1995）提出的植被-不透水层-土壤端元模型（Vegetation - Impervious surface–Soil, V-I-S），V-I-S模型在很多研究中得到应用。而在非城市地区，一般采用植被-土壤-阴影（或干植被）端元模型。

端元的数量和类型确定后，下一步是确定端元的光谱。端元光谱的确定有两种方式：(1) 使用光谱仪在地面或实验室测量到的“参考端元”；(2) 在遥感图像上得到的“图像端元”。参考端元虽然可以精确测量，但由于各种因素（包括不同传感器、大气影响、辐射条件及物候等）造成的噪声，可能会导致其与图像上像元光谱的不一致。要将二者匹配起来需要进行复杂的校正，而且参考端元的微弱噪声就可能引起最后计算得到的端元比例有很大误差。相对而言，直接从图像上寻找端元更加直接方便，因而得到广泛研究。图像端元选取的方法大致可以分为两类：交互式端元提取和自动端元提取。

2.1交互式端元提取
从图像上选取端元的各种方法，都基于这样一种思想：在特征空间中，所有的混合像元都存在于由端元连接而成的多边形（或多面体）内。这样的混合像元才能满足端元面积比例为正值并且总和为1的条件。因此，最简单的交互式提取方法就是在特征空间中（通常是前两个或三个主成分构成的特征空间）目视寻找多边形的顶点作为端元。为了减少目视选取的主观性，一些定量化方法被引进作为选取端元时的参考，这些方法的特点是仍然需要人工参与，所以被称为交互式提取。

(1)PPI指数：像元纯度指数（Pixel Purity Index）是最成功的方法之一。首先对图像进行MNF变换以实现数据降维，接着在由MNF的前几个主要成分组成的特征空间中，随机生成穿过数据云的测试向量，然后将数据点投影到测试向量上。投影在测试向量两端的数据点有较大的可能性属于端元，用一个阈值选出在这个测试向量两端的极值点。继续生成新的随机向量，重复上述步骤，记录图像中每个像元作为极值点的频度，即为PPI指数。PPI指数越高意味着像元的纯度也越高。

(2)MEST算法：MEST（Manual Endmember Selection Tool）算法通过主成分分析来确定混合物中端元的数目。对于三维及三维以内的数据而言，可以直接通过目视选取多边形的顶点来确定端元，但对于更高维数的数据，显示上存在困难。MEST提供了一种在高维空间中寻找端元光谱的方法。

(3)CAR和EAR指数：由于端元内的光谱差异，一类端元往往对应很多条光谱，这两个指数用于解决从多条光谱中选择最具代表性的光谱问题。CAR（Class Average RMSE）计算A类端元光谱用B类端元光谱来分解产生的残差,显然残差越小，A、B两类端元的混淆性越大。EAR（Endmember average RMSE）计算A端元内某一条光谱用A端元内其他光谱来分解产生的残差。显然EAR越低表明这条光谱的代表性越好，如果很高则证明这条光谱可能是离群点，没有代表性。因此，利用CAR和EAR指数可以为端元的合并或分离提供依据，进而指导端元选取。

2.2自动端元提取
自动端元提取指不需要人工参与的端元提取方法。相对于交互式提取方法，自动方法的提取结果不受个人的主观性影响，但是由于自动提取方法一般采用纯数学判据，有时可能产生不具有物理意义的端元。常见的自动提取方法有：N-FINDR、IEA、CCA和ORASIS。

3.反演
最简单的线性混合模型反演过程是无约束的最小二乘法，该方法假设观测误差w为高斯白噪声，此时估计的比例有解析表达式：F=(X’X)-1X’r。考虑到端元比例的物理意义，需要满足约束条件Σfi=1，则可通过拉格朗日乘子法求解部分约束最小二乘法（abundance sum-to-one constraint，ASC），也有解析解。Heinz等提出全约束最小二乘（fully constrained least squares，FCLS）的求解方法。FCLS在原方程组里添加方程Σfi=1，并给予其很大的权重，再用ANC对新的方程组求解，则求得的解近似满足两个约束条件。由于该方法能够同时满足两个约束条件以及高效的运算效率，因而得到广泛认可。但是也有学者认为，准确的混合像元模型的分解结果应当自动满足两个约束条件，强加的约束条件可能会导致新的误差。

三、混合像元分解模型

混合像元分解模型分为：线性模型、概率模型、模糊分析模型、几何光学模型和随机几何模型。

四、线性光谱分离

Linear Spectral Unmixing可以根据物质的波谱特征，获取多光谱或高光谱图像中物质的丰度信息，即混合像元分解过程。假设图像中每个像元的反射率为像元中每种物质的反射率或者端元波谱的线性组合。例如：像元中的25%为物质A，25%为物质B，50%为物质C，则该像元的波谱就是三种物质波谱的一个加权平均值，等于0.25A+0.25B+0.5C，线性波谱分离解决了像元中每个端元波谱的权重问题。

线性波谱分离结果是一系列端元波谱的灰度图像（丰度图像），图像的像元值表示端元波谱在这个像元波谱中占的比重。比如端元波谱A的丰度图像中一个像元值为0.45，则表示这个像元中端元波谱A占了45%。丰度图像中也可能出现负值和大于1的值，这可能是选择的端元波谱没有明显的特征，或者在分析中缺少一种或者多种端元波谱。

利用线性光谱分解，根据物质的光谱特性，确定多光谱或高光谱图像中描述物质的相对丰度。

假设图像各像素处的反射率为像素内各物质(或端元)反射率的线性组合。因此，对于给定的光谱(输入数据)和端元光谱，利用线性光谱分解可以求解每个像素中每个端元的丰度值。

端元数必须小于谱带数，图像中所有的端元数必须被使用。光谱解混结果高度依赖于输入端元。

五、ENVI实现混合像元分解步骤

1.MNF变换
（1）原理
最小噪声分离变换（Minimum Noise Fraction Rotation，MNF Rotation）工具用于判定图像数据内在的维数（即波段数），分离数据中的噪声，减少随后处理中的计算需求量。MNF本质上是两次层叠的主成分变换。第一次变换（基于估计的噪声协方差矩阵）用于分离和重新调节数据中的噪声，这步操作使变换后的噪声数据只有最小的方差且没有波段间的相关。第二步是对噪声白化数据（Noise-whitened）的标准主成分变换。

MNF变换具有PCA变换的性质,是一种正交变换,变换后得到的向量中的各元素互不相关,第一分量集中了大量的信息,随着维数的增加,影像质量逐渐下降,按照信噪比从大到小排列,而不像PCA变换按照方差由大到小排列,从而克服了噪声对影像质量的影响。正因为变换过程中的噪声具有单位方差,且波段间不相关,所以它比PCA变换更加优越。

（2）作用
作用：分离数据中的噪声，减少之后计算的需求量，即数据降维。

（3）操作步骤
在ENVI主菜单中，选择/Transform/MNF Rotation/Forward MNF Estimate Noise Statistics 。

2.计算PPI像元纯度指数
（1）原理
纯像元指数(Pure Pixel Index, PPI)，将投影到大量随机向量端点次数最多的几个像元作为端元。

像元纯度指数算法把每个像元作为一个n维向量, 所有像元就组成了一个向量空间V. 在这个向量空间中基并不唯一, 必然存在全部由位于边界位置的向量组成的基，可以用它们的线性组合来表示所有其他的向量。这些处在边界位置的向量在投影到随机产生的单位向量上时, 出现在随机单位向量边缘的概率最大, 如果有大量的随机单位向量, 这种概率将以频率的形式表现出来。

当图像光谱被视为n维光谱空间中的点时，端元光谱应位于数据云的边缘。光谱点被投影到每个测试向量上，并且在最小和最大投影值的阈值距离内的光谱被标记为极值(如这里的二维图所示)。

（2）操作步骤
在ENVI主菜单中，选择 Spectral-> Pixel Purity Index->[FAST] New Output Band。在打开的Pixel Purity Index Input File对话框中，选择MNF变换结果。（根据情况选择波段，单击Spectral Subset按钮，选择波段）

设置Pixel Purity Index Parameters：

Number of Iterations：迭代次数
Threshold Factor：在数据单元中输入极限像素选择的阈值因子值。例如，值2将来自极端像素(包括高像素和低像素)的所有大于两个DN值的像素标记为极端。这个阈值选择投影向量两端的像素。阈值应该大约是数据中的噪声水平的两到三倍。例如，陆地卫星TM数据的噪声通常小于1 DN，因此阈值为2或3就很好。当使用MNF数据时，将噪声归一化，阈值2或3工作得很好。较大的阈值会导致PPI找到更极端的像素，但它们不太可能是纯端元。

设置阈值，显示感兴趣区：

3.构建n维可视化窗口
（1）定义
n维可视化工具在一个n维空间的数据集中定位、识别和聚类最纯净的像素和最极端的光谱响应(端元)。

（2）操作步骤
在ENVI主菜单中，选择Spectral ->n-Dimensional Visualizer，在n-D Visualizer Input File对话框中选择MNF变换结果，单击OK。

选择1-2-3波段，构建3维散点图。

设置适当的速度（Speed），单击Start按钮，在n-D Visualizer窗口中的点云随机旋转，当在n-D Visualizer窗口中的点云有部分聚集在一块时，单击Stop按钮。
在n-D Visualizer窗口中，用鼠标左键勾画“白点”集中区域，选择的点被标示颜色。
在n-D Controls面板中，选择Class->Items 1:20->White（用于删除点），单击Start按钮，当看到有部分选择的点云分散时候，单击Stop按钮，在n-D Visualizer窗口中选择分散的点，自动会将选择的点删除。借助<-，->，New按钮可以一帧帧从不同视角浏览以辅助删除分散点。

输出端元波谱：输出Spectral Library（或者ASCII），将端元波谱保存为波谱库文件或者文本文件。

4.波谱分离
可参考：http://blog.sina.com.cn/s/blog_764b1e9d0101dnjz.html

（1）在Toolbox中，打开/Classification/Endmember Collection工具，在文件对话框中选择影像；

（2）在Endmember Collection面板中，选择Import >from ASCII file。将显示的端元波谱选中，点击OK；

（3）点击Apply按钮，在弹出的对话框选择输出路径，点击OK执行混合像元分解。

分解结果：

总结：
预处理后的影像经过MNF变换，达到数据降维的目的；然后计算PPI像元纯度指数，设置阈值；利用LSMM线性光谱分离模型或fcls_spectral_unmixing扩展工具进行波谱分离，最后得到不同类型的丰度图。

参考链接：
ENVI扩展工具：完全约束最小二乘法混合像元分解：http://blog.sina.com.cn/s/blog_764b1e9d0101dnjz.html

定量/高光谱遥感之——混合像元分解：
https://blog.csdn.net/lijie45655/article/details/49184529

混合像元分解中的端元波谱获取方法：
http://blog.sina.com.cn/s/blog_764b1e9d0100qtv9.html

高级光谱分析：
http://blog.sina.com.cn/s/blog_764b1e9d0102v564.html

GEE进行混合像元分解：
https://zhuanlan.zhihu.com/p/129881408

Linear Spectral Unmixing：
https://www.l3harrisgeospatial.com/docs/linearspectralunmixing.html#:~:text=Linear%20Spectral%20Unmixing%20has%20two%20constraint%20options%3A%20unconstrained,values%20and%20are%20not%20constrained%20to%20sum-to-unity%20%28one%29.

Spectral Tools References：
https://www.l3harrisgeospatial.com/docs/spectraltoolsreferences.html

Collect Endmember Spectra：
https://www.l3harrisgeospatial.com/docs/collectingendmemberspectra.html#Endmembe2

Spectral Unmixing：
https://www.sciencedirect.com/topics/earth-and-planetary-sciences/spectral-unmixing
~~共同学习，共同进步！