网刊加载中。。。

使用Chrome浏览器效果最佳,继续浏览,你可能不会看到最佳的展示效果,

确定继续浏览么?

复制成功,请在其他浏览器进行阅读

高分五号高光谱数据融合方法比较

  • 张立福 1
  • 赵晓阳 1,2
  • 孙雪剑 1
  • 黄海 3,4,2
  • 彭明媛 1,2
  • 岑奕 1
  • 涂宽 5
1. 中国科学院空天信息创新研究院 遥感卫星应用国家工程实验室, 北京 100101; 2. 中国科学院大学, 北京 100049; 3. 中国科学院月球与深空探测重点实验室, 北京 100012; 4. 中国科学院国家天文台, 北京 100012; 5. 二十一世纪空间技术应用股份有限公司, 北京 100096

最近更新:2022-10-31

DOI:10.11834/jrs.20229318

  • 全文
  • 图表
  • 参考文献
  • 作者
  • 出版信息
EN
目录contents

摘要

数据融合是解决高光谱卫星在时空分辨率等指标上受限的有效途径,探讨不同方法在GF-5高光谱数据上的融合效果,对GF-5高光谱数据的信息挖掘与推广应用有着重要意义。本文本着算法简单易用、适于推广的原则,采用GS(Gram-Schmidt)葛兰—施密特正交变换融合算法、GSA(GS Adaptive)自适应GS融合算法、CNMF(Coupled Non-negative Matrix Factorization)耦合非负矩阵分解融合算法、CRISP-W(Color Resolution Improvement Software Package with Wavelet transform)基于小波变换和CRISP-B(Color Resolution Improvement Software Package with Butterworth)基于巴特沃斯滤波器的分辨率提升融合算法、GLP(Generalized Laplacian Pyramid)广义拉普拉斯金字塔融合算法共6种融合方法,分别对BJ-2、GF-2、GF-1、GF-1C、GF-1D国产卫星多光谱数据与GF-5高光谱数据进行融合实验。通过目视分析、指标评价(相关系数、通用图像质量指标、峰值信噪比、光谱角、全局综合误差)、分类应用、时间成本4种方式对融合结果进行综合比较分析。结果表明,相融合的一组图像系列相同、空间分辨率相差越小,融合结果越好。CRISP-B、CRISP-W、GLP在提升空间分辨率、光谱保真度方面能达到较好的平衡,空间重建方面,GLP稍优且更稳定,CRISP-B、CRISP-W则在光谱信息保持方面稳定性更强且效果更好。数据源会对融合方法产生一定的影响,在光谱特征信息提取、分析等对光谱保真度要求高的工作中,GLP更适合同源数据(如GF-5与GF-1/1C/1D/2)融合,而在多源数据间(如GF-5与BJ-2)进行融合时,则优先选择CRISP-W。CNMF存在一定程度的色彩畸变,且运行时间较长。GSA、GS融合效果最差,其中,GSA不论是光谱保持能力还是空间分辨率提升能力均较GS更稳定。在小样本高光谱图像分类应用中,CRISP-B融合结果分类效果稳定,分类精度较高。GSA融合结果空间细节丰富,虽光谱失真较为严重,但同时增大了地物光谱分离度,仍适用于准确勾勒建筑物、道路等地物。本研究为GF-5高光谱数据与其他国产卫星多光谱数据融合方法的选择提供参考,有助于高分五号高光谱数据的应用与推广。

1 引 言

高分五号卫星(GF-5)于2018-05-09成功发射。这是世界上第一颗可以综合观测陆地和大气的全谱段高光谱卫星,在中国高分辨率对地观测系统重大专项中承担重要角色(

孙允珠 等,2018)。GF-5高光谱数据的空间分辨率为30 m(刘银年,2018),相比高分一号(GF-1)、高分二号(GF-2)等多光谱卫星,其空间分辨率在许多应用中还存在不足(Li等,2018)。高光谱数据与全色/多光谱数据通过空谱融合,可以有效提升高光谱数据的空间分辨率、全色/多光谱数据的光谱分辨能力(Tong等,2014)。对此,国内外学者展开了广泛研究,面向空间维提升(张立福 等,2019),融合算法主要分为成分替换法、多分辨率分析法。成分替换法首先进行投影变换,然后将原始图像的空间信息成分替换为高空间分辨率图像,最后通过逆变换获得空间维提升的遥感数据(Ghassemian,2016)。如Laben和Brower(2000)最早将葛兰—施密特正交变换GS(Gram-Schmidt)算法应用到了遥感图像融合处理上,该算法已作为模块封装在ENVI中,得到了广泛应用。为解决成分替换法对光谱信息保持能力欠佳的问题,学者发展了自适应GSA(GS Adaptive)、GIHSA(GIHS Adaptive)等(Aiazzi等,2007),在一定程度上减少了融合过程中的光谱失真。在多分辨率分析法中,原始图像首先被分解成了不同分辨率的一系列图像,然后将它们在不同分辨率水平上进行融合,最后通过逆变换获得融合数据。如利用小波变换(Gomez等,2001)、广义拉普拉斯金字塔法GLP(Generalized Laplacian Pyramid)(Schmitt和Zhu,2016)。相对于成分替换法来说,多分辨率分析法能较好地保持光谱信息,但它同时要求严格配准数据,不然会出现空间上的失真(Yuan等,2018)。面向光谱维提升的融合方法主要分为概率统计法、光谱解混法、稀疏解混法、人工智能法等。概率统计法如色彩分辨率提升CRISP(Color Resolution Improvement Software Package)(Winter等,2007)高光谱融合算法,基于概率统计理论,结合高空间分辨率和高光谱数据间的物理相关性,实现高光谱高空间数据的最优拟合。光谱解混是近年来针对高光谱影像分析的有力手段,如耦合非负矩阵分解CNMF(Coupled Non-negative Matrix Factorization)算法,它首先对多光谱和高光谱数据进行混合像元分解,然后将多光谱数据的端元丰度信息与高光谱数据的端元光谱信息进行混合,从而获得高空间高光谱数据(Yokoya等,2012)。基于稀疏表达法进行融合时,首先将原始数据分解为稀疏系数、字典矩阵,然后加入稀疏约束条件,进而求解出稀疏系数,从而获得融合重建数据(张立福 等,2019),但该方法相对复杂,运行时间较长。基于人工神经网络、卷积神经网络等,人工智能可以学习系统输入输出的非线性关系以进行较好的融合,但其需要调整的参数较多,且对具体的实验数据较为依赖。在提升高光谱数据的空间分辨率、保持其光谱信息等方面,上述算法各具优势。Yokoya等(2017)研究了最新的10种融合方法对高光谱和的多光谱数据的融合表现,该研究为高光谱数据与多光谱数据的融合提供了一个很好的方向。然而,该研究只使用了模拟的高光谱数据,因此结论在实际应用中受到限制,如GF-5成像仪硬件设计不同、空间分辨率不同、有实际采集几何影响等。Ren等(2020)对GF-5高光谱数据与GF-1、GF-2、Sential-2A多光谱数据进行了融合,全面评价了9种融合方法的表现。但其涉及的国产多光谱卫星数据源较少,如未涉及北京二号(BJ-2)、GF-1C、GF-1D等。

为进一步探讨高光谱融合方法在GF-5高光谱数据上的适用性,本着算法简单易用,适于推广的原则,本文在成分替换类、多分辨率分析类、概率统计类、光谱分解类融合算法中,具体选择GS、GSA、GLP、CRISP-W、CRISP-B、CNMF共6种常用融合方法,分别对6组国产卫星多光谱数据与GF-5高光谱数据进行融合,通过目视分析、指标评价、分类应用、时间成本4种方式,综合比较6种融合方法的空间分辨率增强能力、光谱特征保持能力以及稳定性,以期为GF-5高光谱数据在融合方面的应用研究提供参考。

2 方 法

2.1 融合方法

2.1.1 GS、GSA融合方法

GS融合算法基本流程如图1所示,首先对多光谱图像进行光谱重采样,得到模拟原始全色图像。利用滤波函数对低空间高光谱图像进行光谱重采样,得到一个模拟低空间分辨率全色图像。然后将其作为第一波段与原始拥有N波段的高光谱图像组合,进行GS变换以得到各个正交分量。再将模拟全色图像与GS变换后的第一分量进行直方图匹配,使其均值和方差近似一致,然后进行分量替换以及GS逆变换,接着去除N+1个波段融合图像中的多余波段,从而获得最终的融合图像。

GSA与GS的区别在于生成模拟低分辨率全色图像的方式不同,GS利用滤波函数生成,而GSA则基于高光谱和多光谱波段之间的相关性,将波段分组融合。数学模型可写为

HSFik=HSUik+gikMSi-Iik (1)
gik=covIik,HSUikvarIik (2)
Iik=k=1nwikHSUik (3)

式中,下标iki=1,,m; k=1,, n代表第i组数据中的第k个波段,mn分别表示图像的行列数,HSF为融合后高光谱图像,HSU为上采样高光谱图像,MSi代表多光谱的第i个波段。gikwik分别为正向变换系数与变换系数,I代表强度分量。

图1  GS融合算法基本流程

Fig. 1  Basic process of GS fusion algorithm

2.1.2 GLP融合方法

GLP融合方法根据高光谱图像与多光谱图像之间的相关系数将高光谱图像的波段分为不同的组,并对每组HS图像进行增强。它使用高斯调制传递函数滤波器对多光谱图像进行低通滤波。从原始多光谱数据中减去滤波后的图像,得到高空间细节图像,然后利用全局增益系数将提取的细节图像注入高光谱图像。

HSFik=HSUik+aikMSi-MSiL (4)

式中,MSiL代表第i个过滤的多光谱波段。aiki=1,,m;k= 1,,n代表第i组数据中的第k个高光谱波段的增益系数。

2.1.3 CRISP-B、CRISP-W融合方法

CRISP融合方法是一个包括高光谱数据线性重构和图像滤波等技术在内的组合算法模型。高光谱图像可以用矩阵表示为

H=hK1T,hK2T,,hKRT (5)

式中,hKiT表示图像中第i个像元的一条有K个波段的光谱列向量,R表示图像的总像素数。同理,拥有L个波段的多光谱图像可以表示为

M=mL1T,mL2T,,mLRT (6)

假设多光谱图像为高光谱图像在光谱上的线性采样,可得:

M=FH+ε (7)

式中,F是拟合矩阵,一般表示为多光谱数据的光谱响应矩阵。ε则代表高斯随机误差。而CRISP模型将(7)逆向变换,表达为

GM=H+r (8)

式中,G是转换矩阵,大小为K×Lr代表高斯随机误差。基于最小二乘法,G可以近似为

G=HMTMMT-1 (9)

一条重构的高光谱向量可以通过下式估算:

h1i'h2i'hKi'=g1,1g1,2g1,Lg2,1g2,2g2,LgK,1gK,2gK,L·m1im2imLi (10)
hKi'T=GmLiT (11)

式中,重构高光谱图像的第i个像元向量用hKi'T表示,原始多光谱图像的第i个像元向量则用mLiT表示。得到重构的高光谱图像H',其空间纹理特征与原始多光谱图像完全一致。该过程可以表达为

H'=GM (12)

最小二乘法求解上述过程会损失部分光谱特征,因此当CRISP得到重构高光谱后,还引入了频域滤波(巴特沃斯和小波变换两种滤波器,分别记作CRISP-B和CRISP-W)。

巴特沃斯基本原理是先通过离散余弦变换将图像由空间域转换到频率域,然后通过巴特沃斯低通和高通滤波函数保留原始高光谱图像H和重构高光谱图像H'的低频和高频信息,将二者相加后再反变换到空间域,便可得到最终的融合图像。

小波变换是基于傅里叶变化的一种新变换分析方法。它可以将遥感图像从不同层级进行分解,从而得到图像的低频部分和高频部分。每一更深层级的小波变换,都可以将上一层级的低频部分,继续分解为图像的低频分量以及垂直、水平、对角线方向上的高频分量。CRISP-W算法是将重构高光谱图像H'和原始高光谱图像H进行小波变换,然后将H'的所有高频分量和H的低频分量重新组合,通过小波逆变换得到最终结果。

2.1.4 CNMF融合方法

根据CNMF融合算法,高空间分辨率的高光谱图像Z(ZRλh×Lm)可以空间退化为低空间分辨率的高光谱图像X(XRλh×Lh),同时也可以光谱退化为高空间分辨率多光谱图像Y(YRλm×Lm),且满足λh>λm,Lh<Lm,因此这三者可以表示为

X=ZS+Es (13)
Y=RZ+Er (14)

式中,SRLm×Lh为空间退化矩阵;RRλm×λh为光谱退化矩阵;EsEr为残差项,RS为非负的且假设它们已知。根据线性混合模型,Z可以表示为

Z=WH+N (15)

式中,端元矩阵用WRλh×D表示,端元数据为D,残差为N,丰度矩阵为HRD×Lm,且W0,H0,将式(15)代入式(13)式(14)中,可得:

XWHh (16)
YWmH (17)

经过图像空间的退化关系可得:

HhHS (18)
WmRW (19)

式中,HhRD×Lh为高光谱图像丰度矩阵,WmRλm×D为多光谱端元矩阵。

WWXHhTWHhHhT (20)
WWXHhTWHhHhT (20)
HhHhWTXWTWHh (21)
WmWmYHTWmHHT (22)
HHWmTYWmTWmH (23)

将乘法更新法则应用于XY的NMF解混:

(1)对X进行第一次分解:初始化W,用式(21)更新HhW固定;WHh均通过式(20)式(21)进行更新,直到下一次收敛。

(2)对Y进行NMF分解:通过式(19)初始化Wm,用式(23)更新HWm固定;然后通过式(22)和(23)优化WmH

(3)对X进行后续NMF分解:通过式(18)初始化Hh,用式(20)更新WHh固定;然后通过式(20)和(21)优化WHh

(4)重复步骤(2)(3),直到收敛,从而完成融合:

Z=WH (24)

2.2 评价方法

评价融合结果一般可通过主观和客观的方式进行。主观评价主要是通过目视判读,从图像纹理、色彩变化上定性评估融合的效果。其简单、直接,但具有主观性、不全面性。定量化评估,如使用不同评价指标进行评估、通过地物分类等信息提取类应用进行评估,则为客观评价。本文首先对融合结果进行了目视分析,然后计算了融合图像与参考图像(原始的高光谱图像)之间的相关系数、通用图像质量指标、峰值信噪比、光谱角、全局综合误差指标(

Sun等,2015张立福等,2019),最后计算了融合数据分类结果的总体精度、Kappa系数,并比较了运行时间成本,通过这4种方式来综合分析不同融合算法的优缺点。

2.2.1 指标评价

(1)相关系数。相关系数CC(Correlation Coefficient)表示如下:

CCk=i=1mj=1nRk(i,j)-M(R)¯kFk(i,j)-M(F)¯ki=1mj=1nRk(i,j)-M(R)¯k2i=1mj=1nFk(i,j)-M(F)¯k2 (25)

式中,CCk为第k波段的相关系数。Rk(i,j)Fk(i,j)分别表示第k波段真实图像、融合图像第(i,j)位置像元值,MRk¯MFk¯分别表示二者在k波段的灰度均值。融合图像与参考图像越相似,空间信息变化越小,两者相关系数越接近1(

Alparone等,2007)。

(2)通用图像质量指标。通用图像质量指标UIQI(Universal Image Quality Index)(

Wang和Bovik,2002)表示如下:

UIQIk=SD(RF)kSD(R)kSD(F)k2M(R)¯kM(F)¯kM(R)¯k2+M(F)¯k22SD(R)kSD(F)kSD(R)k2+SD(F)k2 (26)

式中,第一项计算了两幅图像的相关系数,第二项衡量了两者的平均亮度相似度,第三项则是测量了两幅图像对比度的相似度。SDRFk表示两幅图像对应波段的协方差。UIQI越趋近于1,融合效果越好(

Li等,2012)。

(3)峰值信噪比。峰值信噪比PSNR(Peak Signal to Noise Ratio)评估图像的空间信息重构质量,即较原始图像,融合图像的有效信息是否得到增强(

Huynh-Thu和Ghanbari,2008)。一般来说,融合图像空间信息质量越好,PSNR越大。表达式如下:

PSNRk=10log10m, n, (max(k))2/i=1mj=1nRk(i,j)-Fk(i,j)2 (27)

式中,maxk表示第k波段出现的最大值,mn分别表示图像的行列数。

(4)光谱角。光谱角SAM(Spectral Angle Mapper)用来描述两个光谱间的相似度,两幅图像第i,j个像元的光谱角定义为(

Yuhas等,1992):

θ(i,j)=arccosR(i,j),F(i,j)R(i,j)F(i,j) (28)

式中,Ri,jFi,j分别表示两幅图像第i,j个像元的两条光谱向量。一般SAM越小,匹配程度越高,完全匹配SAM为0(

Zhang等,2009)。

(5)全局综合误差指标。全局综合误差指标ERGAS(Erreur Relative Globale Adimensionnelle de Synthèse)(

Wald,2000)主要考量图像辐射值上的畸变(Zurita-Milla等,2008),表达式为

ERGAS=100d1Nk=1NMSEk(M(R)¯k)2 (29)

式中,d为空间降采样因子,N表示融合数据的波段数。ERGAS理想值为0,值越小越好。

2.2.2 计算效率评价

记录运行时间以定量评估不同融合方法的计算效率。GS融合实验在ENVI 5.3中交互实现,GSA、CNMF、GLP、CRISP-B、CRISP-W融合实验均基于Matlab 2019实现。运行平台为Win10系统、Intel Core i5处理器、8 GB内存的计算机。

2.2.3 分类精度评价

通过检查所有融合结果对像素级分类任务的影响,对融合数据进行应用驱动评估。结合Google Earth真实图像、地面调查资料,勾画精细的地面真实地物类型Ground truth,从中产生3%随机点用作训练样本,剩下的97%用作测试样本。基于训练样本,运用随机森林方法对融合结果进行分类,比较不同融合结果中测试样本的总体分类精度、Kappa系数,探讨样本数量小的情况下,适合GF-5高光谱融合数据分类的最佳融合方法。

(1)总体分类精度。总体分类精度OA(Overall Accuracy),表示融合图像的分类结果与地面对应区域的实际类型相一致的概率。

(2)Kappa系数。Kappa系数KC(Kappa coefficient)是一种评价分类结果与实际地物之间吻合度或精度的指标,可表示为

KC=Ni=1cnii-i=1c(ni+n+i)N2-i=1c(ni+n+i) (30)

式中,c是总类别数,即混淆矩阵中的总列数;nii是正确分类的数目,即混淆矩阵第i行、第i列上的像元数量;ni+n+i分别是第i行、第i列的总像元数;N是用于精度评价的总像元数。

3 实验与结果分析

3.1 实验

本文选取鄱阳湖、黄河入海口、洞庭湖区域的GF-1/1C/1D、GF-2、BJ-2的多光谱数据,以及相近时间、相同区域的GF-5高光谱数据,分别进行融合实验。数据如表1所示。

实验数据包含了裸土、建筑、道路、水体、农田、养殖区、水田等典型地物,除E组外,地物均较为破碎。图像数据均为一级产品,像元为辐亮度值。利用ENVI分别对每一组(A-F)数据进行预处理。首先利用RPC Orthorectification模块对多光谱、高光谱数据进行正射校正,并利用三次卷积法,将A-B组中的GF-5高光谱数据重采样为4 m、C-F组的GF-5高光谱数据重采样为8 m。然后以多光谱数据为基准图像、高光谱数据为待配准图像,选取同名控制点,对每组图像进行配准。最后,选取400像素×400像素的相同区域进行融合实验。

表1  实验数据
Table 1  Experiment data
编号区域中心经度中心纬度多光谱时间高光谱时间
A 鄱阳湖 116.39°E 28.77°N BJ-2 (4 m) 2018-10-02 GF-5(30 m) 2018-10-05
B 洞庭湖 112.63°E 28.67°N GF-2 (4 m) 2018-10-05 2018-10-05
C 鄱阳湖 116.31°E 28.75°N GF-1 (8 m) 2018-10-06 2018-10-07
D 鄱阳湖 116.44°E 28.72°N GF-1C (8 m) 2018-10-08 2018-10-07
E 黄河入海口 118.94°E 37.60°N GF-1D (8 m) 2018-11-03 2018-11-01
F 洞庭湖 112.63°E 28.67°N GF-1D (8 m) 2018-10-06 2018-10-05

其中,GS滤波函数采用光谱响应函数;CNMF预设解混的端元数量由虚拟维数方法自动确定;CRISP截止频率为20,其中CRISP-W小波滤波器的层数为3;GLP利用最小二乘法计算回归系数矩阵。

3.2 目视分析

为便于目视判读及定性评价,融合结果用波段组合相同的真彩色图像显示(GF-1/1C/1D,R:Band3(680 nm),G:Band2(576 nm),B:Band1(502 nm);GF-2,R:Band3(660 nm),G:Band2(555 nm),B:Band1(485 nm);BJ-2,R:Band3(635 nm),G:Band2(550 nm),B:Band1(475 nm);GF-5及融合图像,R:Band59(638.4 nm),G:Band38(548.47 nm)),B:Band17(458.59 nm),拉伸方式也保持一致(图2)。

总体来说,GS不能有效提高空间分辨率,融合结果最为模糊,在地物破碎、类型复杂的A、B、C、D、F中尤为明显。CNMF存在较为明显的色彩畸变,在6组实验中均有显现。GSA、GLP均能有效提高空间分辨率,地物的纹理较为清晰,但在色彩保真度方面表现略有不同,在B、C、E、F中,GSA融合结果色彩与多光谱图像相近,而GLP融合结果色彩则与高光谱数据更接近,这在一定程度上说明,GLP的光谱保真度较GSA更好。CRISP-B和CRISP-W方法的融合结果较为接近,色调自然,在6组实验中都表现出不错的光谱保持能力,但空间分辨率提升能力不及GLP。

综上所述,地物边界规则类型较为单一时,六种方法差别不大。地物破碎类型复杂时,GLP、CRISP-B、CRISP-W的光谱保持能力较好,GLP、GSA空间分辨率提升能力较好。GS不能有效提高空间分辨率,CNMF融合结果存在色彩畸变。

3.3 指标评价

表2表6展示了基于不同实验数据、不同融合方法所生成结果的相关指标值,CC、UIQI、PSNR均为波段平均值,SAM为像元平均值。纵向对比,可以看出C、D、E、F组较A、B组,CC、UIQI更接近1,PSNR更大、SAM更小。即融合后的图像与参考图像更相关,亮度、对比度更相近,空间重建质量、光谱保持能力更好。这是因为A、B两组实验均为4 m空间分辨率的多光谱图像和30 m空间分辨率的高光谱图像相融合,而C、D、E、F 4组实验均为8 m多光谱图像与30 m高光谱图像融合。ERGAS表现较为异常,这是因为E组数据有大范围水域,图像波段像素平均值较低,使得ERGAS值偏大,所以在实验组别间进行比较时,不应过分依赖ERGAS指标。不难发现,在多光谱数据空间分辨率为4 m的融合实验中,GF-2较BJ-2数据更优;而在8 m的实验中,GF-1D数据表现更好。综上,相融合的一组图像空间分辨率相差越小,融合效果越好。在将GF-5高光谱数据与4 m多光谱数据融合时,可优先选择GF-2,与8 m多光谱数据融合时,可优先选择GF-1D。

接着进行横向比较,即不同融合方法之间的比较。在空间重建方面,GLP最佳,CRISP-B、CRISP-W表现也较好,尤其是在融合GF-5与BJ-2时,3种融合结果CC均达到了0.84、UIQI均达到了0.81,PSNR均在24.5以上。这说明,图像来源对CRISP-B、CRISP-W、GLP的空间重建能力影响小。在光谱保持方面,CNMF、GLP、CRISP-B、CRISP-W均表现出不错的效果,CRISP-W较CRISP-B稍优,两者的SAM平均值较高,且标准差很低,表现稳定。GLP和CNMF的SAM平均值最高,但标准差较大,这主要是受A组多源融合实验的影响。因此,为较好的保持光谱信息,在同源数据间进行融合时,可优先选择GLP,在多源数据间进行融合时,则优先选择CRISP-W。GS、GSA的空间重建能力、光谱保持能力均较差,在实验A中的表现尤为明显。虽然GS方法在4组(BCDF)实验中表现较优,如4组CC、UIQI均值均约为0.95,PSNR的均值达到了最高,但这是由于GS方法未能有效去除高光谱数据的噪声值造成的,该部分内容将在讨论部分详细说明。此外,GSA不论是光谱保持能力还是空间分辨率提升能力均较GS更稳定。ERGAS是一个综合指标,在利用存在一定配准误差的图像作为参考图像时,具有一定的优势。综上所述,结合综合指标ERGAS,可以得出结论,CRISP-W方法效果最佳,且最稳定,CRISP-B、GLP效果稍次之,CNMF居中,GS、GSA效果最差。

表2  融合结果的评价指标(CC)
Table 2  Evalution indexes of fusion results (CC)
实验区域编号融合方法Mean
GSGSACNMFC-BC-WGLP
A 0.59 0.74 0.75 0.84 0.84 0.87 0.77
B 0.92 0.83 0.86 0.88 0.89 0.94 0.89
F 0.94 0.89 0.87 0.93 0.93 0.92 0.91
E 0.92 0.92 0.92 0.96 0.96 0.95 0.94
D 0.97 0.94 0.91 0.96 0.96 0.96 0.95
C 0.97 0.95 0.94 0.95 0.95 0.96 0.95
Mean 0.89 0.88 0.87 0.92 0.92 0.93
Std 0.13 0.07 0.06 0.04 0.04 0.03

注:   C-B:CRISP-B;C-W:CRISP-W;加粗为最优值。

表3  融合结果的评价指标(UIQI)
Table 3  Evalution indexes of fusion results (UIQI)
实验区域编号融合方法Mean
GSGSACNMFC-BC-WGLP
A 0.58 0.70 0.74 0.81 0.82 0.85 0.75
B 0.91 0.81 0.85 0.87 0.88 0.93 0.88
F 0.94 0.88 0.86 0.92 0.93 0.91 0.91
E 0.92 0.91 0.92 0.96 0.96 0.95 0.94
D 0.97 0.94 0.91 0.95 0.96 0.96 0.95
C 0.97 0.95 0.93 0.94 0.95 0.96 0.95
Mean 0.88 0.87 0.87 0.91 0.91 0.93
Std 0.14 0.09 0.06 0.05 0.05 0.04

注:   C-B:CRISP-B;C-W:CRISP-W;加粗为最优值。

表4  融合结果的评价指标(PSNR)
Table 4  Evalution indexes of fusion results (PSNR)
实验区域编号融合方法Mean
GSGSACNMFC-BC-WGLP
A 21.57 21.18 24.54 25.01 23.39 25.17 23.48
B 32.21 23.92 25.84 26.55 26.28 29.00 27.30
E 27.61 26.79 29.74 30.06 27.58 29.03 28.47
F 38.78 27.36 29.29 29.74 28.23 28.76 30.36
C 34.19 29.72 29.44 30.33 29.66 31.37 30.79
D 32.91 30.23 31.16 31.85 29.54 31.55 31.21
Mean 31.21 26.53 28.34 28.93 27.45 29.15
Std 5.42 3.17 2.34 2.36 2.15 2.11

注:   C-B:CRISP-B;C-W:CRISP-W;加粗为最优值。

图2  不同融合方法、数据的融合结果(左侧A-F为实验编号,括号中注明的传感器为多光谱数据来源)

Fig. 2  Fusion results of different methods and data(On the left, A-F number the experienments and the sensors noted in brackets are the multispectral data sources)

表7可以看出,GSA、C-W、GLP、C-B、CNMF运行时间依次增大,CNMF运行时间过长不可取,其余5种方法运行成本均在可接受范围内。

3.4 分类应用

图像景观类型对评估结果有很大的影响。景观大小和自然特征的多样性,导致空间信息和光谱信息高度可变。表8展示了A、B、C、D、F研究区的标签类别及用于训练、测试的具体样本数量。每类样本数量相等,避免了地物占比不同影响评判结果。E研究区地物种类较少,为保证方法对比的有效性,未对E进行分类。此外,精细的、以地物边界为主的Ground truth,也可有效鉴别融合结果分类好坏。

表5  融合结果的评价指标(SAM)
Table 5  Evalution indexes of fusion results (SAM)
实验区域编号融合方法Mean
GSGSACNMFC-BC-WGLP
A 4.49 4.31 2.85 2.77 2.33 2.17 3.15
B 3.20 3.74 1.79 1.78 2.36 2.11 2.50
F 2.48 2.36 1.49 2.06 2.13 1.95 2.08
E 2.56 2.20 1.98 1.85 1.83 1.75 2.03
C 2.20 1.67 1.94 1.72 2.13 1.88 1.92
D 1.94 1.06 1.81 1.68 1.83 1.63 1.66
Mean 2.81 2.56 1.98 1.98 2.10 1.91
Std 0.84 1.13 0.42 0.38 0.21 0.19

注:   C-B:CRISP-B;C-W:CRISP-W;加粗为最优值。

表6  融合结果的评价指标(ERGAS)
Table 6  Evalution indexes of fusion results (ERGAS)
实验区域编号融合方法Mean
GSGSACNMFC-BC-WGLP
E 4.14 4.10 3.49 3.50 3.26 3.10 3.60
F 2.31 2.29 1.80 1.91 1.80 1.68 1.97
C 2.07 2.12 2.00 1.71 2.04 1.84 1.96
A 2.90 2.29 1.61 1.51 1.40 1.31 1.84
D 1.63 1.88 1.85 1.66 1.68 1.54 1.71
B 1.67 1.55 1.05 0.87 1.15 1.05 1.22
Mean 2.45 2.37 1.97 1.86 1.89 1.75
Std 0.87 0.81 0.75 0.80 0.68 0.65

注:   C-B:CRISP-B;C-W:CRISP-W;加粗为最优值。

表7  不同融合方法的运行时间比较
Table 7  Runtimes comparison of different fusion methods /s
实验区域编号GSGSAC-WGLPC-BCNMF
A ENVI交互 5.4 8.6 13.5 16.4 417.3
B 5.0 8.6 13.4 16.1 422.9
C 5.0 8.7 13.3 16.4 426.0
D 5.0 8.6 13.4 16.1 474.2
E 5.0 8.8 13.6 16.7 466.0
F 5.0 8.7 14.7 16.6 442.4
Mean 5.1 8.7 13.7 16.4 441.5

注:   C-B:CRISP-B;C-W:CRISP-W。

表8  融合结果的评价指标 ( 像素 )
Table 8  Evalution indexes of fusion results
实验区域编号裸土建筑道路水体农田养殖区水田
TrainTestTrainTestTrainTestTrainTestTrainTestTrainTestTrainTest
A 30 970 30 970 30 970 30 970 30 970
B 30 970 30 970 30 970 30 970 30 970
C 30 970 30 970 30 970 30 970 30 970 30 970 30 970
D 30 970 30 970 30 970 30 970 30 970 30 970 30 970
F 30 970 30 970 30 970 30 970 30 970 30 970

表9展示了不同研究区原始多光谱、高光谱图像、各融合结果的分类总体精度、Kappa系数。GS的分类结果最差,这是因为GS的空间分辨率提升能力最差,从而影响分类,尤其是精细分类。CNMF、GLP、CRISP-W效果居中,其中,C、D、F组较A、B组效果好,这可能是由于C、D、F组相融合的一组图像空间分辨率相差小,融合结果较好,因此分类结果也相对较好。在A、B两组中,GSA表现出了较强的优越性,虽然GSA融合结果的光谱保持性较差,但成分替换增大了地物的光谱发散,而在目视分析中也发现,GSA的空间信息与多光谱图像更接近,空间分辨率得到了有效的提高,进而提高了精细分类的精度。CRISP-B在所有融合结果的分类中都表现出了不错的效果,在样本数量小的情况下,较原始GF-5分类结果的总体精度提高了约4个百分点。图3展示了Ground truth,高光谱、多光谱图像的分类结果,以及效果较好的CRISP-W、GSA、CRISP-B融合图像分类结果。可以看出在A、B两组,GSA对道路的分割很清晰,而CRISP-B、CRISP-W分出的道路较粗,有的还呈片状,CRISP-B分出的绿地、裸土等基本呈块状。在C、D两组中,GSA的分类结果较为破碎。F组GSA对建筑、道路的分类精度高,使得总体的分类效果好。

表9  总体精度和Kappa系数
Table 9  Overall Occuracy (OA) and Kappa coefficient
参量实验区域编号高光谱图像多光谱图像GSCNMFGLPCRISP-WGSACRISP-B
OA/% A 81.32 87.96 87.84 87.46 85.57 86.04 88.72* 88.39*
B 81.09 82.08 78.14 79.98 81.24 81.05 88.29*** 83.59*
C 88.26 86.41 86.52 87.78 88.06 88.92* 84.36 90.06***
D 88.47 90.66 90.03 92.64* 92.53* 92.75* 90.24 93.43**
F 88.76 90.55 88.47 91.56 90.86 92.11* 94.06*** 93.25**
Kappa A 0.767 0.850 0.848 0.843 0.820 0.826 0.859* 0.855*
B 0.764 0.776 0.727 0.750 0.766 0.763 0.854*** 0.795*
C 0.863 0.841 0.843 0.857 0.861 0.871* 0.818 0.884**
D 0.866 0.891 0.884 0.914* 0.913* 0.916* 0.886 0.923**
F 0.865 0.887 0.862 0.899 0.890 0.905* 0.929*** 0.919**

注:   *越多,表示分类效果越好。

综上所述,对于区分度较大的水体来说,不同融合方法的分类精度无显著差别,而对于地物破碎、边界不易识别的建筑、道路等,GSA融合结果的分类效果最好,但GSA同时会使片状地物如裸土、农田等产生部分破碎。CRISP-B融合结果分类效果稳定,总体分类精度较高。

4 讨 论

(1)CRISP模型分别采用了小波和巴特沃斯两种滤波器来进行频域滤波过程。小波滤波器的优点是处理速度相对较快,但当图像配准精度较低时,融合图像会出现明显的“斑块”。而巴特沃斯滤波器则受到配准精度的影响相对较小,表现出较强的稳定性。本文实验数据配准精度较高,CRISP-W效果优于CRISP-B,可见融合前严格配准的必要性。然而,无论是小波滤波器还是巴特沃斯滤波器,不同的融合图像达到最佳融合效果所需要的滤波参数是不一样的,这就造成滤波过程对基函数、滤波截止频率、分解层数选择的困难。选择的参数不够理想,会造成融合图像空间和光谱上的进一步信息损失。产生这些问题的主要原因是CRISP模型中的图像重构过程完全是利用数理统计理论,从数学上进行统计优化,缺乏明确的物理意义,容易造成结果的不稳定性。此外,在CRISP模型中,对高空间分辨率高光谱图像的模拟过程与后续的滤波过程是完全独立的两个部分,如果能直接重构得到高空间分辨率的高光谱数据,而不再需要滤波融合过程,将是对CRISP算法的重要改进。

图3  分类结果(A-D、F为实验编号,括号中注明的传感器为多光谱数据来源)

Fig. 3  Classification result(A-D and F number the experienments. The sensors noted in brackets are the multispectral data sources)

(2)GSA融合算法在波段变换过程中,会对光谱信息造成损失,产生光谱畸变现象。好的融合算法不仅能使图像原始光谱信息得到较好的保持,还能使其空间纹理信息增加。但各种融合方法的融合质量,即光谱、空间保真性与遥感应用之间并无绝对相关性。本实验结果表明,GSA光谱失真严重,但同时增大了地物光谱分离度,尤其适用于道路、建筑的分类。下一步可针对植被监测、陆地生产力反演、环境评价、矿产勘测、目标探测等不同应用目的对融合算法做进一步分析。

(3)观察原始高光谱图像,发现除E组外,其他5组图像的前部分波段均有不同程度的噪声、条带,质量较差。从图4也可以看出(此处讨论光谱保持效果较好且最稳定的),E组各波段CC趋势最为平缓。其他5组实验中,CRISP-B、CRISP-W的前部分波段的CC值都有升高的趋势,这是因为基于CRISP-B、CRISP-W融合后的图像噪声、条带有很大的改善,与原始图像有较大的区别(图5),因此CC值较低,符合实际情况。但GS方法却保持了平缓,即与原始噪声大的图像非常相似,没有CRISP-B、CRISP-W去除噪声的作用。A、B、C、F实验中,GS方法的CC在约第80波段(中心波长728.06 nm)之后都有一个陡降的趋势,而其他方法反而有一个上升的趋势并保持在一个比较平稳的状态,这在一定程度上说明CRISP-B、CRISP-W在红外波段保持效果相对GS更有优势。

(a)  A(BJ-2)组实验

(a)  Test A(BJ-2)

(b)  B(GF-2)组实验

(b)  Test B(GF-2)

(c)  C(GF-1)组实验

(c)  Test C(GF-1)

(d)  D(GF-1C)组实验

(d)  Test D(GF-1C)

(e)  E(GF-1D)组实验

(e)  Test E(GF-1D)

(f)  F(GF-1D)组实验

(f)  Test F(GF-1D)

  

图4 各波段CC趋势(A-F分别为括号中注明的数据来源的对比结果)

Fig. 4 CC trend of different bands(A-F are the comparison results of the data sources indicated in brackets respectively)

(a)  原始影像

(a)  Original image

(b)  GS方法结果

(b)  The fusion result of GS

(c)  GRISP-B方法结果

(c)  The fusion result of GRISP-B

(d)  GRISP-W方法结果

(d)  The fusion result of GRISP-W

图5 A组原始、融合后第一波段结果图

Fig. 5 The first band image of origin and fusion result of test A

(4)CNMF模型目标函数具有非凸性,模型可能会收敛到局部极小值,使得分解矩阵不唯一,这可能是本文中CNMF融合效果不佳的原因。此外,CNMF需要预设解混的端元数量,不同的预设方法会导致不同的结果,可进一步探讨预设端元方法对融合结果的影响。

5 结 论

本文选取6种简单、易于推广的高光谱融合方法,对多种国产卫星多光谱数据与GF-5高光谱数据进行融合,通过目视分析、5种经典评价指标、分类应用精度、运行时间成本来综合比较分析了融合方法的优缺点。GLP的光谱信息保持、空间分辨率提升能力均较好,但在融合多源图像时,光谱保持能力不及CRISP-W。CRISP-B、CRISP-W能有效完成大部分融合任务,尤其适用于光谱提取、分析等对光谱保真度要求高的工作。其在融合多源图像时表现稳定,能有效去除原始GF-5图像的条带噪声,CRISP-B的分类效果较好且稳定。CNMF融合效果居中,但融合结果存在色彩畸变,且运行时间过长。GSA、GS融合效果最差,其中,GSA较GS更稳定。GSA融合结果空间细节丰富,虽有较为严重的光谱失真,但同时增大了地物光谱分离度,适用于道路、建筑等地物,即光谱、空间保真性与分类应用之间并无绝对相关性。此外,相融合的一组图像系列相同、空间分辨率相差越小,融合结果越好。本研究考虑了不同传感器-同场景、同传感器—不同场景、数据来源系列相同/不同的融合,为GF-5高光谱数据在融合方面的应用研究提供方法选择依据。

参考文献(References)

Aiazzi B, Baronti S and Selva M. 2007. Improving component substitution pansharpening through multivariate regression of MS + Pan data. IEEE Transactions on Geoscience and Remote Sensing, 45(10): 3230-3239 [DOI: 10.1109/TGRS.2007.901007] [百度学术] 

Alparone L, Wald L, Chanussot J, Thomas C, Gamba P and Bruce L M. 2007. Comparison of pansharpening algorithms: outcome of the 2006 GRS-S data-fusion contest. IEEE Transactions on Geoscience and Remote Sensing, 45(10): 3012-3021 [DOI: 10.1109/TGRS.2007.904923] [百度学术] 

Ghassemian H. 2016. A review of remote sensing image fusion methods. Information Fusion, 32: 75-89 [DOI: 10.1016/j.inffus.2016.03.003] [百度学术] 

Gomez R B, Jazaeri A and Kafatos M. 2001. Wavelet-based hyperspectral and multispectral image fusion//Proceedings of SPIE 4383, Geo-Spatial Image and Data Exploitation II. Orlando, FL, United States: SPIE: 36-42 [DOI: 10.1117/12.428249] [百度学术] 

Huynh-Thu Q and Ghanbari M. 2008. Scope of validity of PSNR in image/video quality assessment. Electronics Letters, 44(13): 800-801 [DOI: 10.1049/el:20080522] [百度学术] 

Laben C A and Brower B V. 2000. Process for enhancing the spatial resolution of multispectral imagery using pan-sharpening. U.S., 6011875 [百度学术] 

Li D, Hao M F, Zhang J Q, Hu B and Lu Q Y. 2012. A universal hypercomplex color image quality index//Proceedings of 2012 IEEE International Instrumentation and Measurement Technology Conference. Graz, Austria: IEEE: 985-990 [DOI: 10.1109/I2MTC.2012.6229639] [百度学术] 

Li X T, Lu J X, Song X N, Sun Y Y, Li L, Lei T J and Qu W. 2018. Application of the GF satellite data in flood disaster monitoring//Proceedings of 2018 IEEE International Geoscience and Remote Sensing Symposium. Valencia, Spain: IEEE: 7293-7296 [DOI: 10.1109/IGARSS.2018.8519548] [百度学术] 

Liu Y N. 2018. Visible-shortwave infrared hyperspectral imager of GF-5 satellite. Spacecraft Recovery and Remote Sensing, 39(3): 25-28 [百度学术] 

刘银年. 2018. “高分五号”卫星可见短波红外高光谱相机的研制. 航天返回与遥感, 39(3): 25-28 [DOI: 10.3969/j.issn.1009-8518.2018.03.003] [百度学术] 

Ren K, Sun W W, Meng X C, Yang G and Du Q. 2020. Fusing China GF-5 hyperspectral data with GF-1, GF-2 and sentinel-2A multispectral data: which methods should be used? Remote Sensing, 12(5): 882 [DOI: 10.3390/rs12050882] [百度学术] 

Schmitt M and Zhu X X. 2016. Data fusion and remote sensing: an ever-growing relationship. IEEE Geoscience and Remote Sensing Magazine, 4(4): 6-23 [DOI: 10.1109/MGRS.2016.2561021] [百度学术] 

Sun X J, Zhang L F, Yang H, Wu T X, Cen Y and Guo Y. 2015. Enhancement of spectral resolution for remotely sensed multispectral image. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 8(5): 2198-2211 [DOI: 10.1109/jstars.2014.2356512] [百度学术] 

Sun Y Z, Jiang G W, Li Y D, Yang Y, Dai H S, He J, Ye Q H, Cao Q, Dong C Z, Zhao S H and Wang W H. 2018. GF-5 Satellite: overview and application prospects. Spacecraft Recovery and Remote Sensing, 39(3): 1-13 [百度学术] 

孙允珠, 蒋光伟, 李云端, 杨勇, 代海山, 何军, 叶擎昊, 曹琼, 董长哲, 赵少华, 王维和. 2018. “高分五号”卫星概况及应用前景展望. 航天返回与遥感, 39(3): 1-13 [DOI: 10.3969/j.issn.1009-8518.2018.03.001] [百度学术] 

Tong Q X, Xue Y Q and Zhang L F. 2014. Progress in hyperspectral remote sensing science and technology in China over the past three decades. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 7(1): 70-91 [DOI: 10.1109/jstars.2013.2267204] [百度学术] 

Wald L. 2000. Quality of high resolution synthesised images: is there a simple criterion?//Proceedings of the 3rd Conference “Fusion of Earth Data: Merging Point Measurements, Raster Maps and Remotely Sensed Images”. Sophia Antipolis: Thierry Ranchin, Lucien Wald Editors: 99-103 [百度学术] 

Wang Z and Bovik A C. 2002. A universal image quality index. IEEE Signal Processing Letters, 9(3): 81-84 [DOI: 10.1109/97.995823] [百度学术] 

Winter M E, Winter E M, Beaven S G and Ratkowski A J. 2007. Hyperspectral image sharpening using multispectral data//Proceedings of 2007 IEEE Aerospace Conference. Big Sky, MT: IEEE: 1-9 [DOI: 10.1109/AERO.2007.353060] [百度学术] 

Yokoya N, Yairi T and Iwasaki A. 2012. Coupled nonnegative matrix factorization unmixing for hyperspectral and multispectral data fusion. IEEE Transactions on Geoscience and Remote Sensing, 50(2): 528-537 [DOI: 10.1109/tgrs.2011.2161320] [百度学术] 

Yokoya N, Grohnfeldt C and Chanussot J. 2017. Hyperspectral and multispectral data fusion: a comparative review of the recent literature. IEEE Geoscience and Remote Sensing Magazine, 5(2): 29-56 [DOI: 10.1109/MGRS.2016.2637824] [百度学术] 

Yuan Q Q, Wei Y C, Meng X C, Shen H F and Zhang L P. 2018. A multiscale and multidepth convolutional neural network for remote sensing imagery pan-sharpening. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 11(3): 978-989 [DOI: 10.1109/jstars.2018.2794888] [百度学术] 

Yuhas R H, Goetz A F H and Boardman J W. 1992. Discrimination among semi-arid landscape endmembers using the Spectral Angle Mapper (SAM) algorithm. Proceedings of the 1992 Summaries of the Third Annual JPL Airborne Geoscience Workshop. Pasadena: JPL, 1992: 147-149. [百度学术] 

Zhang L F, Peng M Y, Sun X J, Cen Y and Tong Q X. 2019. Progress and bibliometric analysis of remote sensing data fusion methods (19922018). Journal of Remote Sensing, 23(4): 603-619 [百度学术] 

张立福, 彭明媛, 孙雪剑, 岑奕, 童庆禧. 2019. 遥感数据融合研究进展与文献定量分析(19922018). 遥感学报, 23(4): 603-619 [DOI: 10.11834/jrs.20199073] [百度学术] 

Zhang Y F, De Backer S and Scheunders P. 2009. Noise-resistant wavelet-based bayesian fusion of multispectral and hyperspectral images. IEEE Transactions on Geoscience and Remote Sensing, 47(11): 3834-3843 [DOI: 10.1109/tgrs.2009.2017737] [百度学术] 

Zurita-Milla R, Clevers J G P W and Schaepman M E. 2008. Unmixing-based Landsat TM and MERIS FR data fusion. IEEE Geoscience and Remote Sensing Letters, 5(3): 453-457 [DOI: 10.1109/lgrs.2008.919685] [百度学术] 

欢迎关注学报微信

遥感学报交流群