为了克服线性光谱混合模型的上述不足, 许多学者利用非线性模型对高光谱图像的光谱混合方式进行了描述。

鉴于线性混合是非线性混合在光学多次散射被忽略情况下的特例^[11] , 本文定义了广义的非线性光谱混合模型, 其表达式为

y=f(M, a)+n, (1)

式中: y=[y₁, …, y_L]^T为L波段的高光谱图像; f是一未知的非线性函数向量; M=[m₁, …, m_R]为端元光谱矩阵; a=[a₁, …, a_R]^T为端元的丰度向量(R为端元的个数); n为0均值且方差为σ ²的高斯噪声向量。

由于f函数是未知的, 因此基于非线性光谱混合模型的光谱解混是个病态问题, 需要给出合适且确定的非线性混合模型才能进行解混。为此, 本文提出了基于二次散射的非线性混合模型即二次散射模型。

光子在进入传感器前一般会在地物间发生多次碰撞和散射, 为了使计算简化, 本文忽略了3次以上的散射, 只考虑光子在进入传感器前与地物发生2次撞击的情况, 提出二次散射模型, 表达式为

y=c(Ma)+d(Ma)²+n, (2)

式中: M为端元光谱矩阵; a为端元丰度向量, 其中

0≤ a_i≤ 1, ∑i=1Ra_i=1; c+d=1(c≥ 0, d≥ 0), d为二次散射系数(表示光子照射到地物上发生二次散射的概率), c为与d有约束关系的系数(表示光子照射到地物上被直接反射的概率)。此模型概括了线性混合和二次散射混合情况, 当d=0时, 该模型为线性混合模型。

假设某一混合像元内只存在A和B 2种地物, 且地物的空间位置及比例都未知; 设2种地物的光谱向量分别为m₁和m₂, 丰度向量分别为a₁和a₂, 即M=[m₁, m₂], a=[a₁, a₂]^T。

光谱混合方式分线性和二次散射2大类, 如图1所示。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 地物散射情况Fig.1 Ground object scattering case

从图1可以看出, 当光子入射到某种地物上后, 将以概率c直接反射到传感器, 以概率d(d=1-c)散射到其他地物上。以概率d散射的这部分光子又分别以概率a₁和a₂碰撞到地物A或B后再反射到传感器中。该混合像元得到的光谱为

y=ca₁m₁+da₁ma₁m₁+da₁m₁a₂m₂+ca₂m₂+da₂m₂a₁m₁+da₂m₂a₂m₂+n。 (3)

式中: 前6项分别对应图1(a)— (f)6种情况; 第7项为误差项。则式(3)可进一步化简为

y=ca₁m₁+ca₂m₂+2da₁m₁a₂m₂+d(a₁m₁)²+d(a₂m₂)²+n=c(a₁m₁+a₂m₂)+d(a₁m₁+a₂m₂)²+n=c(Ma)+d(Ma)²+n。 (4)

光谱混合从形式上分为聚合式、整合式和致密式3种, 如图2所示。多次散射大多发生在致密式分布的地物中。

二次散射模型是半确定性模型, 需要首先确定二次散射系数d。二次散射系数能直观反映地物分布的复杂程度: 地物分布越简单, 二次散射系数越低(极限情况d=0, 为线性混合); 相反, 地物分布越复杂, 非线性程度越强, 则二次散射系数越高。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 光谱混合的3种形式Fig.2 3 forms of spectrum mixture

假设某混合区由2种纯地物组成, 如果某个像元只存在其中一种地物, 可认为二次散射系数d=0(线性情况)。假设像元内有2种地物且均等分布, 则地物的交界地带最多, 可认为此种情况下二次散射系数最大(图3(a)); 如果2种地物不均等分布, 比如地物A占10%, 地物B占90%, 则认为此种情况地物间存在交界, 但交界地带相对较少, 二次散射系数较低(图3 (b))。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 地物分布复杂度Fig.3 Complexity of ground object distribution

地物丰度方差能表征地物的复杂程度, 当像元中只存在一种纯地物, 则丰度方差最大, 对应的地物复杂度最低, 即二次散射系数最小(d=0); 当像元中地物分布均等时, 则地物丰度方差为0, 但地物分布最复杂, 二次散射系数最大。可见, 二次散射系数与地物丰度方差成反比。

如果忽略地物组分的颗粒度、粗糙度及透明度等微观因素, 则二次散射系数可定义为

d(i, j)=d_max- dmax·var(i, j)varmax, (5)

式中: d(i, j)∈ [0, d_max], d_max为最大的二次散射系数(本文设为0.5); var(i, j)为地物丰度方差, (i, j)为图像坐标; var_max为各混合像元中地物丰度方差的最大值。

2.1.1 模拟数据的生成

从ENVI USGS Spectral Libraries中找出天然碱、泻利盐和黝帘石3种地物光谱数据^[12]。该光谱数据的波长区间为395~2 560 nm, 光谱分辨率为10 nm, 波段数为210个。分别根据光谱线性混合模型和二次散射模型生成大小为100像元× 100像元且包含210个波段的高光谱图像, 并增加不同信噪比的高斯白噪声。3种纯地物的光谱曲线如图4所示。

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 3种纯地物的光谱曲线Fig.4 Spectral curves of 3 pure mineral

根据模拟数据生成的3种纯地物丰度图像如图5所示。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 模拟生成的3种纯地物丰度图Fig.5 Simulated abundance image of 3 minerals

将每种纯地物丰度图像变成一个1× 10 000的行向量, 然后将3个行向量组成一个3× 10 000的地物丰度矩阵。3个地物的光谱数据矩阵为210行× 3列, 每一列对应一种纯地物的光谱向量, 每种纯地物光谱向量有210个采样点。

首先, 对模拟数据进行光谱线型混合: 将光谱矩阵与地物丰度矩阵相乘, 得到210× 10 000的混合矩阵。该矩阵表示高光谱图像每1行为1个波段, 共210个波段。图6为第100和150波段的图像, 不同颜色代表不同反射率值。

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 线性混合模拟生成的高光谱图像Fig.6 Simulated hyperspectral image based on linear mixture

根据式(1)对模拟数据进行光谱非线性混合。图7是二次散射系数d=0.3时模拟得到的第100和150波段的图像。不同颜色代表不同反射率值。

图7

Fig.7

	Figure Option ViewDownloadNew Window
	图7 非线性混合模拟生成的高光谱图像Fig.7 Simulated hyperspectral image based on nonlinear mixture

在实验中, 对模拟得到的数据加入不同信噪比的噪声, 信噪比计算公式为

R_SNR=10 lg var(s)var(n), (6)

式中var(s)和var(n)分别代表信号和噪声方差。

图8为非线性光谱混合后再加入50 dB噪声后的第100波段图像。

图8

Fig.8

	Figure Option ViewDownloadNew Window
	图8 非线性光谱混合后加入50 dB噪声的第100波段图像Fig.8 100th band simulated hyperspectral image based on nonlinear mixture(RSNR=50 dB)

2.1.2 实验方法及结果

将模拟得到的数据分别用全约束最小二乘法(fully constraind least squares, FCLS)^[13]和非线性变换后的最小二乘法进行解混, 得到估计的地物丰度数据。利用式(7)计算估计值的均方误差^[14], 即

MSE=1N∑n=1Ny^n-yn2(7)

式中: y_n表示各纯地物的实际丰度值向量; y^n表示各纯地物的估计丰度值向量; N为采样点个数, 在本实验中N=60 000。

表1为在相同信噪比、不同二次散射系数情况下, 基于线性和非线性模型解混得到的均方误差。

表1

Tab.1

表1(Tab.1)

表1 相同信噪比(RSNR=50)、不同散射系数下的解混均方误差 Tab.1 MSE of solved abundance at same RSNR(RSNR=50) but different scatter coefficient 二次散 射系数d d=0 d=0.1 d=0.2 d=0.3 d=0.4 d=0.5 线性解混 0.001 0.060 0.117 0.170 0.214 0.245 非线性解混 0.001 0.001 0.001 0.001 0.001 0.001 提高率/% 0 98.0 99.1 99.4 99.5 99.6

表1 相同信噪比(RSNR=50)、不同散射系数下的解混均方误差 Tab.1 MSE of solved abundance at same RSNR(RSNR=50) but different scatter coefficient

由表1可以看出, 随着二次散射系数的增大, 即地物分布复杂情况的加重, 基于线性混合模型解混得到的均方误差值增大。可见, 地物分布越复杂, 线性混合模型的适用性越差; 而基于非线性混合模型的解混效果好得多。同时可以发现, 随着二次散射系数的提高, 非线性解混效果也变得更好。

表2为在相同二次散射系数(d=0.3)、不同信噪比情况下, 基于线性和非线性模型解混得到的均方误差。

表2

Tab.2

表2(Tab.2)

表2 相同二次散射系数(d=0.3)、不同信噪比下的解混均方误差 Tab.2 MSE of solved abundance at same d(d=0.3) but different RSNR RSNR 线性解混 非线性解混 提高百分比/% 10 0.205 0.137 33.4 20 0.187 0.034 81.8 30 0.175 0.011 93.7 40 0.171 0.004 98.0 50 0.170 0.001 99.4

表2 相同二次散射系数(d=0.3)、不同信噪比下的解混均方误差 Tab.2 MSE of solved abundance at same d(d=0.3) but different RSNR

由表2可以看出, 不管是基于线性混合模型还是基于非线性混合模型, 信噪比R_SNR越低, 解混效果越差。与线性解混结果相比, 非线性的解混效果更优, 提高的百分比随着信噪比的增大而增大, 表明非线性模型抗噪性更好。

图9为在二次散射系数d=0.3、信噪比R_SNR=40情况下, 基于线性混合模型解混得到的地物A丰度图像及其误差图像。

图9

Fig.9

	Figure Option ViewDownloadNew Window
	图9 基于线性混合模型解混得到的地物A丰度图像(左)及其误差图像(右)(d=0.3, RSNR=40)Fig.9 Solved abundance image(left)of mineral A and it’ s error image(right)based on linear mixture(d=0.3, RSNR=40)

图10为在二次散射系数d=0.3、信噪比为R_SNR=40情况下, 基于非线性混合模型解混得到的地物A丰度图像及其误差图像。

图10

Fig.10

	Figure Option ViewDownloadNew Window
	图10 基于非线性混合模型解混得到的地物A丰度图像(左)及其误差图像(右)(d=0.3, RSNR=40)Fig.10 Solved abundance image(left)of mineral A and it’ s error image(right)based on nonlinear mixture(d=0.3, RSNR=40)

基于图9, 10可以发现, 基于非线性混合模型的解混误差明显小于基于线性混合模型的解混误差。

2.2.1 数据的选取

选取美国内华达州南部Cuprite地区的AVIRIS高光谱数据^[15]。图像大小为250像元× 191像元, 包含224个波段, 波长区间为0.389~2.467μ m。除去低信噪比和水气吸收波段(波段1~2, 104~113, 148~167, 221~224), 最后可用数据为188个波段数据。图11为第100波段高光谱图像。

图11

Fig.11

	Figure Option ViewDownloadNew Window
	图11 Cuprite地区AVIRIS第100波段图像Fig.11 100th band image of AVIRIS of Cuprite district

Cuprite区主要分布6种矿物^[16], 分别为DesertVarnish(沙漠岩漆)、Alunite(明矾石)、Montmorillonite(蒙脱石)、Chalcedony(玉髓)、Kaolinite(高岭石)#1和Kaolinite(高岭石)#2。

2.2.2 实验方法及结果

利用顶点成分分析算法^[17]提取6种矿物端元, 基于线性混合模型、利用FCLS 方法^[13]得到的6种矿物丰度图像如图12所示。

图12-1

Fig.12-1

	Figure Option ViewDownloadNew Window
	图12-1 Cuprite地区6种矿物线性解混结果Fig.12-1 Solved abundance results of 6 minerals based on linear mixture in Cuprite district

图12-2

Fig.12-2

	Figure Option ViewDownloadNew Window
	图12-2 Cuprite地区6种矿物线性解混结果Fig.12-2 Solved abundance results of 6 minerals based on linear mixture in Cuprite district

利用式(5)确定基于该解混结果的地物复杂度, 设定最大二次散射系数d_max=0.5, 得到的地物复杂度分布如图13所示。

基于非线性混合模型(式(1)), 使用上面得到的二次散射系数对每个混合像元进行解混, 得到的解混结果如图14所示。

对比图12, 14可以看出, 基于二次散射模型的解混效果更好, 图像上道路影像特征清晰(图14); 而基于线性模型解混的道路影像特征则模糊(图12)。主要原因是道路与周围地物是主要边界地带, 而边界地带容易发生二次散射。

图13

Fig.13

	Figure Option ViewDownloadNew Window
	图13 Cuprite地区地物分布复杂度图Fig.13 Complexity of ground object distribution in Cuprite district

图14

Fig.14

	Figure Option ViewDownloadNew Window
	图14 Cuprite地区6种矿物非线性解混结果Fig.14 Solved abundance results of 6 minerals based on nonlinear mixture in Cuprite district

图15显示了在不同散射系数下(d=0~0.5), 蒙脱石的分布情况。

图15

Fig.15

	Figure Option ViewDownloadNew Window
	图15 不同二次散射系数下的蒙脱石分布情况Fig.15 Solved abundance of Montmorillonite in different secondary scatter coefficient

从图15可以发现, 基于二次散射模型解混的道路影像特征更加清晰连续, 而基于线性混合模型(d=0)的道路影像特征则比较模糊, 并出现道路影像局部中断现象。