基于最小體積約束的非負(fù)矩陣分解模型的高光譜解混算法探究

余肖玲，黃光鑫

(成都理工大學(xué) 管理科學(xué)學(xué)院，四川 成都 610059)

0 引 言

在測量光譜圖像的過程中，受自然條件和環(huán)境等因素的影響，大多數(shù)地面的真實(shí)測量數(shù)據(jù)幾乎不可能獲得，從而導(dǎo)致圖像中的端元無法自動獲取.對于這種無法從圖像中提取端元的情況，若能實(shí)現(xiàn)對高光譜遙感數(shù)據(jù)進(jìn)行數(shù)據(jù)驅(qū)動下的盲分解，則非監(jiān)督的混合像元分解就顯得重要.非監(jiān)督的混合像元分解，能做到的就是在完全不了解端元信息的情況下，直接以遙感圖像本身為目標(biāo)，根據(jù)混合像元的光譜模型及約束條件等信息，把混合像元分解為端元光譜矩陣及其在像元中所占的比例矩陣，即系數(shù)矩陣或豐度矩陣［1］.在現(xiàn)有資料研究中，混合像元模型主要包括2 類:線性光譜混合模型和非線性光譜混合模型.近年來，隨著對光譜混合特性的進(jìn)一步研究，線性光譜混合模型越來越引起更多學(xué)者的重視［2－9］.對此，本研究首先對線性光譜混合模型進(jìn)行描述，然后在線性混合模型的基礎(chǔ)上，再對混合像元分解算法進(jìn)行討論.

1 線性光譜混合模型

線性光譜混合模型［10］的表達(dá)式為，

式中，x 為m(m 為影像波段數(shù))維混合像元光譜，是遙感影像圖已獲得的實(shí)際測量值，A 為m ×p(p 為端元數(shù)目)的端元矩陣，其中每一列代表一個端元的光譜向量，向量s 為像元中各端元的豐度，ε 為m 維模型誤差或高斯隨機(jī)噪聲.

若把遙感圖像中的全部像元均考慮在內(nèi)，式(1)可寫為，

式中，X ∈Rm×n，其列向量表示各個像元的光譜;S∈Rp×n構(gòu)成系數(shù)矩陣或豐度矩陣.

在線性光譜混合模型中，系數(shù)矩陣S 需滿足2個限制條件［6］，即:1)非負(fù)性.

2)豐度和為1.

2 非負(fù)矩陣分解模型

約束非負(fù)矩陣分解(Constrained Nonnegative Matrix Factorization，CNMF)算法通過加入非負(fù)限制可將一個非負(fù)矩陣分解(Nonnegative Matrix Factorization，NMF)分為2 個非負(fù)矩陣A 和S 的乘積，

X = A × S

其中，A 為m × p 階矩陣，S 為p × n 階矩陣，p ∈min(m，n).

目標(biāo)函數(shù)式為，

由式(5)可知，NMF 算法理論上就是一個最優(yōu)化問題，其目的就是目標(biāo)函數(shù)在基于歐氏距離的要求下尋找到最小值，并滿足A 和S 滿足非負(fù)性的約束［3］.式(5)中，σ 和τ 是規(guī)則化參數(shù)，σ‖A‖2/2 和τ‖S‖2/2 為引入的平滑性約束，目的是調(diào)節(jié)逼近誤差和平衡約束的消長關(guān)系.在尋求最小值的過程中，A 和S 遵循的迭代準(zhǔn)則為，

在式(6)的指引下，式(5)是單調(diào)非增的.當(dāng)A和S 值穩(wěn)定時，目標(biāo)函數(shù)值不再發(fā)生變化.在高光譜解混中，可以有意識地把X 設(shè)想成m 維混合像元光譜，利用上述迭代規(guī)則尋求滿足非負(fù)性的端元矩陣和豐度矩陣.需要注意的是，為了滿足式(4)的成立，每一步迭代后仍需對S 按，

做歸一化處理［5］.

3 基于最小體積約束的非負(fù)矩陣分解方法

基于最小體積約束的非負(fù)矩陣分解方法(Nonnegative Matrix Factorization Method Based on Minimum Volume Constraint，MVC-NMF)在最小二乘和凸面幾何結(jié)合的基礎(chǔ)上，把體積約束附加到NMF 中，其相應(yīng)的目標(biāo)函數(shù)式為，

用來計算單體體積的懲罰項.φ ∈R.

把交替最小平方法引入到上述優(yōu)化問題中，就可以利用，

把初始問題分解為2 個小問題.在計算過程中，保持矩陣A 和S 任何一個不變來交替更新另一個，并采用梯度投影學(xué)習(xí)規(guī)則強(qiáng)化每個小問題中的矩陣非負(fù)限制.若出現(xiàn)不滿足約束限制的負(fù)值，則在保持非負(fù)變量不變的同時，利用函數(shù)max(0，x)將負(fù)值設(shè)定為0，則有，

將目標(biāo)函數(shù)f(A，S)對矩陣S 求偏導(dǎo)數(shù)，

Sf(A，S)= AT(AS－X)，

令，

則目標(biāo)函數(shù)對矩陣A 的偏導(dǎo)數(shù)為，

由此，MVC-NMF 算法的具體步驟為:

①輸入混合光譜數(shù)據(jù)矩陣X ∈Rm×n(m 為波段數(shù)，n 為像元數(shù));

②用虛擬維度［9］的方法估計端元數(shù)P;

③從X 中任意選取P 列充當(dāng)初始源矩陣A ∈Rm×P;

④生成初始系數(shù)矩陣S ∈RP×n，并根據(jù)式(7)對S 進(jìn)行歸一化處理;

⑤設(shè)定最大迭代次數(shù)MaxIter;

⑥根據(jù)式(10)分別更新A 和S 并根據(jù)式(7)對矩陣S 再做歸一化處理;

⑦設(shè)置最大迭代次數(shù)MaxIter;

⑧若MaxIter 已達(dá)到最大值或目標(biāo)值不再收斂，則算法結(jié)束;否則重復(fù)步驟⑥.

4 實(shí) 驗

為了檢驗MVC-NMF 算法對混合像元分解的性能，在實(shí)驗中，采用USGS 光譜庫中的五條光譜來構(gòu)建模擬數(shù)據(jù)，添加噪聲所用的信噪比為30 dB，光譜曲線如圖1 所示，同時，還將本算法與平滑LO 模約束的MVC-NMF(MVC-NMF with Smoothed LO norm Constraint，SC-MVC-NMF)算法的實(shí)驗結(jié)果進(jìn)行對比.

1)端元數(shù)不同時的實(shí)驗對比.為了對比在不同端元數(shù)目情況下2 種算法的結(jié)果，實(shí)驗中給模擬圖像添加了30 dB 的高斯噪聲且采用的最大迭代次數(shù)為60 次.MVC-NMF 算法和SC-MVC-NMF 算法均采用VCA+FCLS 進(jìn)行初始化.將端元數(shù)設(shè)定為5、6、7時，2 種算法提取結(jié)果的運(yùn)行時間及相似性測度的結(jié)果如表1 所示.

圖1 構(gòu)建模擬圖像的5 條光譜

表1 端元數(shù)不同時的實(shí)驗對比

由表1 可知，在端元數(shù)目不同的情況下，MVCNMF 算法的精度略低于SC-MVC-NMF 算法，但其運(yùn)行時間明顯較快.

2)噪聲的穩(wěn)健性對比.為了對比2 種算法的抗噪性，利用上述5 條光譜構(gòu)建模擬圖像，并在模擬數(shù)據(jù)中分別添加SNR =10、20、30、40 dB 的高斯噪聲且設(shè)定最大迭代次數(shù)為60.MVC-NMF 算法和SCMVC-NMF 算法均利用VCA +FCLS 進(jìn)行初始化.在不同的信噪比下，MVC-NMF 和SC-MVC-NMF 算法獲取的端元光譜、豐度矩陣的相似性測度如表2 所示.

表2 噪聲的穩(wěn)健性對比

從表2 可以看出，在信噪比減小的情況下，2 種算法的結(jié)果都逐漸變差，盡管SC-MVC-NMF 算法運(yùn)行時間略有增加，但其混合像元分解的結(jié)果始終優(yōu)于本研究算法提取的結(jié)果.

3)實(shí)驗結(jié)果表明，噪聲的穩(wěn)健性方面，在端元數(shù)不同的情況下，本研究算法收斂速度較快，但其混合像元解混的結(jié)果略差，因此，附加條件下MVCNMF 的改進(jìn)成為下一步研究的重點(diǎn).

5 結(jié) 論

本研究提出的MVC-NMF 算法對于高光譜圖像混合像元分解問題有較大的實(shí)用價值，然而MVCNMF 算法在計算過程中表現(xiàn)出來的收斂效果略差的缺點(diǎn)，卻影響了高光譜混合像元分解的效率，且很難找到最優(yōu)解，這些劣勢限制了非負(fù)矩陣分解模型的應(yīng)用.因此，在以后的研究過程中，應(yīng)進(jìn)一步探索更為快速高效的高光譜混合像元分解的NMF 算法.

［1］賈森.非監(jiān)督的高光譜圖像解混技術(shù)研究［D］.杭州:浙江大學(xué)，2006.

［2］Maurich D C.Minimum volume transforms for remotely sensed data［J］.IEEE Transactions on Geoscience and Remote Sensing，1994，32(3):542－552.

［3］Lin C J.Projected gradient methods for non—negative matrix factorization［J］.Neural Computaion，2007，19(10):2756－2779.

［4］Pauca U，Piper J，Plemmons R.Nonnegative matrix factorization for spectral data analysis［J］.Linear Algebra and Application，2006，416(1):29－47.

［5］薛綺.基于線性混合模型的高光譜圖像端元提取方法研究［D］.長沙:國防科學(xué)技術(shù)大學(xué)，2004.

［6］Heinz D C，Chang C I.Fully constrained least squares linear spectral mixture analysis method for material quantification in Hyperspectral imagery［J］.IEEE Geoscience and Remote Sensing Letters，2001，39(3):529－545.

［7］Miao L，Qi H.Endmember extraction from highly mixed data using mininlum volun ℃constrained nonnegative malrix factorization［J］.IEEE Transactions on Geoscience and Remote Sensing，2007，45(3):765－777.

［8］Chang C I，Du Q.Estimation of number of spectrally distinct signal sources in hyperspectral imagery［J］.IEEE Transcation on Geoscience and Remote Sensing，2004，42(3):608－619.

［9］劉雪松.基于非負(fù)矩陣分解的高光譜遙感圖像混合像元分解研究［D］.上海:復(fù)旦大學(xué)，2011.

［10］趙春暉，成寶芝，楊偉超.利用約束非負(fù)矩陣分解的高光譜解混算法［J］.哈爾濱工程大學(xué)學(xué)報，2012，32(3):96－99.

［11］施蓓琦，劉春，孫偉偉，等.應(yīng)用稀疏非負(fù)矩陣分解聚類實(shí)現(xiàn)高光譜影像波段的優(yōu)化選擇［J］.測繪學(xué)報，2013，46(3):351－358.