基于微分搜索的高光谱图像非线性解混方法
【专利摘要】本发明属于高光谱图像处理技术领域，尤其涉及一种基于微分搜索的高光谱图像非线性解混方法，按照下述步骤进行，对实际高光谱图像进行端元提取；确定搜索个体的维数和位置编码；按照目标函数计算各搜索个体适应度值；计算经停位置StopoverSitei；比较每个搜索个体的当前位置Xi和其经停位置StopoverSitei的适应度值；根据条件决定是否进行计算。该方法将非线性解混问题转化为最优化问题，克服了传统梯度类优化算法对初始值要求高、易陷入局部收敛的局限性，较之线性解混算法和基于梯度优化的非线性解混算法，具有解混精度更高、稳定性更好的优点。
【专利说明】
基于微分搜索的高光谱图像非线性解混方法
技术领域
[0001] 本发明属于高光谱图像处理技术领域，尤其涉及一种基于微分搜索的高光谱图像 非线性解混方法。
【背景技术】
[0002] 高光谱图像解混是高光谱图像分析处理中的一项重要工作，其目的是为了解决高 光谱图像拍摄过程中空间分辨率不足的问题，已被广泛应用于遥感学、材料学和显微光谱 学等多学科领域。目前，大部分解混算法认为混合像元是端元光谱成分的线性混合。然而， 当高光谱图像场景中存在泥土、植被和水等地物时，光线会在不同地物之间发生反射和相 互作用，基于线性混合模型(Linear Mixing Model，LMM)的解混算法效果往往并不理想。应 该尝试在非线性混合模型(Nonlinear Mixing Model ,NLMM)下研究解混算法。
[0003] 目前，已有非线性混合模型主要基于辐射传输理论，但针对该模型的解混算法复 杂度很高，限制了高光谱图像解混技术的实际应用。为了避免高复杂度的数学计算，一些学 者在保证实际物理意义的基础上提出了简化而有效的非线性混合模型。其中，新近提出的 GBM是一个广义模型，它对未知混合类型的实际场景解混具有更广泛的适应性和更佳的性 能，开始成为国内外研究高光谱图像解混的重点模型。
[0004] 针对该模型解混的关键在于进行有效的丰度估计。但目前针对该类模型的丰度估 计方法主要是Bayesian估计和梯度法。然而，Bayesian估计方法存在参数复杂，计算量大等 缺点;基于梯度法的解混方法存在易陷入局部收敛、解混精度低等问题。这些问题使得该模 型下的高光谱图像解混效果仍然不尽如人意。
【发明内容】
[0005] 为了克服上述现有技术的不足，本发明的目的是提供一种基于微分搜索的高光谱 图像非线性解混方法，该方法将非线性解混问题转化为最优化问题，克服了传统梯度类优 化算法对初始值要求高、易陷入局部收敛的局限性，较之线性解混算法和基于梯度优化的 非线性解混算法，具有解混精度更高、稳定性更好的优点。
[0006] 为了实现上述目的，本发明采用如下技术方案：
[0007] -种基于微分搜索的高光谱图像非线性解混方法，其特征在于，按照下述步骤进 行：
[0008] 步骤1利用VCA算法对实际高光谱图像进行端元提取；
[0009] 步骤2根据高光谱图像中的端元数目R，按照高光谱图像GBM混合模型确定搜索个 体的维数和位置编码；
[0010] 步骤3按照微分搜索算法的初始化原理在搜索空间随机产生搜索个体，设定搜索 范围上限upj = l和搜索范围下限Iowj = O;
[0011] 步骤4构建目标函数，按照目标函数计算各搜索个体当前位置X1的适应度值；
[0012] 步骤5根据每个搜索个体的当前位置Xi计算经停位置StopoverSitei;再按照目标 函数计算经停位置StopoverSitei的适应度值；
[0013] 步骤6比较每个搜索个体的当前位置Xi和其经停位置StopoverSitei的适应度值， 如果经停位置StopoverSi tei的适应度值小于当前位置Xi，则由经停位置StopoverSitei代 替当前位置X1作为该搜索个体的当前位置，否则仍然保留搜索个体的当前位置X1;
[0014] 步骤7对当前位置X1中的丰度矢量(^,之,L 部分进行归一化处理；
[0015] 步骤8如果步骤7中归一化处理的结果已经达到最大进化代数，则输出当前搜索群体中
否则，转到步骤4进行下一代的进化搜索过程；
[0016] 步骤9如果已经对高光谱图像中的所有像素进行了解混，则停止计算；否则，返回 步骤3，针对下一个像素执行解混工作。
[0017] 所述步骤2中，高光谱图像GBM混合模型的表达式为：
[0018]
[0019]式中，y=[yi，y2, Λ，yL]T为每个像素点的观测值，L为波段数;R为端元数目，M为L \尺维的端元光谱矩阵，其每一列1111{=[1111,1{，1112,1<，八，1]11^] 1'(1^=1，2，八，1〇代表一种端元光 谱;a为实际图像的丰度矢量，且a = [ai，a2，L，ai，aj，L aR]T，ai，aj为像素点的丰度值，满足丰 度非负约束和丰度和为一约束条件；γ U是控制端元间相互作用程度的非线性参数，且〇< γ i,j<l ;η=[ηι，η2, Λ，η?为附加噪声。
[0020] 所述步骤4中，构建的目标函数为基于重构误差理论的目标函数，表达式为：
[0021]
[0022]式中，y为高光谱图像的观测值，Ibm为高光谱图像的重构值;丰度i为每个像素中 各端元成分所占比例值，丰度S中各元素满足非负特性约束(ANC)以及和为1约束(ASC);同 时，非线性参数矢量^满足γ i,j<i的要求。
[0023] 所述步骤5中，按照下式计算新的经停位置StopoverSitei，同时对经停位置 StopoverSitei进行边界控制：
[0024] StopoverSi tei = Xi+Rg(donor-Xi)
[0025] 式中，donor= [Xrancbm_shUffiing(i)]为从当前N个搜索个体中随机选取出的一个个体 所处的位置，利用(Cbnor-X 1)的交叉变异操作产生新的搜索方向；心为随机数，由gamma随机 数产生器和均匀分布随机数产生器共同运算得到;&为每个搜索个体的当前位置。
[0026] 随机数Rg的产生公式为:Rg = randg[2 · randi] · (rand2_rand3)〇
[0027] 本发明的有益效果是：（1)本发明是在新近提出的GBM(广义双线性模型）的基础上 研究的非线性解混算法，有效克服了线性解混算法在处理存在非线性效应高光谱图像时的 局限性。（2)本文发明采用群智能优化方法完成解混工作，有效克服了传统梯度类优化方法 对初始值要求高、易陷入局部收敛的缺点，具有更高的解混精度和稳定性。（3)本文发明物 理意义明确、算法研究过程无须进行复杂的公式推导。因此，可以利用本发明的算法框架针 对其它混合模型进行推广而得到新的解混算法。
【附图说明】
[0028]图1为本发明方法的解混流程图；
[0029]图2为真实场景的Moffett Field高光谱图像；
[0030]图3为真实场景的Japser Ridge高光谱图像；
[0031]图4为从Moffett Field高光谱图像中提取的端元光谱曲线；其中图4(a)为从 Moffett Field高光谱图像中提取的vegetation端元光谱曲线；图4(b)为从Moffett Field 高光谱图像中提取的water端元光谱曲线；图4(c)为从Moffett Field高光谱图像中提取的 soil端元光谱曲线；
[0032]图5为从Japser Ridge高光谱图像中提取的端元光谱曲线；其中图5(a)为从 Japser Ridge高光谱图像中提取的tree端元光谱曲线；5(b)为从Japser Ridge高光谱图像 中提取的water端元光谱曲线；图5(c)为从Japser Ridge高光谱图像中提取的soil端元光 谱曲线；图5(d)为从Japser Ridge高光谱图像中提取的road端元光谱曲线；
[0033]图6为使用本发明方法得到的Moffett Field高光谱图像解混丰度图；其中，图6 (a)为使用本发明方法得到的vegetation端元的Moffett Field高光谱图像解混丰度图；图 6(b)为使用本发明方法得到的water端元的Moffett Field高光谱图像解混丰度图；图6(c) 为使用本发明方法得到的soil端元的Moffett Field高光谱图像解混丰度图；
[0034]图7为使用本发明方法得到的Japser Ridge高光谱图像解混丰度图，其中图7(a) 为使用本发明方法得到的tree端元的Japser Ridge高光谱图像解混丰度图；图7(b)为使用 本发明方法得到的water端元的Japser Ridge高光谱图像解混丰度图；图7(c)为使用本发 明方法得到的soil端元的Japser Ridge高光谱图像解混丰度图；图7(d)为使用本发明方法 得到的road端元的Japser Ridge高光谱图像解混丰度图。
【具体实施方式】
[0035]下面结合附图对本发明的【具体实施方式】做出说明。
[0036] 结合图1、图2和图3对本发明的高光谱图像非线性解混的具体步骤如下(解混的流 程图如图1所不）：
[0037] 步骤1输入高光谱仪采集到的实际高光谱图像，利用VCA算法针对图像进行端元提 取。
[0038]真实场景数据选用Moffett Field数据（如图2所示）和Japser Ridge数据（如图3 所示)两种高光谱图像数据。
[0039] (I)Moffett Field数据
[0040] 该数据源自美国加州San Francisco Bay南端的Moffett Field高光谱图像。该图 像具有189个波段，波长范围为400腦-2500腦，光谱分辨率为10腦。实验选取50\50的子图 像进行算法解混性能评估，该场景主要由vegetation、water和soil三种端元成分组成。
[0041] (2)Japser Ridge数据
[0042]该数据源自美国加州地区的Japser Ridge自然保护区的高光谱图像。该图像具有 224个波段，波长范围为380nm-2500nm，光谱分辨率为9.46nm。在去除了 1-3、108-112、154-166和220-224共26个波段的数据后(为去除水汽和大气影响），剩余198个有效波段数据。实 验选取100 X 100的子图像进行算法解混性能的评估，该场景主要由Tree、Water、Soi 1和 Road四种端元成分组成。
[0043] 步骤1利用VCA算法对实际高光谱图像进行端元提取；
[0044] 步骤2根据高光谱图像中的端元数目R，按照高光谱图像GBM混合模型确定搜索个 体的维数和位置编码；
[0045] 其中所述步骤2中，高光谱图像GBM混合模型的表达式为：
[0046]
[0047]式中，y=[yi，y2, Λ，yL]T为每个像素点的观测值，L为波段数;R为端元数目，M为L \尺维的端元光谱矩阵，其每一列1111{=[1111,1{，1112,1<，八，1]11^] 1'(1^=1，2，八，1〇代表一种端元光 谱;a为实际图像的丰度矢量，且a = [ai，a2，L，ai，aj，L aR]T，ai，aj为像素点的丰度值，满足丰 度非负约束和丰度和为一约束条件；γ U是控制端元间相互作用程度的非线性参数，且〇< γ i,j<l ;η=[ηι，η2, Λ，η?为附加噪声。
[0048] 步骤3按照微分搜索(Differential Search,DS)算法的初始化原理在搜索空间随 机产生搜索个体，设定搜索范围上限Upj = I和搜索范围下限Iowj = O;
[0049] DS算法及其初始化原理为现有技术，具体实现按照论文中进行，该论文为 Civicioglu P. Transforming geocentric cartesian coordinates to geodetic coordinates by using differential search algorithm[J].Computers and Geosciences，2012，46:229-247·
[0050] 步骤4构建目标函数，按照目标函数计算各搜索个体当前位置X1的适应度值；目标 函数如下式所示：
[0051]
[0052] 式中，y为高光谱图像的观测值，i，c,BM为高光谱图像的重构值;丰度0为每个像素中 各端元成分所占比例值，求解过程必须保证丰度6中各元素的非负特性约束(ANC)以及和为 1约束(ASC);同时，非线性参数矢量#还要满足0< γ i,j<l的要求。
[0053] 步骤5根据每个搜索个体的当前位置Xi计算经停位置StopoverSitei;再按照步骤4 的目标函数计算经停位置StopoverSitei的适应度值;经停位置StopoverSitei的计算式如 下，同时对经停位置StopoverSitei进行边界控制。
[0054] StopoverSitei = Xi+Rg(donor-Xi)
[0055] 式中，donor= [Xrancbm_shUffiing(i)]为从当前N个搜索个体中随机选取出的一个个体 所处的位置，进而利用(donor-Xi)的交叉变异操作产生新的搜索方向。R g为随机数，由ga_a 随机数产生器和均勾分布随机数产生器共同运算得到，Rg = randg[2 · randi] · (rand2_ rand3)，其作用是保证各搜索个体向着新产生的方向（Cbnor-X1)以随机的跨度探索。
[0056] 步骤6比较每个搜索个体的当前位置Xi和其经停位置StopoverSitei的适应度值， 如果StopoverSi tei的适应度值小于当前位置Xi，则由经停位置StopoverSitei代替当前位 置X1作为该搜索个体的当前位置，否则仍然保留搜索个体的当前位置X 1。
[0057]步骤7对当前位置X1中的丰度矢量部分进行归一化处理。
[0058]步骤8如果上述归一化处理的结果已经达到最大进化代数，则输出当前搜索群体中最优 搜索个体的位置，4,；Pu，fu,L，;?li?，;> 2,3，;?24,L ,;P2i?,L 从而得到丰度 矢量4 =[名為1為]τ和非线性参数矢量H八2，八."[，L,^，;^,4，L ,心，L A 否则，转到步骤4进行下一代的进化搜索过程。
[0059] 步骤9如果已经对高光谱图像中的所有像素进行了解混，则停止计算；否则，返回 步骤3，针对下一个像素执行解混工作。
[0060] 使用本发明的解混方法进行实验，选用：（1)美国加州San Francisco Bay南端的 Moffett Field高光谱图像，如图2所示。该图像具有189个波段，波长范围为400nm-2500nm， 光谱分辨率为l〇nm。实验选取50X50的子图像。（2)选用美国加州地区的Japser Ridge自然 保护区的高光谱图像如图3所示。该图像具有224个波段，波长范围为380nm-2500nm，光谱分 辨率为9.46nm。在去除了 1-3、108-112、154-166和220-224共26个波段的数据后(为去除水 汽和大气影响），剩余198个有效波段数据。实验选取100 X 100的子图像。
[0061] 为了实现全自动解混，端元提取方法选用VCA算法。图4(a)、图4(b)和图4(c)所示 为利用VCA算法分别从Moffett Field高光谱图像数据中提取出vegetation、water和soil 三种端元的光谱曲线。图5(a)、图5(b)、图5(c)和图5(d)所示为利用VCA算法分别从Japser Ridge高光谱图像数据中提取出Tree、Water、Soil和Road四种端元的光谱曲线。
[0062]图6(a)、图6(b)和图6(c)为采用本发明解混Moffett FieId高光谱图像得到的 vegetation、water和soil的丰度图。图7(a)、图7(b)、图7(c)和图7(d)为采用本发明解混 Japser Ridge高光谱图像得到的Tree、Water、Soil和Road丰度图。
[0063] 通过对比分析可知，本发明的解混方法解混出的丰度与真实地物成分的分布是一 致的。
[0064] 为了定量评价本发明对真实高光谱遥感数据的解混性能，采用RE和SAM两个指标 本发明的效果进行评价。
[0065] (I)RE
[0066]
[0067] 式中，L为波段数目，I X J为高光谱图像像素数。y(p)和少(/?)分别为高光谱图像第 P个像素点的真实观测数据和本发明算法重构出的数据。
[0068] (2)SAM
[0069]
[0070]
[0071] 将本发明方法与Bays、FCLS和GBM-GDA方法进行比较，得到的RE、SAM列在表1中，从 表1中可以看出，本发明方法的性能优于Bays和FCLS两种线性解混方法，并且还优于基于梯 度优化的非线性解混方法GBM-GDA，具有最小的RE和SAM值。
[0072]表1不同方法真实数据解混性能 l·
l·
[0075]以上对本发明的一个实例进行了详细说明，但所述内容仅为本发明的较佳实施 例，不能被认为用于限定本发明的实施范围。凡依本发明申请范围所作的均等变化与改进 等，均应仍归属于本发明的专利涵盖范围之内。
【主权项】
1. 一种基于微分搜索的高光谱图像非线性解混方法，其特征在于，按照下述步骤进行： 步骤1利用VCA算法对实际高光谱图像进行端元提取； 步骤2根据高光谱图像中的端元数目R，按照高光谱图像GBM混合模型确定搜索个体的 维数和位置编码； 步骤3按照微分搜索算法在搜索空间随机产生搜索个体，设定搜索范围上限uW = l和搜 索范围下限l〇wj = 0; 步骤4构建目标函数，按照目标函数计算各搜索个体当前位置t的适应度值； 步骤5根据每个搜索个体的当前位置Xi计算经停位置StopoverSitei;再按照目标函数 计算经停位置StopoverSitei的适应度值； 步骤6比较每个搜索个体的当前位置Xi和其经停位置StopoverSitei的适应度值，如果 经停位置StopoverSi tei的适应度值小于当前位置Xi，贝lj由经停位置StopoverSitei代替当 前位置t作为该搜索个体的当前位置，否则仍然保留搜索个体的当前位置χ 1; 步骤7对当前位置&中的丰度矢量,??λ))部分进行归一化处理； 步骤8如果步骤7中归一化处理的结果已经达到最大进化代数，则输出当前搜索群体中最优 搜索个体的位置:(4,毛，L，毛,八2，八3，L 六,4，L ,?^,L从而得到丰度矢 量6 = 毛，l 和非线性参数矢量:ΗΑ,；!，Λ,3，ι;,六』，U2, 4，l，u; 否则，转到步骤4进行下一代的进化搜索过程； 步骤9如果已经对高光谱图像中的所有像素进行了解混，则停止计算；否则，返回步骤 3，针对下一个像素执行解混工作。2. 如权利要求1所述的基于微分搜索的高光谱图像非线性解混方法，其特征在于:所述 步骤2中，高光谱图像GBM混合模型的表达式为：式中，y = [yi，y2, Λ，yiJT为每个像素点的观测值，L为波段数;R为端元数目，Μ为LXR维 的端元光谱矩阵，其每一列mk=[mi,k，m2,k，λ，mL,k]T(k=l，2, Λ，R)代表一种端元光谱;a为 实际图像的丰度矢量，且a= [ai，a2，L，ai，aj，L aR]T，ai，aj为像素点的丰度值，满足丰度非负 约束和丰度和为一约束条件；γ U是控制端元间相互作用程度的非线性参数，且〇< γ 1 ;η=[ηι，Π2, Λ，η?为附加噪声。3. 如权利要求1所述的基于微分搜索的高光谱图像非线性解混方法，其特征在于:所述 步骤4中，构建的目标函数为基于重构误差理论的目标函数，表达式为：式中，y为每个像素点的观测值，为高光谱图像的重构值;丰度?为每个像素中各端 元成分所占比例值，丰度?中各元素满足非负特性约束以及和为1约束；同时，非线性参数矢 量歹满足〇< γ i,j<l的要求。4. 如权利要求1所述的基于微分搜索的高光谱图像非线性解混方法，其特征在于:所述 步骤5中，按照下式计算新的经停位置StopoverSitei，同时对经停位置StopoverSi tei进行 边界控制： StopoverSitei = Xi+Rg(donor-Xi) 式中，d〇n〇r= [Xrancbm_shUffling(i)]为从当前N个搜索个体中随机选取出的一个个体所处 的位置，利用(donor-Xi)的交叉变异操作产生新的搜索方向;Rg为随机数，由gamma随机数产 生器和均匀分布随机数产生器共同运算得到;&为每个搜索个体的当前位置。5.如权利要求4所述的基于微分搜索的高光谱图像非线性解混方法，其特征在于:随机 数Rg的产生公式为：Rg = randg[2 · randi] · (rancb-rancb)。
【文档编号】G06T1/00GK106056524SQ201610352111
【公开日】2016年10月26日
【申请日】2016年5月25日
【发明人】陈雷, 张立毅, 费腾, 张勇, 孙云山
【申请人】天津商业大学