袁博，胡冰

(1.南陽(yáng)理工學(xué)院 計(jì)算機(jī)與信息工程學(xué)院，河南 南陽(yáng) 473004；2.南陽(yáng)理工學(xué)院 學(xué)報(bào)編輯部，河南 南陽(yáng) 473004)

0 引言

研究發(fā)現(xiàn)，光譜往往在很大程度上表征了地物的本征特性[1]，光譜分辨率的提高有助于對(duì)地物的精確識(shí)別和分類[2]。然而，高光譜遙感數(shù)據(jù)的數(shù)據(jù)量大、冗余信息多，高光譜影像處理算法普遍具有約束條件多、運(yùn)算量大、計(jì)算復(fù)雜等特點(diǎn)，為了降低計(jì)算復(fù)雜度并改善處理精度，不僅需要從原始數(shù)據(jù)中提取地物光譜信息，還應(yīng)該在保留原始光譜結(jié)構(gòu)與有效光譜信息的前提下，盡可能減少參與運(yùn)算的數(shù)據(jù)量，降低光譜信息冗余[3]。

高光譜的數(shù)據(jù)降維方法主要可分為維數(shù)約簡(jiǎn)和波段選擇兩大類。維數(shù)約簡(jiǎn)就是通過(guò)線性或非線性數(shù)學(xué)變換，將高維數(shù)據(jù)轉(zhuǎn)化為可解的低維數(shù)據(jù)[4]。波段選擇則是從多達(dá)數(shù)百個(gè)波段的高光譜影像中，選擇具有較好分類識(shí)別能力的波段組合，代替原始全波段數(shù)據(jù)參與分析與應(yīng)用[5]。

維數(shù)約簡(jiǎn)通常涉及數(shù)學(xué)運(yùn)算，容易破壞數(shù)據(jù)原始光譜結(jié)構(gòu)，不利于地物光譜信息提取，波段選擇則不存在此問(wèn)題。波段選擇方法的研究成果很多，包括基于稀疏的波段選擇[6](sparse based band selection，SpaBS)、稀疏矩陣分解的波段選擇法[7](sparse nonnegative matrix factorization，SNMF)、自適應(yīng)波段選擇[8](adaptive band selection，ABS)、加權(quán)概率原型分析的高光譜影像波段選擇[9]等。

不同于上述方法，本文提出一種相對(duì)熵結(jié)合互信息的光譜綜合信息量分析波段選擇(Kullback-Leibler divergence and mutual information based band selection，KLMIBS)方法來(lái)研究高光譜影像的波段選擇問(wèn)題。通過(guò)構(gòu)建高光譜影像光譜信息的總體分布模型進(jìn)行波段選擇，并利用一個(gè)高光譜影像公開(kāi)數(shù)據(jù)集進(jìn)行實(shí)驗(yàn)與精度分析，驗(yàn)證本文提出的KLMIBS方法的有效性。

1 光譜信息分布模型

高光譜影像波段選擇的結(jié)果中，除了要求不改變?cè)脊庾V數(shù)據(jù)結(jié)構(gòu)，保持原始波段物理含義不變之外，還應(yīng)該盡可能多地保留地物的有效光譜信息[10]。因此，需要對(duì)整幅影像的光譜信息總體分布情況進(jìn)行分析與建模。KLMIBS方法選擇相對(duì)熵(relative entropy，RE)[11]結(jié)合互信息[12](mutual information，MI)構(gòu)建光譜信息總體分布的定量描述模型，指導(dǎo)高光譜影像的波段選擇。

1.1 待選擇的波段數(shù)量

研究波段選擇問(wèn)題，首先需要確定待選擇的波段數(shù)量。高光譜影像處理算法通常要求選出的波段數(shù)量必須大于或等于端元數(shù)量[13]。端元數(shù)量可通過(guò)主成分分析[14](principal component analysis，PCA)、最小噪聲分離[15](minimum noise fraction，MNF)、奇異值分解[16](singular value decomposition，SVD)、噪聲子空間投影[17](noise subspace projection，NSP)和基于最小誤差高光譜信號(hào)識(shí)別[18](hyperspectral signal identification by minimum error，HySime)等方法得到。

其中，HFC和NWHFC方法都是基于相關(guān)矩陣和協(xié)方差矩陣特征值之間的協(xié)方差趨近于0的假設(shè)條件，利用虛擬維度估計(jì)端元數(shù)量的。不同之處在于，NWHFC在執(zhí)行HFC之前，對(duì)樣本相關(guān)矩陣和樣本協(xié)方差矩陣進(jìn)行了噪聲白化，去除了噪聲方差的相關(guān)性，使得虛擬維度的估計(jì)更加精確。HFC和NWHFC的主要問(wèn)題在于，估計(jì)的準(zhǔn)確度比較一般，且樣本數(shù)量不是很大時(shí)，其假設(shè)條件不成立。NSP僅采用白化后的樣本協(xié)方差矩陣進(jìn)行虛擬維數(shù)的估計(jì)，不受樣本數(shù)量的影響。在噪聲估計(jì)正確的前提下，尤其是樣本個(gè)數(shù)不多時(shí)，該方法優(yōu)于HFC和NWHFC。

對(duì)于高光譜遙感影像的端元數(shù)量估計(jì)問(wèn)題，PCA丟棄的主分量中可能包含真實(shí)端元；HFC和NWHFC的假設(shè)條件受樣本數(shù)量影響很大；NSP 則需要精確估計(jì)噪聲和給定誤警率。HySime盡管計(jì)算過(guò)程較復(fù)雜，但不需要任何參數(shù)且估計(jì)準(zhǔn)確度較高，綜合性能最優(yōu)。因此，選擇HySime算法估計(jì)高光譜影像的端元數(shù)量。

1.2 高光譜影像光譜差異信息模型

相對(duì)熵又稱KL散度(Kullback-Leibler divergence)或信息散度(information divergence)，是2個(gè)概率分布P和Q差別的非對(duì)稱性的度量，用來(lái)度量使用基于Q的編碼來(lái)編碼來(lái)自P的樣本，平均所需的額外比特個(gè)數(shù)。為避免混淆，下文將相對(duì)熵、KL散度和信息散度統(tǒng)稱為KL散度。

設(shè)P(x)和Q(x)是X取值的2個(gè)離散概率分布，則P對(duì)Q的KL散度為：

(1)

KL散度僅當(dāng)概率P和Q各自總和均為1，且對(duì)于任何x都滿足Q(x)>0及P(x)>0時(shí)，才有定義。

對(duì)于連續(xù)隨機(jī)變量，其概率分布P和Q可按積分方式定義為式(2)：

(2)

式中：p和q分別表示分布P和Q的密度。

從式(1)和式(2)可以看出，KL散度總是非負(fù)，當(dāng)且僅當(dāng)分布P和Q完全相同時(shí)，KL散度為0；P和Q的相似度越高，KL散度越??；相似度越低，KL散度越大?？梢詫⒎植糚和Q分別與高光譜圖像中2個(gè)波段的光譜數(shù)據(jù)進(jìn)行類比，則2波段數(shù)據(jù)的相似度越低，信息差異(即有效信息量)越大，則對(duì)應(yīng)的KL散度數(shù)值越大，反之則越小。因此，選擇利用KL散度參數(shù)描述波段間的有效光譜信息。

DKL(P‖Q)只能描述2個(gè)波段P和Q之間的光譜信息差異，而波段選擇過(guò)程需要對(duì)高光譜影像全部波段間光譜信息差異的整體分布進(jìn)行定量描述。因此，定義了一個(gè)KL散度矩陣MKL。對(duì)于高光譜數(shù)據(jù)X∈RN×M(R代表非負(fù)實(shí)數(shù)集，N是全部像元的數(shù)量，M是所有波段的數(shù)量)，MKL是一個(gè)M×M的方陣，如式(3)所示。

(3)

式中:Xi表示高光譜數(shù)據(jù)中所有像元在第i波段的數(shù)據(jù)；DKL(Xi‖Xj)表示所有像元從波段i到波段j的KL散度，可用公式(1)計(jì)算。

1.3 高光譜影像光譜冗余信息模型

為了更全面地描述光譜信息量整體分布情況，不僅需要建立有效光譜信息的整體分布模型，還需要對(duì)冗余光譜信息進(jìn)行整體分布的定量化描述，建立光譜冗余信息的整體分布定量模型。為此，論文引入了信息學(xué)中的互信息MI(mutual information)參數(shù)。

互信息是用來(lái)度量隨機(jī)變量之間獨(dú)立性的基本準(zhǔn)則，可以表示成KL散度的形式。多個(gè)隨機(jī)變量之間的互信息定義為其聯(lián)合概率密度函數(shù)與各邊緣密度函數(shù)乘積之間的KL散度，即：

(4)

將式(4)右邊繼續(xù)展開(kāi)：

(5)

由于KL散度的非負(fù)性可知：

(6)

互信息量I(ai；bj)在聯(lián)合概率空間P(XY)中的統(tǒng)計(jì)平均值，稱為平均互信息(average mutual information，AMI)，其計(jì)算公式如(7)所示。

(7)

利用平均互信息量是為了克服互信息量的隨機(jī)性，使其成為一個(gè)確定的量。關(guān)于平均互信息量的物理意義，以I(X；Y)為例，即在對(duì)Y一無(wú)所知的情況下，X的先驗(yàn)不定度與收到Y(jié)后關(guān)于X的后驗(yàn)不確定度之差?？梢杂闷骄バ畔⒘縼?lái)定量表示譜間冗余信息，為波段選擇提供依據(jù)。

基于上述分析，本文定義了一個(gè)平均互信息矩陣，記為MI，作為高光譜影像整體光譜冗余信息模型。MI是一個(gè)M行、M列的方陣，如式(8)所示。

(8)

式中:I(Xi；Xj)表示波段i和波段j之間所有像元反射率的平均互信息，該數(shù)值可依據(jù)式(7)進(jìn)行計(jì)算。

1.4 高光譜影像光譜綜合信息模型

為了更加全面地表述所有波段之間的光譜信息分布，同時(shí)抑制單參數(shù)的負(fù)面影響，算法將KL散度矩陣MKL和平均互信息量矩陣MI以不同的權(quán)重系數(shù)進(jìn)行結(jié)合，生成一個(gè)綜合光譜信息矩陣S，作為高光譜影像的整體光譜信息模型，并將其作為波段選擇過(guò)程的主要參考。

通過(guò)對(duì)KL散度和平均互信息參數(shù)的綜合，可以完整地反映光譜信息量的譜間相關(guān)特征。式(9)給出了圖像中任意2個(gè)譜段數(shù)據(jù)X和Y之間的譜間相關(guān)模型。

S(X，Y)=DKL(P‖Q)-λ·MI(X；Y)

(9)

式中：S(X，Y)為波段X和波段Y之間的譜間相關(guān)性；DKL(P‖Q)為波段X的概率分布P和波段Y的概率分布Q之間的KL散度；MI(X；Y)為波段X和波段Y之間的平均互信息；λ為權(quán)重系數(shù)，由DKL(P‖Q)和MI(X；Y)的均值的相對(duì)幅度大小決定。

容易知道，與MKL和MI相同，S也是一個(gè)M行、M列的方陣，如式(10)所示。

(10)

S中較大的數(shù)值代表較大的光譜信息差異(也可以理解為有效光譜信息，區(qū)別于冗余光譜信息)和對(duì)應(yīng)波段之間較小的光譜信息依賴程度(也可以理解為較小的光譜信息冗余)。S(i，j)表示波段i和波段j之間所有像元的綜合光譜信息差異(或可以理解為有效光譜信息)。S(i，j)的計(jì)算如式(11)所示。

S(i，j)=cKLDKL(Xi‖Xj)-cIMI(Xi；Xj)

(11)

式中：cKL和cI分別對(duì)應(yīng)KL散度矩陣MKL和平均互信息量矩陣MI的權(quán)重系數(shù)。由于KL散度和平均互信息量對(duì)于波段選擇過(guò)程的參考意義基本上同等重要，為了使對(duì)應(yīng)矩陣MKL和MI對(duì)于S的貢獻(xiàn)度基本相當(dāng)，需要參照歸一化思想，對(duì)MKL和MI按照其內(nèi)部元素的平均幅值的比例關(guān)系分別賦予相應(yīng)權(quán)重。具體來(lái)說(shuō)，即權(quán)重系數(shù)cKL和cI由MKL和MI中所有元素絕對(duì)值的均值的比值確定，如式(12)所示。

(12)

由于只關(guān)心cKL和cI之間的相對(duì)比值，而非二者的絕對(duì)數(shù)值，因此S的計(jì)算公式可以轉(zhuǎn)換為式(13)。

(13)

1.5 KLMIBS處理流程

希望利用矩陣S同時(shí)選出全部待選波段是不切實(shí)際的。本文算法主要基于如下技術(shù)思路：每次僅從剩余波段中選擇一個(gè)波段，該波段在所有待選波段中，包含了相對(duì)于已選擇波段集合最多的有效光譜信息。

對(duì)于矩陣S中的隨機(jī)一列：si={S(1，i)，S(2，i)，…，S(M，i)}，它代表波段i相對(duì)于其余全部波段的有效光譜信息(差異化光譜信息)。

S是MKL和MI的綜合，S中較大數(shù)值的元素，意味著對(duì)應(yīng)2個(gè)波段之間同時(shí)具有較大的KL散度值和較小的平均互信息，也就是說(shuō)，對(duì)應(yīng)2個(gè)波段之間的光譜信息差異較大，光譜信息冗余較小，這正是波段選擇算法的理想依據(jù)。

KLMIBS波段選擇方法流程如下：

①利用Hysime估計(jì)高光譜圖像的端元數(shù)量，將端元數(shù)量估計(jì)結(jié)果確定為波段選擇數(shù)量T。

②計(jì)算綜合信息量矩陣S?；谑?1)和式(3)計(jì)算KL散度矩陣MKL，基于式(7)和式(8)計(jì)算平均互信息矩陣MI，基于式(12)和式(13)計(jì)算S。

③選擇第一個(gè)波段。選擇x1∈[1，M]，M為波段數(shù)量，使得sx1={S(1，x1)，…，S(M，x1)}在所有的si={S(1，i)，…，S(M，i)}，i=1，2，…，M中取最大值，則序號(hào)為x1的波段即為第一個(gè)被選擇波段。

④依次選擇剩余T-1個(gè)波段。假設(shè)全部已選擇波段的序號(hào)集合為K={K1，…，KU}，其中U為已選波段數(shù)量，選擇xi∈[1，M]，i∈[2，T]，xi?K，計(jì)算SK=S(K1，xi)+…+S(KU，xi)，使得SK在所有的xi中取最大值，則xi即為第i個(gè)被選擇波段序號(hào)。按照④中步驟依次選擇剩余的T-1個(gè)波段，并將選中的每一個(gè)波段序號(hào)加入集合K，則最終生成的波段序號(hào)集合K即為算法的波段選擇結(jié)果。

2 真實(shí)數(shù)據(jù)波段選擇實(shí)驗(yàn)

真實(shí)實(shí)驗(yàn)數(shù)據(jù)為高光譜數(shù)據(jù)圖像收集實(shí)驗(yàn)儀器(hyperspectral digital imagery collection experiment，HYDICE)在美國(guó)華盛頓特區(qū)采集到的高光譜數(shù)據(jù)集，如圖1所示。該高光譜影像的行數(shù)與列數(shù)均為400。該數(shù)據(jù)原始波段數(shù)量210，波段范圍0.4～2.4 μm，波段寬度為10 nm，涵蓋了可見(jiàn)光和近紅外譜段范圍；在去除了0.9～1.4 μm范圍內(nèi)的大氣吸收波段后，剩余191個(gè)有效波段。

圖1 美國(guó)華盛頓特區(qū)HYDICE高光譜數(shù)據(jù)

根據(jù)KLMIBS波段選擇方法流程，首先利用HySime方法得到該影像子區(qū)的端元數(shù)量估計(jì)值為10，確定待選擇的波段數(shù)量為10。在191個(gè)有效波段中，KLMIBS最終選擇的10個(gè)波段的下標(biāo)集合如下：{5，13，26，39，47，56，71，82，136，191}。

圖2給出了波段選擇算法結(jié)果對(duì)應(yīng)的綜合光譜信息矩陣S的三維視圖(圖2(a))，以及10個(gè)選擇波段的光譜三維視圖(圖2(b))。如圖2(a)所示，圖中峰谷起伏程度越大，代表各波段之間的光譜差異信息(有效信息)越多，光譜相關(guān)信息(冗余信息)越少；說(shuō)明KLMIBS方法中的矩陣S能夠?qū)蓛刹ǘ沃g的光譜信息差異有效凸顯出來(lái)，同時(shí)能夠?qū)λ胁ǘ沃g的光譜有效信息整體分布情況進(jìn)行定量描述。如圖2(b)所示，10個(gè)選定波段的集合中各波段之間的光譜幅值差異十分顯著，說(shuō)明選出的波段集合包含的光譜信息冗余小，有效光譜信息量大。

圖2 KLMIBS波段選擇實(shí)驗(yàn)結(jié)果

為進(jìn)一步定量分析和驗(yàn)證KLMIBS方法的性能，將KLMIBS與前文所述的SpaBS、SNMF和ABS 3種代表性波段選擇方法的波段選擇結(jié)果在真實(shí)高光譜遙感圖像中的像元解混精度進(jìn)行對(duì)比。

像元解混算法統(tǒng)一采用最小體積約束的非負(fù)矩陣分解(minimum volume constrained NMF，MVCNMF)。MVCNMF是目前光譜解混研究中成果最多的基于凸面幾何學(xué)理論的方法中很有代表性的一種，通過(guò)將體積限制加入到NMF中，將最小二乘分析和凸面幾何結(jié)合起來(lái)，具有較高的解混精度。MVCNMF提出的代價(jià)函數(shù)包括兩部分，一部分為估量觀測(cè)數(shù)據(jù)與端元和豐度重建數(shù)據(jù)之間的近似誤差，另一部分由最小體積限制組成。

利用MVCNMF的端元光譜估計(jì)值與地面實(shí)測(cè)真值之間的光譜夾角距離[13](spectral angle distance，SAD)作為解混結(jié)果的精度評(píng)價(jià)參數(shù)，單位為弧度。SAD越小，代表像元解混精度越高，解混效果越好。

KLMIBS與SpaBS、SNMF和ABS的波段選擇結(jié)果的MVCNMF解混精度對(duì)比如表1所示。

表1 不同波段選擇算法的MVCNMF解混精度

如表1所示，與其他3種方法相比，KLMIBS對(duì)應(yīng)的SAD最小，代表其選擇結(jié)果的MVCNMF解混精度最優(yōu)，說(shuō)明該方法能夠很好地保留有效光譜信息，去除冗余光譜信息，波段選擇效果優(yōu)于多數(shù)同類算法。

3 結(jié)束語(yǔ)

論文研究結(jié)果表明，在波段選擇過(guò)程中對(duì)高光譜影像進(jìn)行光譜信息量分析，通過(guò)引入KL散度與(平均)互信息2種信息學(xué)參數(shù)，分別對(duì)波段間的光譜差異信息與光譜冗余信息進(jìn)行建模，并生成反映整幅影像有效光譜信息整體分布情況的綜合模型，能夠在數(shù)量非常有限的波段集合中，相對(duì)于其余幾種代表性波段選擇方法保留更多的光譜有效信息，有助于后續(xù)高光譜影像分析與處理。

相對(duì)于前人在高光譜影像波段選擇方面的研究成果，論文主要在以下兩個(gè)方面進(jìn)行了改進(jìn)與創(chuàng)新：一是通過(guò)賦予不同權(quán)重系數(shù)，融合2種代表性信息學(xué)參數(shù)，對(duì)高光譜影像的有效光譜信息整體分布狀況進(jìn)行建模和定量描述；二是提出的KLMIBS方法屬于非監(jiān)督方法，能夠顯著降低對(duì)先驗(yàn)知識(shí)的依賴程度，相對(duì)于部分對(duì)先驗(yàn)知識(shí)依賴程度較高的波段選擇方法，在一定程度上改善了實(shí)用性，在后續(xù)處理中，其波段選擇結(jié)果也具有令人滿意的精度。

在今后工作中，該方法的適用范圍與穩(wěn)定性還需進(jìn)一步分析與驗(yàn)證；基于KLMIBS的技術(shù)思路，其他信息學(xué)或統(tǒng)計(jì)學(xué)參數(shù)組合，以及其他權(quán)重系數(shù)計(jì)算方法，也都值得進(jìn)一步嘗試與檢驗(yàn)。