背景技术
遥感技术(Remote Sensing,RS),指的是开始于二十世纪六十年代,从远距离空间(航天)或者外太空空间(航空)对地球表面进行的综合性的观测技术。遥感,即从远距离或外太空间接观测监控,不接触目标,从光学角度获取对目标、现象以及区域的相关信息,从而进行数据的融合分析以及推断,最终达到获取所需的目标信息的一种手段、技术和科学。高光谱图像遥感技术(Hyper-spectral Graphic Remote Sensing),是一种立于高光谱技术的技术基础上的遥感探测技术,是一种融合新型光谱探测技术、微信号探测技术、光学精密机械、信号高速处理技术、计算机处理信息技术在一起的先进性、综合性、科学性技术。该技术同时由于影响到地球学、环境保护学、野生生物学、信息技术学、空间地理学等众多科学领域的关键技术,其技术的进步发展备受国内外学者的强烈瞩目,目前已广泛应用于植被绿化、土壤分析、精细农业管理、大气环境监测、水环境监测、勘探矿产资源分布等方面,充分展示出了高光谱遥感技术的潜力与优势。
在全球遥感界中通常会这样认为:多光谱遥感(Multi-spectral RemoteSensing)指的是光谱辨析率在1~9.9×10-1λ的光谱范围内的遥感,在这种范围内的遥感传感器内只有分布在可见光区和近中红外光光谱区的很少的波段数;高光谱遥感(Hyper-spectral Remote Sensing)指的是光谱辨析率在1~9.9×10-2λ的光谱范围内的遥感技术;而超高光谱遥感(Ultra-spectral Remote Sensing)则指的是经过更高技术的进步后,光谱辨析率达到1~9.9×10-1λ的光谱范围内的遥感技术。遥感技术,其发展历程是在全色(即黑白色)摄影与彩色影像阶段之后,多光谱遥感迅速发展于二十世纪中后半叶时期,取得了较大的进步,并已运用于环境与资源探测领域,但其分辨率仍然处于10-1λ数量级光谱的波长范围内,只有相对较少的采样点,随着科技的发展,人们对地球资源与环境认识越来越深入,其运用精度越来越不能满足需求,对分辨率更深入的发展要求更加迫切,主要体现在空间分辨率和光谱分辨率的深入与提高。图片中在空间分辨率指的是在遥感图片中相邻的两个地形地物之间能够被区分出的最小长度,也就是通常所说的图片的清楚程度,可用于对图片数据的直观解释。而图片中在光谱分辨率指的是光谱检测中的光谱反射出的能够区分出的最小光谱间隔的波长,更确切的讲,就是检测区分光谱的性能。相对于空间分辨率的直观观察,光谱观测技术能够更加直接和有效的反应地物的结构与性质,尤其是在远洋探测监测、动植物研究分类、农业精细化规划、勘探矿产调查和现代化军事运用方面具有更好的效果和更快的速率。随着科技的日新月异,航天航空工业的飞速发展,遥感平台和光学传感器的不断更新进步升级,提高遥感光谱分辨率的已经处在迅猛的发展势头中,高光谱遥感技术的发展已经成为当今遥感技术研究学者们的研究热点。由于种类性质不一样的物质在一定的波长的光谱照射下,具有各自不同比例的吸收和反射的特性,通过对比其反射光谱(或吸收)光谱的之间的差异,可以导出其物质的组成成分和物理结构上之间的差异。如何才能从数据复杂的图像中提取出物质的特有属性(特征)并且清晰了解与周边物质及在整体形态中的关系,成为研究中主要解决的问题。
高光谱遥感影像(Hyper-spectral Remote Sensing Images),是指在电磁波谱中的中远红外光谱、近红外光谱、可见光光谱和紫外光谱区域内,在光谱成像仪的作用下,得到的很多光谱分布连续且光谱区域范围很小的影像信息数据(如图1所示)。随着成像光谱仪的迅速发展与进步,所能获取到的高光谱遥感影像相比于传统的二维成像技术的遥感图像,其是具有几十以至上百个波段叠加而成的光谱图,其中的每一个基本的像元结构都是从几十以至上百个连续波段通道中所获取的得到的光谱图像,其正对应的光谱实物的反射特性,最后都会得到一条较为完整的光谱曲线。高光谱遥感图像已不再是二维成像技术,将独特的光谱理论与遥感成像技术有机的结合起来,形成的连续的地物光谱曲线,使得利用高光谱技术能够成功反演地形地物的细节。
高光谱遥感图像技术是一种新型的综合性对地观测技术,该技术的发展已在众多领域得到广泛的运用。高光谱图像由于其光谱分辨率高,在成像中能够记录几十以至几百个波段的光谱信息,能够很好的识别分析地形地物的种类与性质,其良好的效果与广泛的用途使得高光谱技术得到了众多学者们的格外关注,引起了对高光谱图像技术的研究热潮。
相对于纯像元而言,混合像元是高光谱遥感图像中的每个数据单元是由不同地物的光谱信息组成的混合信息。在现实中,由于地形地物分布不均、交综复杂以及成像光谱仪的分辨率有限,导致高光谱遥感图像中一定会大量的存在着不均匀像元(即混合像元),尤其是在地物种类繁多且分布复杂的地形地域,若是将图片中像元看成是纯净地物,那么会造成较大的误差和分析的精度下降。这使得在处理分析过程中对混合像元的研究成为高光谱遥感图像的主流研究方向,成为一个新的研究热点。混合像元的存在阻碍了对高光谱图像的研究利用,对混合像元的分解成为了高光谱图像技术的基础。
1999年,D.D.Lee与H.S.Seung在美国著名科学杂志《Nature》上发表了一篇极其重要的文章,其创造性地提出了一种带有非负性约束的矩阵分解的方法:非负矩阵分解算法(Non-negative Matrix Factorization,NMF)。非负矩阵分解,是一种在非负约束下的矩阵分解算法,由于其算法的便捷性和对广泛物理意义的表征性使得非负矩阵分解算法得到广泛的运用。该算法是在非负性的约束下,将原矩阵分解成一个基矩阵和系数矩阵的乘积,其中基矩阵和系数矩阵中所有元素满足非负的约束。光谱信息的非负性也适合于非负矩阵的要求,分解后的结果是端元向量矩阵与其丰度矩阵的乘积,在物理意义的吻合使得非负矩阵分解算法在混合像元分解的研究中得到了广泛的运用。但是非负矩阵分解算法具有迭代收敛慢和容易陷入局部最优解的局限性。
2002年P.O.Hoyer在IEEE(Institute of Electrical and ElectronicsEngineers)上发表了重要文章,提出了一种将稀疏编码与NMF算法相融合的方法,该方法是在NMF算法加入了一种稀疏度约束。W.Liu等人在此基础上提出了稀疏NMF(Sparse NMF,SNMF)算法。SNMF算法在NMF的目标函数中加入了关于表征丰度稀疏性的约束项,而对于自然信号,在科学与工程领域中而言,稀疏性是一种比较常见的属性。
加权NMF(Weighted NMF,WNMF)算法是众多的改进算法中的一类。它是在NMF目标函数中加入了对角矩阵U,该矩阵中对角线上的元素值大小表示着各个训练样本向量的重要程度,从而对各类样本的重要程度进行了区分,这就会有以利于防止出现相似的训练样本向量会产生冗余信息。加权式的改进措施可以适用于任何一个优秀的方法中。
在NMF分解中,产生的基矩阵W通常具有显示训练样本特征的物理意义,而这特征也就是训练样本的局部性,进而局部NMF(Local NMF,LNMF)算法就应运而生,该算法是在目标函数里加入了一定的约束因子,如此确保分解结果尽可能的正交化,从而使得分解产生的结果能够更多的体现出训练样本的局部特征。
多年来,在基本NMF算法的基础上,很多改进的NMF算法被提出,SNMF、WNMF、LNMF等,不同的改进算法有着不同的运用方向,经实际运用证明,各改进算法均具有更优异的不同的属性和更加令人满意的收敛速度,对于NMF算法的实际运用取得了很大的增强效果。
发明内容
发明目的:在高光谱遥感图像的研究中,由于地形地物的复杂性,即使较高的分辨率也不能使得图像的像元为纯净的像元,而是由像元空间区域内的各种端元的混合反射光谱,即混合像元。混合像元的存在阻碍了对高光谱图像的研究利用,非负矩阵分解是一种在非负约束下的有效的矩阵分解算法。为了利于高光谱遥感图像中混合像元分解,使得分解所得的基矩阵和系数矩阵能够得到期望的高稀疏性和低相关性的要求,对非负矩阵分解过程的目标函数中的基矩阵和系数矩阵施加稀疏性约束并对系数矩阵施加相关性约束。
技术方案:一种应用于高光谱图像处理的非负矩阵分解方法,该方法考虑矩阵的稀疏特性和相关性,并一定程度上影响了混合像元分解的精度。具体如下:此方法最终将把给定的非负矩阵Vm×n分解成基矩阵Wm×r和系数矩阵Hr×n的乘积,即Vm×n≈Wm×rHr×n,首先选择非负矩阵V,随机初始化W和H,然后在目标函数中对系数矩阵H施加最小相关性约束,对基矩阵W和系数矩阵H施加稀疏性约束,然后按照迭代公式进行迭代至矩阵W和H收敛。
具体包括如下步骤:
步骤1:选择非负矩阵V,随机初始化W和H;
步骤2:通过对目标函数F添加了对系数矩阵H进行了最小相关性约束和对基矩阵W和系数矩阵H都施加了2-范数约束,得到式子(5),利用基本NMF算法中的下降法给出一种更新规则,得到新的迭代公式(6),按式子(6)进行迭代;
步骤3:在每一步迭代后将W和H的负元素置零,并对W按列进行归一化处理;
步骤4:循环进行步骤二~三直至收敛,此时得到基矩阵W和系数矩阵H。
为了更好的理解本发明所涉及的技术和方法,在此对本发明涉及的理论进行介绍。
1、基本NMF算法
假设给定原矩阵Vm×n满足非负,则NMF算法的目的是求解出两个非负低秩子矩阵Wm×r和Hr×n,从而使得
Vm×n≈Wm×rHr×n (7)
其中Wm×r称为基矩阵或特征矩阵,Hr×n称为系数矩阵。一般情况下r<<min(m,n),从而达到降低数据存储维数的目的,同时使得子矩阵Wm×r和Hr×n具有更加明确的解释意义。
NMF算法是一种低秩逼近算法,通常其目标函数可以是多种多样的,D.D.Lee与H.S.Seung在文章中主要对以下两种目标函数进行了研究推导:
一是运用原矩阵V和基矩阵W、系数矩阵H之间的欧氏距离(Euclidian Distance),将其平方作为目标函数,推导出迭代公式,然后通过逐步的迭代,尽量使得欧氏距离最小化从而达到最佳逼近的目的,其目标函数为:
当F(W,H)=0时,V=WH。目标函数式(8)关于W和H分别都是凸的,但对W和H联合起来时并非是凸的,即同时关于W和H的时候目标函数是非凸的,若在全局范围内求解最小值则是非常困难的,所以NMF求出的结果只能是局部最优解。利用梯度下降法给出了一种更新规则,首先根据目标函数式(8)对W和H分别求梯度得到:
再通过逐步迭代,其迭代规则为:
上式中当学习速率ηW、ηH按下式进行取值时:
将式(9)和式(11)分别代入式(10),从而可以得到一个乘性迭代的规则:
将其写成矩阵的形式为:
二是运用V和WH之间的散度(KL-散度),也就是概率意义上的距离作为目标函数来达到最佳逼近:
该解法如同欧氏距离,可得到以下乘性迭代规则
NMF算法是一种低秩迭代算法,其实现的基础是目标函数的收敛性,下面进行收敛性证明:
首先引入相关的定义,假定目标函数是F,则设G(h,ht)为F(h)的辅助函数,其满足条件:G(h,ht)≥F(h),当h=ht时相等,即G(ht,ht)=F(ht),其中h为一变量。G(h,ht)与F(h)的关系如图2所示。
引理1:假定G为辅助函数,则在以下的迭代规则,F是非增的:
F(hmin)≤…≤F(ht+1)≤F(ht)≤…≤F(h1)≤F(h0) (17)
对于利用目标函数对矩阵W和H的迭代过程是否收敛,我们只要能够构造出上面所说的辅助函数,即可证明迭代过程收敛,下面以矩阵H为例,W的证明与H是一样的,不再赘述。
引理2:假定K(ht)为一个对角矩阵,Kab(ht)=δ(WTWht)a/hta,其中δ函数是一个冲激函数,则
其中F(h)为目标函数对矩阵H的迭代,即要证明该G(h,ht)函数为F(h)的构造函数只需证明G(h,ht)≥F(h)与G(ht,ht)=F(ht)即可。显然,把h=ht带入式(18)中G(ht,ht)=F(ht),下面证明G(h,ht)≥F(h):
对F(h)在h=ht处进行泰勒公式展开即
要证明G(h,ht)≥F(h),只需证明G(h,ht)-F(h)≥0,即
(h-ht)T(K(ht)-WTW)(h-ht)≥0 (20)
为了证明其非负性,构造矩阵Mij(ht)=hti(K(ht)-WTW)htj,只要证明Mij(ht)是半正定的即可证明式(20),即只要证明vTMv≥0,即
所以,NMF算法的主要步骤和要求为:
(1)初始化W和H;
(2)选择学习速率和参数,保证迭代过程始终非负;
(3)按照迭代公式进行迭代;
(4)将W和H强制非负,且将H进行归一化处理;
(5)反复进行3、4至收敛。
2、稀疏NMF算法
2002年P.O.Hoyer在IEEE(Institute of Electrical and ElectronicsEngineers)上发表了重要文章,提出了一种将稀疏编码与NMF算法相融合的方法,该方法是在NMF算法加入了一种稀疏度约束。W.Liu等人在此基础上提出了稀疏NMF(Sparse NMF,SNMF)算法。所谓稀疏性,指的是矩阵中的元素数据分布的一种特性,具有稀疏性的矩阵叫做稀疏矩阵,是指其矩阵中的元素中存在大量的零元素且占有比例很高的矩阵。对于自然信号,在科学与工程领域中而言,稀疏性是一种比较常见的属性。
SNMF算法是在NMF的目标函数中加入了关于表征丰度稀疏性的约束项,从而使得目标函数变为
约束项是以系数矩阵H中所有元素的代数和来表征对其稀疏度的约束,并在迭代中使该约束项的数值达到最小,以使得结果尽可能的接近原矩阵,λ称为稀疏因子。
3、稀疏性约束
一般情况下,高光谱遥感图像中的各个像元中的端元分布不可能在整张图像中分布,一般会比较集中分布于一片或若干片较为集中的区域。所以,在对高光谱遥感图像非负分解的过程中施加稀疏性的约束,是符合实际情况,并且有利于图像处理的后续工作的进行。
稀疏度,在矩阵论中对向量来说有一种利用1-范数与2-范数之间的关系作为度量的方式,其具体方式为:
其中n为向量x的维度,||x||1=Σ|xi|为向量x的1-范数,为向量x的2-范数。在对范数的学习中有这样的关系式:
所以可以得到0≤S(x)≤1,当且仅当向量x只有一个元素时S(x)=1,且当向量x中所有元素相等且不为零时S(x)=0。
可以看出,我们想得到较大稀疏度的条件是使得1-范数与2-范数之比较大,即该式的值较小,因为在NMF算法的迭代过程中,对迭代的结果都要作归一化处理,即1-范数的值等于1(||x||1=1),所以我们要最大化2-范数,可以使得分解后的结果具有较大的稀疏度。
在NMF算法中,在添加稀疏约束时不能仅仅是简单的从向量变成矩阵的2-范数,矩阵的2-范数的定义为若仅仅是最大化基矩阵W的2-范数,则有可能导致产生分解结果中W中有一个列向量的2-范数很大,但其他的却很小,从而导致结果不能达到期望的稀疏性的要求,所以在施加约束是考虑将最大化基矩阵W所有列向量的2-范数的平方之和作为目标,即
同理对系数矩阵H施加为:
4、相关性约束
相关性,指的是度量两个或两个以上的具备相互之间有所关联的性质的元素变量之间的相互关联的大小指数。在高光谱图像中,各个像元的特征曲线之前的相关程度主要存在于地物之间的差异,高光谱由于其光谱分辨率比较高,因而其相关性应该比较小,但是由于地物的复杂性,混合像元在去光谱特征具有一定的相似性,表现为具有一点的相关性,而各个端元之间的相关性是很小的。
在数理统计理论中,相关系数是指表现两个或两个以上的信号之间的相关性的大小。相关系数越大表示相关性就越大,且在不相关的一组信号数据间所表现出的相关系数的绝对值必定会小于它们之间的混合信号数据之间所表现出的相关系数的绝对值。在NMF算法中,我们希望得到的结果是具有较小的相关性,若在分离时得到的信号之间的相关系数的绝对值越小就表明分离出来的端元之间通常都会包含一定量的相似的物质,在光谱特征上就变现为具有一定的相关性。这种相关性通常会与像元之间的平滑性较为接近,表现为混合像元的光谱号间的相关性越小,因此可以把分离信号相关系数最小作为对目标函数的迭代过程中的一个约束条件。在信息学中,相关系数的定义为:
当i≠j,若<hi,hj>=0,则R(H)取最小值。
具体实施方式
下面结合具体实施例,进一步阐明本发明,应理解这些实施例仅用于说明本发明而不用于限制本发明的范围,在阅读了本发明之后,本领域技术人员对本发明的各种等价形式的修改均落于本申请所附权利要求所限定的范围。
应用于高光谱图像处理的非负矩阵分解方法,包括如下步骤:
步骤1:选择非负矩阵V,随机初始化W和H;
步骤2:按式子(6)进行迭代;
步骤3:在每一步迭代后将W和H的负元素置零,并对W按列进行归一化处理;
步骤4:循环进行步骤二~三直至收敛,此时得到基矩阵W和系数矩阵H。
仿真实验结果分析
1.实验图像
通过仿真实验对算法的性能进行分析和评价。为了验证基于约束的NMF算法的有效性,本文对针对某一标准高光谱遥感影像进行试验。所使用的遥感影像采用USGS数据库图像。
以下4张高光谱遥感图像是从USGS官网上下载的南京地区周围的4个区域的图像。图像格式大小为420×314像素,因为高光谱遥感图像中的像元都混合的端元组成的复杂光谱曲线,需要经过预处理,本实验采用USGS官网中可供下载的已处理过的图片,可直接用于下一步处理。实际图片如图(3)所示。
USGS波普数据库其中主要是分析立体三维的成像光谱遥感数据,其包含了多达498条,其共有有针对性地反应各种自然矿物质(即矿物矿产)资源的444个植物或其他示例样品的共计498个反射波谱段数据。数据库的光谱波段范围为0.2~3.0μm,波谱经过进一步的校正后,其实际波谱范围为0.4~2.5μm,其分辨率在4nm(其波谱在可见光波段0.2~0.8μm)和10nm(其波谱在近红外波段0.8~2.35μm)。
2、实验预处理和参数选择
下载的图片格式为彩色图片,大小为420×314×3,对于计算起来有些不方便,所以首先转化成灰阶图片,其格式大小为420×314,像素值大小为0~1,为double型。转化成的图像如下图4所示。
对该4张图片使用本算法进行试验,对参数的选择主要有降维后的基矩阵W的列数和系数矩阵H的行数r、迭代算法的学习速率ηW和ηH以及世界约束条件时的正则化参数αW和αH。经过多次试验对比,为了使得试验效率更高和结果更好,本实验的参数选择为r=140,学习速率取式(11),正则化参数取αW=0.3,αH=0.5,βH=0.2,迭代次数为1000。
3、实验结果与分析
实验结果效果图如图5~8所示。其中V是原图像,WH是分解后重构的图像,为了作以对比V-WH是指的重构前后的误差值,下面的W和H是分解得到的基矩阵和系数矩阵。在这里我们主要关心的是误差值与分解的结果的优良性,但是从图中我们无法辨别。为此,在试验中除了采用矩阵的稀疏度作为对比的标准外,还将引入一个表征数据准确度的标准:信噪比(Signal-to-Noise Ratio,SNR)。信噪比,又名“讯噪比”,指的是在电子领域中,某个设备或系统中的原信号与在该设备或系统运行过程中所产生的噪声的比例大小。在图像处理的领域中,信噪比有如下定义(带入本实验的变量):
其中V为原图像的数据矩阵,WH为分解后重构的图像的数据矩阵,信噪比SNR的数据量的单位采用分贝(dB)。在试验中我们将本文算法与NMF和SNMF两种算法进行对比,结果如下:
表.1 南京禄口机场附近区域
算法 |
H的稀疏度 |
信噪比 |
NMF |
0.2648 |
23.6317 |
SNMF |
0.6184 |
25.6542 |
本文算法 |
0.6257 |
31.5333 |
表.2 南京青龙山附近区域
算法 |
H的稀疏度 |
信噪比 |
NMF |
0.2896 |
23.4918 |
SNMF |
0.6254 |
25.2826 |
本文算法 |
0.6458 |
31.1413 |
表.3 南京江心洲附近区域
算法 |
H的稀疏度 |
信噪比 |
NMF |
0.3152 |
23.0784 |
SNMF |
0.7185 |
25.3059 |
本文算法 |
0.7148 |
31.1424 |
表.4 南京栖霞区附近区域
算法 |
H的稀疏度 |
信噪比 |
NMF |
0.2800 |
23.6474 |
SNMF |
0.6605 |
25.2716 |
本文算法 |
0.6572 |
31.0521 |
通过以上实验结果可以明显发现本文算法会优于NMF和SNMF算法,与NMF算法相比,在系数矩阵H的稀疏度上,本文算法占有很大优势,信噪比方面有一定提升,与SNMF算法相比,在系数矩阵H的稀疏度上,与本文算法比较接近,但是在信噪比方面,本文则会优秀很多。
在实验过程中,发现对于降维后的维数r并没有较多分析,仅仅只有r<<min(m,n)的限定,而在试验中进行了对比与分析,参数选择如前所述,为了达到实验效果,本次试验采用固定迭代次数,取1000以确保实验效果。实验结果如图9所示。
从结果中可以看出当r<50时,信噪比上升较快,此时的由于低维度r的取值过低,导致原图像信息在分解的结果中信息丢失率比较严重,此时的r不适合于试验取值;当100<r<200时,信噪比随低维度r的变化较为平稳,可以用此段的r的值作为试验取值;当r>200时,信噪比已基本不随低维度r的变化而有较大变化了,且对于降维的效果也影响较大。试验中采用r=140是比较合适的,有一个不错的信噪比,即分离后重构的图像的效果比较好,同时,140相对于420和314都小很多,说明了降维效果也比较好。总体能够达到我们的要求。
表.5 试验收敛时信噪比与迭代次数
实验图片 |
实验结果信噪比 |
实验迭代次数 |
a |
30.0018 |
454 |
b |
30.0021 |
507 |
c |
30.0052 |
442 |
d |
30.0006 |
493 |
从图10可以看出迭代次数在400以下时,信噪比上升的速率很快,也就意味着在这期间由于迭代次数不够导致分解的结果不太理想,重构的图片与原图片差距较大,不能作为实验结果,当达到400至500次时信噪比上升的速率明显放缓,迭代次数的上升对于信噪比的影响不是很大,当迭代次数超过500次时,基本信噪比就不再变化了,已经接近实验比较好的结果了,此时再上升迭代的次数已经没有什么收益了,也不适宜实验的效率。按照前面实验中得出450~500次的迭代次数是比较合理的,基本已经达到要求。
经过实验和分析可以看出,将本文提出的对NMF算法中施加相关性约束和稀疏性约束的一种约束NMF算法运用于实际的图像中时,得出的结果较为满意,具有较强的稀疏性和信号较高的信噪比。