袁 博

南陽(yáng)理工學(xué)院計(jì)算機(jī)與信息工程學(xué)院，河南 南陽(yáng) 473000

空間分辨率限制及地物復(fù)雜多樣性，使得高光譜遙感圖像像元中通常包含多種地物成分，形成所謂混合像元[1-2]。像元的“解混”是實(shí)現(xiàn)地物分類(lèi)與目標(biāo)識(shí)別等高光譜遙感定量化、精細(xì)化應(yīng)用的關(guān)鍵步驟之一。非負(fù)矩陣分解(non-negative matrix factorization，NMF)適合處理高維海量數(shù)據(jù)，能夠減小對(duì)先驗(yàn)知識(shí)的依賴(lài)[3]，同時(shí)與線性混合模型的數(shù)學(xué)形式相似，因此，基于NMF的線性盲解混算法，逐漸成為高光譜解混領(lǐng)域的研究熱點(diǎn)[4]。文獻(xiàn)[5—8]表明，盲解混算法是指不依賴(lài)大量先驗(yàn)知識(shí)，只根據(jù)混合像元數(shù)據(jù)本身即可同時(shí)估計(jì)出端元光譜信息和豐度信息的一類(lèi)非監(jiān)督解混算法。

近年來(lái)，相關(guān)性分析開(kāi)始出現(xiàn)在高光譜像元解混研究中[9-12]。文獻(xiàn)[13]通過(guò)添加空間相關(guān)性約束項(xiàng)和訓(xùn)練數(shù)據(jù)的空間信息，提出了一種新的高光譜解混方法。文獻(xiàn)[14]利用先驗(yàn)概率密度函數(shù)表達(dá)兩個(gè)相鄰區(qū)域的空間相關(guān)程度，提出了一種區(qū)域相關(guān)的NMF解混算法(area-correlated spectral unmixing method based on Bayesian nonnegative matrix factorization,ACBNMF)。文獻(xiàn)[15]通過(guò)聯(lián)合最小化光譜相關(guān)度函數(shù)和NMF誤差函數(shù)，提出了一種最小化光譜相關(guān)度約束的NMF方法(minimum spectral correlation constraint NMF,MSCCNMF)。文獻(xiàn)[16]基于圖像空間特征和聯(lián)合稀疏解混思想，提出一種基于高光譜圖像空間相關(guān)性度量的多任務(wù)聯(lián)合稀疏解混方法。

上述基于相關(guān)性分析的解混算法存在如下問(wèn)題：①缺少針對(duì)高光譜遙感圖像在空間和光譜兩個(gè)維度典型相關(guān)性特征的綜合分析與利用[15-16]；②屬于監(jiān)督類(lèi)算法，需要大量難以獲取的先驗(yàn)知識(shí)，一定程度上限制了算法的適用范圍[13-14]。因此，有必要深入挖掘并充分利用高光譜遙感圖像中蘊(yùn)含的各類(lèi)相關(guān)性特征，探索其在降低先驗(yàn)知識(shí)依賴(lài)程度，改善解混算法性能方面的應(yīng)用潛力。

本文針對(duì)NMF在線性盲分解中存在的目標(biāo)函數(shù)非凸導(dǎo)致的局部極小問(wèn)題，綜合利用馬爾科夫隨機(jī)場(chǎng)(Markov random field,MRF)模型[17]和信號(hào)預(yù)測(cè)度技術(shù)[18]，提出一種基于空間與譜間相關(guān)性模型的NMF線性盲分解算法(spatial and spectral correlation model-based NMF pixel unmixing algorithm，S2CNMF)。

1 標(biāo)準(zhǔn)NMF線性盲分解模型及存在的主要問(wèn)題

1.1 NMF線性盲分解模型

1.1.1 線性光譜混合模型(linear spectral mixture model，LSMM)

LSMM的數(shù)學(xué)形式如式(1)所示

(1)

式中，X是L(L為波段數(shù)量)維的已知高光譜圖像；M是L×P(P為端元數(shù)目)的端元矩陣，每一列代表一個(gè)端元的光譜向量；向量s是像元中各端元的豐度；Ne是L維的高斯隨機(jī)噪聲。如果把高光譜圖像的全部C(假設(shè)C為像元個(gè)數(shù))個(gè)像元均考慮進(jìn)式(1)，則式(1)可擴(kuò)展為

X=MS+Ne

(2)

1.1.2NMF解混模型

已知非負(fù)矩陣X1∈Rm×n和正整數(shù)p(p

X1≈X2X3

(3)

或等價(jià)表示為

X1j≈X2X3j

(4)

(5)

(6)

1.2 NMF線性盲分解存在的主要問(wèn)題

NMF存在局部極小問(wèn)題，如果不采取對(duì)應(yīng)措施，將損害NMF線性盲分解的精度和穩(wěn)定性。

NMF的求解過(guò)程等價(jià)于通過(guò)迭代求目標(biāo)函數(shù)最小值的過(guò)程，理想情況是目標(biāo)函數(shù)為凸函數(shù)。函數(shù)凸性是數(shù)學(xué)分析中的一個(gè)重要概念，凸函數(shù)的重要性質(zhì)是：任何局部極小同時(shí)也是全局最小。只要NMF的目標(biāo)函數(shù)取得局部極小收斂，就說(shuō)明獲取了全局最小值，也即全局最優(yōu)解。文獻(xiàn)[19]給出了標(biāo)準(zhǔn)NMF目標(biāo)函數(shù)(5)的凸性判定過(guò)程與收斂性證明，結(jié)論如下：NMF目標(biāo)函數(shù)分別對(duì)于端元光譜矩陣和豐度矩陣都是凸函數(shù)，但同時(shí)對(duì)于二者是非凸函數(shù)。也就是說(shuō)，式(5)取局部極小值時(shí)，端元光譜矩陣和豐度矩陣的解并不是二者的全局最優(yōu)解。“非凸性”會(huì)使收斂結(jié)果淪為局部極小，增大解的不確定性并降低算法整體精度。

此外，高光譜圖像同時(shí)包含豐富的空間特征和光譜特征，對(duì)于混合像元的組成結(jié)構(gòu)和解混過(guò)程都具有重要影響。但是，標(biāo)準(zhǔn)NMF只將高光譜圖像視為無(wú)規(guī)律的光譜測(cè)量記錄，忽略其空間特征，也可能導(dǎo)致解混結(jié)果出現(xiàn)較大偏差。

由此可見(jiàn)，如果能夠深入挖掘高光譜圖像中蘊(yùn)含的相關(guān)性特征，建立對(duì)應(yīng)約束條件，沿著符合地物真實(shí)情況的方向?qū)MF盲分解過(guò)程進(jìn)行引導(dǎo)，將有助于消除或弱化上述問(wèn)題帶來(lái)的不利影響，改善算法性能。

2 基于MRF改進(jìn)模型的空間相關(guān)性模型

在高光譜遙感圖像中，每類(lèi)地物的空間分布往往都有各自主導(dǎo)區(qū)域，并在各自主導(dǎo)區(qū)域內(nèi)部連續(xù)分布，導(dǎo)致相鄰像元觀測(cè)數(shù)值(DN值、輻亮度或地表反射率等)之間存在近似性和相互依賴(lài)性，是地物蘊(yùn)含于混合像元及其形成過(guò)程的重要物理屬性。論文選擇以MRF模型為載體，在解混算法中顯式且定量地描述相鄰像元的空間相關(guān)特征。

MRF包含“Markov性質(zhì)”和“隨機(jī)場(chǎng)”兩個(gè)要素，可簡(jiǎn)單解釋為具有Markov性質(zhì)的隨機(jī)場(chǎng)。MRF的具體定義與主要性質(zhì)可參考文獻(xiàn)[20]。二維數(shù)字圖像可看作隨機(jī)場(chǎng)，空間相關(guān)性質(zhì)也可類(lèi)比于Markov性質(zhì)(僅相鄰時(shí)刻間的狀態(tài)相關(guān))，則空間相關(guān)特征顯著的高光譜遙感圖像可用MRF近似表示。圖像的空間能量大小與地物的空間變化頻率和幅度成正比。對(duì)于空間相關(guān)特征顯著的高光譜圖像，解混結(jié)果中空間能量越小，越可能接近地面空間能量分布的真實(shí)情況。本文以文獻(xiàn)[21]關(guān)于MRF的圖像分割模型為基礎(chǔ)，建立反映圖像空間相關(guān)特征的能量函數(shù)模型如下

E=EFeature+ERelevance

(7)

模型將能量函數(shù)分為EFeature和ERelevance兩部分，能夠反映像元本身以及與相鄰像元之間兩種特征。EFeature描述像元本身特征，只和觀測(cè)值分布有關(guān)。假設(shè)觀測(cè)圖像含有高斯噪聲，可得能量函數(shù)EFeature

(8)

ERelevance描述像元與鄰域間的相關(guān)性。由于MRF模型在求解最大后驗(yàn)概率估計(jì)的過(guò)程中使用一致平滑假設(shè)，會(huì)在圖像邊緣處造成過(guò)平滑。為克服過(guò)平滑，文獻(xiàn)[22]提出的基于DA-GMRF的無(wú)監(jiān)督圖像分割方法，利用間斷自適應(yīng)(discontinuity-adaptive,DA)思想，基于邊緣信息定義基團(tuán)勢(shì)能，即首先運(yùn)用邊緣算子對(duì)圖像提取邊緣，得到二值邊緣圖，然后定義能量函數(shù)ERelevance為

δ(fi,fi:+r)δ(bi,bi:+r)]

(9)

式中，C為圖像中像元數(shù)；D=2表示一鄰域模型，D=4表示二階鄰域模型；{ψ1,ψ2,…,ψD}是權(quán)重系數(shù)；xi為像素i的灰度值；fi為基團(tuán)中像素i的標(biāo)記值；bi為像素i處的邊緣值。當(dāng)xi=xi:+r時(shí)，δ(xi,xi:+r)=-1，否則δ(xi,xi:+r)=0；當(dāng)fi=fi:+r時(shí)，δ(fi,fi:+r)=-1，否則δ(fi,fi:+r)=1；當(dāng)bi=bi:+r時(shí)，δ(bi,bi:+r)]=1，否則δ(bi,bi:+r)]=0。式(9)考慮了平滑約束在圖像邊緣處的自適應(yīng)性，避免了邊緣處的過(guò)平滑。

文中假設(shè)所有高光譜圖像行列數(shù)均相等，行列數(shù)不相等的圖像可作為多個(gè)行列數(shù)相等子圖像的組合?；谏鲜黾僭O(shè)，下文中所有豐度矩陣S均為方陣，即存在逆矩陣。

設(shè)U為分離矩陣，即豐度矩陣S的逆，Y為端元矩陣M的估計(jì)，由式(3)可得

U=S-1

(10)

Y=AU

(11)

式(11)及后續(xù)公式中，統(tǒng)一用A表示高光譜圖像矩陣。下面以U為求解對(duì)象，通過(guò)NMF和MRF約束交替運(yùn)行的方式計(jì)算U。

為簡(jiǎn)化推導(dǎo)過(guò)程，這里僅考慮一階鄰域系統(tǒng)相關(guān)性。如果把所有權(quán)重系數(shù)ψ設(shè)置為1，用W表示特征向量的均值矩陣，同時(shí)假設(shè)特征向量中各個(gè)特征的標(biāo)準(zhǔn)差σk相等，結(jié)合式(8)和式(9)可以推出E(U)如式(12)所示

0.5Nlnσ2+0.5(UTATRU-2WAU+WWT)/σ2

(12)

能量函數(shù)E(U)取最小值時(shí)，MRF模型描述的“空間能量”最小，空間分布情況最理想，此時(shí)式(12)中E(U)關(guān)于U的偏導(dǎo)數(shù)應(yīng)為零，據(jù)此可最終推導(dǎo)出U的解如式(13)所示

U=(AT·A)-1·(W·A)T

(13)

3 基于復(fù)雜度映射的譜間相關(guān)性模型

譜間相關(guān)性源于傳感器的頻譜交疊，一般而言，光譜分辨率越高，譜間相關(guān)越強(qiáng)。分段平滑即光譜曲線在幾個(gè)特定的光譜范圍內(nèi)具有明顯的平滑特征，在光譜范圍邊界處平滑特征減弱，存在不同程度的突變。該特征反映了譜線的結(jié)構(gòu)性和冗余性，即根據(jù)某波長(zhǎng)位置的數(shù)值可大致預(yù)測(cè)出附近其他波長(zhǎng)位置的數(shù)值。真實(shí)高光譜遙感數(shù)據(jù)中，源于譜間相關(guān)的分段平滑特征非常顯著和普遍。

信號(hào)預(yù)測(cè)度也可稱(chēng)為復(fù)雜度映射(complexity pursuit)，利用信號(hào)的可預(yù)測(cè)性(predictability)度量信號(hào)復(fù)雜度，通過(guò)在混合信號(hào)的變換中尋找具有最小復(fù)雜度的信號(hào)來(lái)達(dá)到恢復(fù)原始信號(hào)的目的，其中預(yù)測(cè)度反映了利用信號(hào)本身某些部分預(yù)測(cè)其他部分的難易程度。

這里仍然用M表示高光譜數(shù)據(jù)的端元矩陣，用F(M)表示高光譜圖像端元光譜信號(hào)的預(yù)測(cè)度，則F(M)的計(jì)算公式為

(14)

(15)

最大化Vp/Up需要同時(shí)滿足兩個(gè)條件：

(1)Mlp存在一定的變化范圍；

(2)Mlp的值在時(shí)間上變化“緩慢”。

由于Mlp代表高光譜圖像的端元光譜信號(hào)，絕大多數(shù)的端元光譜都具有分段平滑特征，即同時(shí)具備變化性和平滑性，很好地滿足了上述兩個(gè)條件。因此，可利用預(yù)測(cè)度函數(shù)定量描述端元光譜的分段平滑特征，作為NMF目標(biāo)函數(shù)的新約束項(xiàng)。加入光譜預(yù)測(cè)度約束后的NMF目標(biāo)函數(shù)為

(16)

(17)

S←S·MTA/(MTMS)

(18)

式中，MT為M的轉(zhuǎn)置矩陣。

此外，為滿足ASC，S需要被歸一化，歸一化計(jì)算公式如下

(19)

4 基于空間與譜間相關(guān)性模型的NMF線性盲分解

為減小與已有約束(NMF目標(biāo)函數(shù)的歐氏距離約束或KL散度約束，非負(fù)約束、全加約束等)的相互干擾，S2CNMF算法中的空間相關(guān)約束采用了與NMF目標(biāo)函數(shù)交替運(yùn)行的獨(dú)立函數(shù)形式，即迭代解混過(guò)程由含有譜間相關(guān)約束的NMF步驟和空間相關(guān)約束步驟交替組成。

基于MRF的空間相關(guān)特征約束，以及基于預(yù)測(cè)度技術(shù)的譜間相關(guān)約束，針對(duì)NMF解混存在的問(wèn)題，分別引入了反映高光譜圖像特征的不同約束條件，但二者作用目標(biāo)和側(cè)重點(diǎn)所有不同。前者由于與NMF目標(biāo)函數(shù)獨(dú)立、交替運(yùn)行，主要側(cè)重于修正NMF流程在端元空間能量分布方面的偏差，提高解混精度；而后者存在于NMF目標(biāo)函數(shù)內(nèi)部，主要側(cè)重于修正NMF目標(biāo)函數(shù)的非凸性，確保并加快NMF目標(biāo)函數(shù)的收斂。

可得S2CNMF算法步驟如下：

(1) 利用基于最小誤差的高光譜信號(hào)識(shí)別法(hyperspectral signal identification by minimum error，HySime)估算端元數(shù)量P；

(2) 初始化端元矩陣M和豐度矩陣S；

(3) 含有譜間相關(guān)約束的NMF迭代：根據(jù)更新規(guī)則式(17)和式(18)分別計(jì)算端元矩陣M和豐度矩陣S的迭代結(jié)果，并利用M和S計(jì)算分離矩陣U的迭代結(jié)果；

(4) 空間相關(guān)約束：歸一化分離矩陣U的每一列，同時(shí)估計(jì)像元特征向量的均值矩陣W，然后計(jì)算U；

(5) 重復(fù)步驟(3)和(4)，繼續(xù)迭代，直到各自的停止準(zhǔn)則同時(shí)滿足，得到一個(gè)估計(jì)的成分，將其轉(zhuǎn)換為矩陣，即獲得一個(gè)端元分布。繼續(xù)迭代，直至滿足閾值條件(預(yù)先設(shè)置的一個(gè)很小的正數(shù)如10-4等，作為迭代停止條件)，得到端元分布的估計(jì)結(jié)果。

其中，端元數(shù)量估計(jì)算法HySime是一種估計(jì)高光譜信號(hào)子空間的方法，計(jì)算過(guò)程復(fù)雜，但不需任何參數(shù)，具有自適應(yīng)性，估計(jì)準(zhǔn)確度較高。文獻(xiàn)[23]給出了該算法的原理與實(shí)現(xiàn)過(guò)程。

算法步驟(3)中NMF迭代的目標(biāo)函數(shù)收斂性證明可參見(jiàn)文獻(xiàn)[24]；步驟(3)收斂保證了分離矩陣U的迭代結(jié)果收斂，則算法步驟(4)中對(duì)U的每一列進(jìn)行歸一化得到的特征向量均值矩陣W也是收斂的；容易證明式(13)關(guān)于W收斂，因此步驟(4)中的目標(biāo)函數(shù)收斂。由于步驟(3)和步驟(4)共同組成了一個(gè)帶有初始邊界值的并行交替迭代法，該算法可抽象為一種Dirichlet-Neuman(D-N)交替迭代法，結(jié)合文獻(xiàn)[25]對(duì)于D-N交替迭代法及其收斂性的分析結(jié)果，可以證明S2CNMF算法的收斂性。

需要進(jìn)一步說(shuō)明的是，本文算法主要針對(duì)空間相關(guān)程度較高，或者空間相關(guān)程度一般但譜間相關(guān)程度較高的高光譜遙感場(chǎng)景；反之也成立。也就是說(shuō)，只有空間、譜間兩種相關(guān)特征同時(shí)不顯著時(shí)，才會(huì)嚴(yán)重影響算法性能，多數(shù)自然場(chǎng)景下的真實(shí)高光譜數(shù)據(jù)一般不會(huì)出現(xiàn)上述情況。

5 試驗(yàn)與分析

5.1 試驗(yàn)1

試驗(yàn)1中的數(shù)據(jù)來(lái)自超光譜數(shù)字圖像收集實(shí)驗(yàn)儀器(hyperspectral digital imagery collection experiment,HYDICE)在美國(guó)華盛頓特區(qū)采集的高光譜數(shù)據(jù)集，如圖1所示，大小為400×400像元。該數(shù)據(jù)原始波段數(shù)量210，波段范圍400～2400 nm，波段寬度為10 nm，涵蓋了可見(jiàn)光和近紅外譜段范圍；在去除900～1400 nm的大氣吸收波段后，該數(shù)據(jù)剩余191波段。該數(shù)據(jù)可在美國(guó)地質(zhì)調(diào)查局(United States Geological Survey)網(wǎng)站下載(https:∥speclab.cr.usgs.gov/spectral.lib04/lib04-HYDICE.html)，同時(shí)可下載的還有該地區(qū)的實(shí)地調(diào)查報(bào)告及配套實(shí)測(cè)數(shù)據(jù)，據(jù)此可知該地區(qū)圖像中主要含有水體、裸土和樹(shù)木3種地物。

圖1 美國(guó)華盛頓特區(qū)HYDICE數(shù)據(jù)Fig.1 HYDICE data of Washington D.C.,USA

鄰接規(guī)則選擇queen(具有公共邊界或公共頂點(diǎn)的像元為相鄰像元)，距離計(jì)算方法選擇歐氏距離，鄰接距離d設(shè)置為3(空間距離不大于3的像元互為相鄰像元)，根據(jù)式(20)計(jì)算整幅圖像的全局Moran’I空間自相關(guān)指數(shù)，結(jié)果為0.502 6，說(shuō)明各像元間的空間相關(guān)特征顯著，這也是大部分真實(shí)高光譜遙感圖像的共同特點(diǎn)

(20)

式中，n為樣本量，即空間位置的個(gè)數(shù)；xi、xj是空間位置i和j的觀察值；wij表示空間位置i和j的鄰近關(guān)系，當(dāng)i和j為鄰近的空間位置時(shí)，wij=1；反之，wij=0。全局Moran指數(shù)I的取值范圍為[-1，1]。

利用HySime估算端元數(shù)量P，結(jié)果為3，這也驗(yàn)證了從地面實(shí)測(cè)報(bào)告得出的主要含有水體、裸土和樹(shù)木3種地物的結(jié)論。因此，將解混試驗(yàn)的端元數(shù)量設(shè)置為3。由于HYDICE為機(jī)載光譜成像儀，空間分辨率高，圖像中每種端元都含有大量純像元。因此，通過(guò)在原圖像中人工選擇參考點(diǎn)的方式收集每種端元(地物)的光譜作為參考值，利用全約束最小二乘法計(jì)算端元的豐度參考值。

通過(guò)均值法重采樣使圖像的空間分辨率降低為原來(lái)的0.2倍?？臻g分辨率的大幅降低，形成了大量混合像元，便于驗(yàn)證S2CNMF算法。由于在地物范圍確定的前提下，豐度分解結(jié)果的尺度是和空間分辨率大小嚴(yán)格對(duì)應(yīng)的，同樣利用均值法重采樣技術(shù)把之前求出的端元豐度參考值進(jìn)行聚合，得到行、列元素?cái)?shù)量均為原始圖像0.2倍的新圖像豐度參考值，作為解混結(jié)果精度分析中的近似真值。降分辨率后的圖像如圖2所示，降分辨率后圖像的豐度參考值如圖3所示。圖3中純白色代表端元在該像元內(nèi)部面積比例為1(100%)，純黑色代表端元在該像元內(nèi)部面積比例為0，其余各階灰度分別對(duì)應(yīng)0～1的不同比例。

圖2 降分辨率后的美國(guó)華盛頓特區(qū)HYDICE數(shù)據(jù)Fig.2 Resolution descended HYDICE data of Washington D.C.,USA

圖3 降分辨率后的美國(guó)華盛頓特區(qū)HYDICE數(shù)據(jù)豐度參考值Fig.3 Abundance reference of resolution descended HYDICE data of Washington D.C.,USA

試驗(yàn)選擇如下比較方法：區(qū)域相關(guān)的NMF解混算法(ACBNMF)、最小化光譜相關(guān)度約束的NMF方法(MSCCNMF)，這兩種都是基于相關(guān)性分析的NMF解混算法。同時(shí)選擇最小體積約束的非負(fù)矩陣分解(MVCNMF)作為比較方法，以評(píng)價(jià)提出的約束條件的有效性。MVCNMF目標(biāo)函數(shù)包括兩部分：①估計(jì)觀測(cè)數(shù)據(jù)與端元和豐度重建數(shù)據(jù)之間的近似誤差；②最小體積限制。

圖4—圖7分別列出了該數(shù)據(jù)的S2CNMF、ACBNMF、MSCCNMF和MVCNMF解混豐度估計(jì)結(jié)果。論文利用式(21)計(jì)算光譜角距離SAD，代表端元光譜分解精度

SAD(EUnmixing,EReference)=

(21)

式中，EUnmixing表示某種端元在光譜解混結(jié)果中的端元光譜估計(jì)值；EReference表示某種端元的參考真值。

圖4 降分辨率后美國(guó)華盛頓特區(qū)HYDICE數(shù)據(jù)S2CNMF豐度估計(jì)Fig.4 S2CNMF’s abundance estimation of resolution descended HYDICE data of Washington D.C.,USA

圖5 降分辨率后美國(guó)華盛頓特區(qū)HYDICE數(shù)據(jù)ACBNMF豐度估計(jì)Fig.5 ACBNMF’s abundance estimation of resolution descended HYDICE data of Washington D.C.,USA

圖6 降分辨率后美國(guó)華盛頓特區(qū)HYDICE數(shù)據(jù)MSCCNMF豐度估計(jì)Fig.6 MSCCNMF’s abundance estimation of resolution descended HYDICE data of Washington D.C.,USA

圖7 降分辨率后美國(guó)華盛頓特區(qū)HYDICE數(shù)據(jù)MVCNMF豐度估計(jì)Fig.7 MVCNMF’s abundance estimation of resolution descended HYDICE data of Washington D.C.,USA

(22)

式中，L為高光譜圖像波段數(shù)量；C為全部像元數(shù)量。

端元光譜和豐度解混的精度分析結(jié)果分別見(jiàn)表1和表2，其中粗體代表每一列的精度最高項(xiàng)。可以看出：4種解混算法均能以較為理想的精度有效分解出3種主要地物類(lèi)型，其中本文提出的S2CNMF算法綜合性能最好，以端元光譜分解結(jié)果為準(zhǔn)，S2CNMF比其他3種算法的分解精度分別提高了12.8%、14.7%和15.0%；以豐度估計(jì)結(jié)果為準(zhǔn)，S2CNMF比其他3種算法的分解精度分別提高了11.4%、13.7%和13.4%。上述結(jié)果表明，S2CNMF算法中空間和譜間相關(guān)性特征約束項(xiàng)的聯(lián)合，對(duì)于克服NMF的局部極小問(wèn)題，進(jìn)一步提高解混精度具有積極意義和較理想能力。

表1 端元光譜分解結(jié)果精度

表2 豐度分布估計(jì)結(jié)果精度

此外，S2CNMF算法屬于非監(jiān)督算法，不需要像頂點(diǎn)成分分析(vertex component analysis,VCA)等基于純像元的解混算法一樣假設(shè)每種地物至少存在一個(gè)純像元，不需要像文獻(xiàn)[9—10]中提出的算法一樣需要預(yù)知圖像空間特征和地物分布規(guī)律，也不需要預(yù)知每種地物的具體種類(lèi)和預(yù)先獲取每種地物的光譜數(shù)據(jù)，只需要使用Hysime等方法預(yù)先估計(jì)出地物(端元)數(shù)量即可，因此對(duì)于解混先驗(yàn)知識(shí)的依賴(lài)程度較小，適用范圍較廣。

最后，在表3給出了上述4種解混算法的運(yùn)行時(shí)間。可以看出，4種算法的運(yùn)算效率從高到低依次為：S2CNMF、ACBNMF、MSCCNMF和MVCNMF，其中S2CNMF和ACBNMF的運(yùn)算效率比較接近，明顯優(yōu)于MSCCNMF和MVCNMF。究其原因，MSCCNMF為了度量多條光譜的相關(guān)程度，需要計(jì)算光譜相關(guān)度函數(shù)，而該函數(shù)是由每?jī)蓷l光譜協(xié)方差函數(shù)的平方和所組成的，由于每次迭代都需要計(jì)算光譜相關(guān)度函數(shù)，故算法整體計(jì)算時(shí)間較長(zhǎng)；MVCNMF算法中的目標(biāo)函數(shù)最小化是一個(gè)組合優(yōu)化問(wèn)題，MVCNMF采用交替非負(fù)最小二乘法(alternating nonnegative least squares)結(jié)合投影梯度學(xué)習(xí)法(projected gradient learning)求解該問(wèn)題并實(shí)現(xiàn)非負(fù)約束，對(duì)應(yīng)的迭代更新準(zhǔn)則在求目標(biāo)函數(shù)對(duì)豐度的偏導(dǎo)數(shù)過(guò)程中，涉及大量高階矩陣的行列式計(jì)算，需要耗費(fèi)大量時(shí)間，導(dǎo)致算法迭代速度較慢，運(yùn)行時(shí)通常耗時(shí)較長(zhǎng)。

表3 算法運(yùn)行時(shí)間對(duì)比

5.2 試驗(yàn)2

試驗(yàn)1中，高光譜成像場(chǎng)景以成片分布的各類(lèi)型農(nóng)田為主，可以預(yù)見(jiàn)其空間相關(guān)性較高。為進(jìn)一步驗(yàn)證算法對(duì)空間相關(guān)特征不顯著場(chǎng)景的解混效果，試驗(yàn)2中的數(shù)據(jù)選擇獲取于1995年7月美國(guó)Nevada州Cuprite采礦區(qū)的AVIRIS高光譜數(shù)據(jù)，如圖8所示。圖8中的A、B、C、K、M等字樣分別代表該地區(qū)5種廣泛分布的礦物Alunite、Buddingtointe、Calcite、Kaolinite和Muscovite的大致分布位置。該圖像大小為400列，350行，空間分辨率為20 m，波長(zhǎng)范圍為0.4～2.5 μm，光譜分辨率為10 nm，共224個(gè)波段，試驗(yàn)中剔除了存在水汽吸收和信噪比較低的波段，使用剩余的188個(gè)波段。該地區(qū)位于美國(guó)Nevada州南部，地表多為裸露礦物，基本無(wú)植被覆蓋。相對(duì)于試驗(yàn)1中的數(shù)據(jù)，該數(shù)據(jù)的空間相關(guān)特征不明顯，幾種主要礦物空間分布的隨機(jī)性較大。

由于AVIRIS為機(jī)載光譜成像儀，空間分辨率(20 m)相對(duì)較高，可合理假設(shè)：對(duì)于每一種端元(地物)，都能找到它具有較高豐度值的混合像元，這些混合像元的光譜數(shù)據(jù)可作為其中居主導(dǎo)地位的端元的光譜參考值。因此，通過(guò)在原圖像中人工選擇參考點(diǎn)的方式，收集每種端元(地物)的光譜作為參考值，利用全約束最小二乘法計(jì)算端元的豐度參考值。

圖8 Cuprite采礦區(qū)AVIRIS高光譜數(shù)據(jù)Fig.8 AVIRIS hyperspectral data of Cuprite mining field

重復(fù)試驗(yàn)1中的試驗(yàn)環(huán)節(jié)，對(duì)比和驗(yàn)證算法針對(duì)不同場(chǎng)景時(shí)的性能表現(xiàn)。需要說(shuō)明的是，試驗(yàn)2只針對(duì)上述5種主要礦物類(lèi)型進(jìn)行解混試驗(yàn)，對(duì)于其余小目標(biāo)代表的地物類(lèi)型的假設(shè)和處理與數(shù)據(jù)1完全相同，具體情況不再贅述。

圖9為本文算法S2CNMF對(duì)應(yīng)的Cuprite采礦區(qū)AVIRIS數(shù)據(jù)豐度估計(jì)結(jié)果。由于篇幅限制，MVCNMF等其他3種算法的結(jié)果圖不再一一列出。

圖9 Cuprite采礦區(qū)AVIRIS數(shù)據(jù)S2CNMF豐度估計(jì)Fig.9 S2CNMF’s abundance results of AVIRIS hyperspectral data of Cuprite mining field

表4和表5分別給出了4種解混算法的端元光譜和豐度解混的精度分析結(jié)果。其中粗體代表每一列的精度最高項(xiàng)。可以看出：4種解混算法均能有效分解出5種主要礦物類(lèi)型，其中本文提出的S2CNMF算法性能仍然最好，但與其余3種算法均比較接近。以端元光譜分解結(jié)果為準(zhǔn)，S2CNMF比其他3種算法的分解精度分別提高了8.85%、10.8%和9.34%；以豐度估計(jì)結(jié)果為準(zhǔn)，S2CNMF比其他3種算法的分解精度分別提高了6.27%、6.85%和8.03%。其余3種算法精度基本相當(dāng)。

表4 Cuprite采礦區(qū)數(shù)據(jù)的端元光譜分解結(jié)果精度

Tab.4 Precision of end member estimation of Cuprite mining field

解混方法S2CNMFACBNMFMSCCNMFMVCNMFAlunite0.21220.22570.23750.2401Buddingtointe0.19870.21380.22080.2267Calcite0.19310.22950.21940.2103Kaolinite0.19850.20560.22360.2238Muscovite0.21240.23790.23550.2175平均0.20280.22250.22740.2237

表5 Cuprite采礦區(qū)數(shù)據(jù)的豐度分布估計(jì)結(jié)果精度

Tab.5 Precision of abundance estimation of Cuprite mining field

解混方法S2CNMFACBNMFMSCCNMFMVCNMFAlunite0.21220.21930.22360.2412Buddingtointe0.19980.22250.21940.2179Calcite0.22010.23570.23080.2337Kaolinite0.21520.23190.22770.2358Muscovite0.20670.21530.23010.2175平均0.21080.22490.22630.2292

上述結(jié)果表明，本文算法在試驗(yàn)數(shù)據(jù)空間相關(guān)程度顯著減小的情況下，相對(duì)于其他幾種NMF代表性算法，解混精度的優(yōu)勢(shì)略有減小。分析其原因，試驗(yàn)數(shù)據(jù)空間相關(guān)程度的降低，會(huì)不可避免地導(dǎo)致算法基于空間相關(guān)特征的約束部分出現(xiàn)一定程度的性能下降；但由于算法還包括了譜段間的相關(guān)性約束，能夠在較大程度上彌補(bǔ)空間相關(guān)約束部分解混性能降低帶來(lái)的負(fù)面影響，在整體上仍保持較為理想的解混性能。

6 結(jié)論與展望

論文結(jié)果表明，在NMF線性盲分解過(guò)程中分別利用MRF模型與復(fù)雜度映射，引入高光譜圖像的空間與譜間相關(guān)特征約束，能夠顯著降低目標(biāo)函數(shù)非凸問(wèn)題帶來(lái)的不利影響，相較于ACBNMF、MSCCNMF和MVCNMF等幾種代表性NMF擴(kuò)展算法，能夠進(jìn)一步改善解混精度。

相對(duì)于前人的研究成果，本文主要從3個(gè)方面進(jìn)行了改進(jìn)：①針對(duì)影響NMF線性盲分解性能的“局部極小”問(wèn)題，在解混模型中同時(shí)構(gòu)建空間與譜間相關(guān)性特征約束，增強(qiáng)了NMF解混模型對(duì)高光譜遙感數(shù)據(jù)特征的表述能力；②引入空間相關(guān)特征約束時(shí)，放棄了添加目標(biāo)函數(shù)新約束項(xiàng)的傳統(tǒng)做法，而是采用與NMF目標(biāo)函數(shù)并行迭代的新形式，有助于簡(jiǎn)化算法模型，避免過(guò)多內(nèi)部函數(shù)項(xiàng)之間的相互影響；③提出的屬于非監(jiān)督算法，能夠降低對(duì)先驗(yàn)知識(shí)的依賴(lài)，相對(duì)于更為依賴(lài)先驗(yàn)知識(shí)的現(xiàn)有監(jiān)督類(lèi)算法，實(shí)用性進(jìn)一步改善。

最后應(yīng)當(dāng)指出，盡管相關(guān)性特征在高光譜圖像中普遍存在，但并非全部都同時(shí)具有顯著的空間與譜間相關(guān)特征。如何進(jìn)一步拓展適用范圍，穩(wěn)定算法性能，將是下一步研究工作的主要內(nèi)容。