利用獨(dú)立性約束非負(fù)矩陣分解的高光譜解混算法

楊秀坤，王東輝

(哈爾濱工程大學(xué) 信息與通信工程學(xué)院，黑龍江 哈爾濱 150001)

高光譜圖像的主要特點(diǎn)是，成像光譜儀將圖像維和光譜維信息融合為一體，其探測(cè)波段幾乎以納米和亞微米級(jí)的光譜分辨率覆蓋了可見光、近紅外、短波紅外、中波紅外和熱紅外等全部光譜區(qū)[1]。由于空間分辨率的限制以及地物的復(fù)雜多樣性，一個(gè)像元內(nèi)往往會(huì)包含多種地物類型，稱為混合像元。混合像元分解是高光譜研究的重要課題。為了提高對(duì)真實(shí)地表覆蓋的描述準(zhǔn)確性，需要對(duì)混合像元進(jìn)行分解，計(jì)算一種地物類型(端元)在該像元中所占的比例(豐度)。由于非負(fù)矩陣分解NMF(nonnegative matrix factorization)[2]不需要假設(shè)原始數(shù)據(jù)的概率密度分布等先驗(yàn)信息，利用非負(fù)約束獲取數(shù)據(jù)，因此，基于統(tǒng)計(jì)學(xué)方法的NMF近年在高光譜解混中受到關(guān)注。因?yàn)槟繕?biāo)函數(shù)具有非凸性，NMF存在局部極小值，采用懲罰函數(shù)法，需要根據(jù)應(yīng)用領(lǐng)域的特點(diǎn)增加輔助約束條件克服該缺點(diǎn)。目前，已經(jīng)提出了平滑性[3]、最小體積(MVCNMF)[4]、分段平滑和豐度稀疏(PSSCNMF)[5]、最小豐度協(xié)方差[6]、端元光譜差距[7]等約束條件。獨(dú)立成分分析ICA(independent component analysis)[8]基于信號(hào)高階統(tǒng)計(jì)特性，將數(shù)據(jù)線性分解為統(tǒng)計(jì)獨(dú)立的成分。混合像元中各端元光譜相互獨(dú)立，因此，可將度量端元光譜數(shù)學(xué)期望的非高斯性的四階累積量和負(fù)熵約束應(yīng)用于非負(fù)矩陣分解，提出了一種端元光譜獨(dú)立性約束的NMF方法，I-NMF(independent nonnegative matrix factorization)，適用于沒有純像元的情形。

1 基于端元光譜獨(dú)立性約束非負(fù)矩陣分解

1.1 線性光譜混合模型

線性光譜混合模型(linear spectral mixture model, LSMM)是混合像元分解的常用方法。這種模型計(jì)算簡單、應(yīng)用最為廣泛。在線性模型中，混合光譜是端元光譜與其比例的線性組合。設(shè)X為像元光譜矢量，S=[s1s2…sN]為端元光譜矩陣，A=[a1a2…aN]T是N維矢量，其各分量元素為對(duì)應(yīng)端元豐度，N為隨機(jī)噪聲。線性光譜混合模型為

(1)

1.2 非負(fù)矩陣分解

NMF算法通過最小化歐式距離目標(biāo)函數(shù)，在已知X的情況下解得S和A的最優(yōu)解：

(2)

迭代公式為

(3)

基于NMF的光譜解混算法屬于端元生成(endmember genera-tion algorithm, EGA)算法，不需要假定存在純像元，在提取端元的同時(shí)獲得相應(yīng)的豐度。

因NMF分解算法式(1)具有非負(fù)性限制，其主要缺陷是容易陷入局部最小值，并且分解結(jié)果具有幅度上的不確定性[7]。I-NMF算法以端元光譜之間的獨(dú)立性為約束條件，使得目標(biāo)函數(shù)綜合了圖像的分解誤差和端元光譜的影響。

1.3 端元光譜數(shù)學(xué)期望的四階累積量約束

統(tǒng)計(jì)論的中心極限定理規(guī)定，一組獨(dú)立隨機(jī)變量和的分布比任何信號(hào)更接近高斯分布，可將非高斯性作為隨機(jī)信號(hào)相互獨(dú)立的度量。四階累積量是一種典型的非高斯性測(cè)度，在理論和計(jì)算上簡單，易于實(shí)現(xiàn)[9]。非高斯隨機(jī)變量有正有負(fù)，采用端元光譜數(shù)學(xué)期望的四階累積量絕對(duì)值作為約束條件:

(4)

(5)

I表示其分量都為1的向量，即IN=[1 1 … 1]∈UN。由于X需進(jìn)行中心化和白化預(yù)處理，則si的數(shù)學(xué)期望為0，均方差為1。迭代中用瞬時(shí)值代替數(shù)學(xué)期望計(jì)算，該約束應(yīng)被最大化，將該約束求導(dǎo)數(shù)，整理可得到

(6)

(7)

1.4 端元光譜數(shù)學(xué)期望負(fù)熵約束

僅考慮式(5)作為非負(fù)矩陣分解的約束條件不能獲得理想解混精度，還需要增加端元光譜的負(fù)熵約束。信息論中，各信息分量的獨(dú)立性可用互信息表示，互信息最小等同于負(fù)熵最大。端元光譜負(fù)熵可用一種非線性函數(shù)組合逼近[10]:

(8)

(9)

1.5 I-NMF算法描述

增廣式(4)，I-NMF的目標(biāo)函數(shù)為

(10)

式中：λ1、λ2為權(quán)值。

約束條件(懲罰函數(shù))J1(S)、J2(S)保證端元光譜具有較大的獨(dú)立性，使其接近真實(shí)的端元光譜。目標(biāo)函數(shù)綜合考慮了最小近似誤差和端元光譜懲罰函數(shù)的限制條件。權(quán)值λ1、λ2調(diào)整分解的重構(gòu)誤差與懲罰函數(shù)J1(S)、J2(S)之間的關(guān)系。式(10)等號(hào)后的第1項(xiàng)保證盡可能小的分解誤差，而懲罰函數(shù)保證所提取的端元光譜接近真實(shí)的端元光譜，從而獲得全域最優(yōu)解。

通過NMF線性混合模型式(1)非負(fù)性約束得到滿足，為滿足式(1)和為1約束，采用文獻(xiàn)[12]的方法：

(11)

I-NMF的投影梯度迭代公式為

(12)

投影函數(shù)[13]選擇為PΩ(x)=max(0,x)

2 仿真實(shí)驗(yàn)

2.1 性能指標(biāo)

光譜角距離(spectral angle distance, SAD)和均方根誤差(root mean square error, RMSE)常被用來計(jì)算光譜和豐度解混估計(jì)值與真實(shí)值的近似程度，并用其來衡量解混效果。對(duì)于第i個(gè)端元，SAD定義為

(13)

RMSE定義為

(14)

2.2 模擬數(shù)據(jù)實(shí)驗(yàn)

從USGS礦物光譜庫中選取5種線性獨(dú)立的端元光譜(明礬石、水銨長石、方解石、高嶺石、白云母)，按Dirichlet分布混合，端元豐度之和進(jìn)行歸一化，加上不同的白噪聲，形成模擬實(shí)驗(yàn)數(shù)據(jù)。

不同像元個(gè)數(shù)解混性能實(shí)驗(yàn)：像元個(gè)數(shù)不同時(shí)，I-NMF與MVCNMF、PSSCNMF算法解混性能比較，如圖2。像元數(shù)目H分別為1 600、3 600、6 400、10 000，信噪比固定為20 dB。

無純像元情況的解混性能實(shí)驗(yàn)：對(duì)I-NMF與MVCNMF、PSSCNMF算法在無純像元時(shí)解混性能比較，如圖3。全部端元的豐度A從1降到0.55(1表示數(shù)據(jù)中存在純像元)，信噪比為30 dB。觀察實(shí)驗(yàn)結(jié)果，I-NMF算法解混性能要優(yōu)于MVCNMF、PSSCNMF算法。

2.3 真實(shí)數(shù)據(jù)實(shí)驗(yàn)

采用美國JPL實(shí)驗(yàn)室的AVIRIS獲得的內(nèi)華達(dá)州Cuprite高光譜圖像，圖4為第183、193、207波段圖像合成的偽彩色圖像。去除信噪比低和水蒸氣吸收波段余下188個(gè)有效波段。文獻(xiàn)[14]給出了該地區(qū)的真實(shí)地物分布。

圖4 Cuprite偽彩色圖像

(a)明礬石 (b)水銨長石 (c)方解石

(d)黃鉀鐵礬 (e)高嶺石 (f)蒙脫石

(g)白云母 (h)綠脫石 (i)榍石

表1 Cuprite數(shù)據(jù)的光譜角比較

3 結(jié)束語

針對(duì)線性混合像元分解的特點(diǎn)，提出一種基于獨(dú)立性約束的非負(fù)矩陣分解的混合像元分解方法。該方法適用于無純像元的情況。通過模擬數(shù)據(jù)與真實(shí)數(shù)據(jù)實(shí)驗(yàn)分析，表明本文提出的方法能夠提高混合像元分解的精度。投影梯度算法的收斂速度較慢，需采用更有效的優(yōu)化算法來提高解混的效率。

參考文獻(xiàn)：

[1]杜培軍，譚琨，夏俊士. 高光譜遙感影像分類與支持向量機(jī)應(yīng)用研究[M]. 北京：科學(xué)出版社，2012:6-9.

[2]LEE D D, SEUNG H S. Learning the parts of objects by non-negative matrix factorization[J]. Nature, 1999, 401:788-791.

[3]劉雪松，王斌，張立明. 基于非負(fù)矩陣分解的高光譜遙感圖像混合像元分解[J]. 紅外與毫米波雷達(dá)學(xué)報(bào)，2011，30(1):27-32.

LIU Xuesong, WANG Bin, ZHANG Liming. Hyperspectral unmixing based on nonnegative matrix factorization[J]. J Infrared Millim Waves, 2011, 30(1):27-32.

[4]MIAO Lidan, QI Hairong. Endmember extraction from highly mixed data using minimum volume constrained nonnegative matrix factorization[J]. IEEE Trans Geosci Remote Sens, 2007, 45(3):765-777.

[5]賈森，錢沄濤，紀(jì)震，等. 基于光譜和空間特性的高光譜解混方法[J]. 深圳大學(xué)學(xué)報(bào):理工版，2009, 26(3):262-267.

JIA Sen, QIAN Yuntao, JI Zhen, et al. Spectral and spatial character-based hyperspectral unmixing[J]. Journal of Shenzhen University:Science and Engineering, 2009, 26(3):262-267.

[6]趙春暉，成寶芝，楊偉超. 利用約束非負(fù)矩陣分解的高光譜解混算法[J]. 哈爾濱工程大學(xué)學(xué)報(bào)，2012, 33(3):377-382.

ZHAO Chunhui，CHENG Baozhi，YANG Weichao. Algorithm for hyperspectral unmixing using constrained nonnegative matrix factorization[J]. Journal of Harbin Engineering University, 2012, 33(3):377-382.

[7]吳波，趙銀娣，周小成. 端元約束下的高光譜混合像元非負(fù)矩陣分解[J]. 計(jì)算機(jī)工程，2008, 34(22):229-233.

WU Bo，ZHAO Yindi，ZHOU Xiaocheng. Unmixing mixture pixels of hyperspectral imagery using endmember constrained nonnegative matrix factorization[J]. Computer Engineering, 2008, 34(22):229-233.

[8]HYVARINEN A, KARHUNEN J, OJA E. Independent component analysis[M]. New York: Wiley, 2001:1-5.

[9]李娜，趙慧潔. 高光譜數(shù)據(jù)非監(jiān)督分類的改進(jìn)獨(dú)立成分分析方法[J]. 國土資源遙感，2011, 89(2):70-74.

LI Na, ZHAO Huijie. An improved independent component analysis method for unsupervised classification of hyperspectral data[J]. Remote Sensing for Land & Resources,2011, 89(2):70-74.

[10]路威，余旭初，馬永剛，等. 高光譜遙感影像海上艦船目標(biāo)檢測(cè)算法的研究[J]. 海洋測(cè)繪，2005, 25(4):8-12.

LU Wei, YU Xuchu, MA Yonggang, et al. Research of sea fleet targets detection approach in hyperspectral image[J]. Hydrographic Surveying and Charting, 2005, 25(4):8-12.

[11]ARNAVUT Z, NARUMALANI S. Lossless compression of multispectral images using permutations[C]//IEEE IGARSS ’96, 1996:463-465.

[12]HEINZ D C, CHANG C I. Fully constrained least squares linear spectral mixture analysis method for material quantification in hyperspectral imagery[J]. IEEE Trans Geosci Remote Sens, 2001, 39(3):529-545.

[13]GILLIS N, GLINEUR F. Using under-approximations for sparse nonnegative matrix factorization[J]. Pattern Recognition, 2010, 43(4):1676-1687.

[14]CIARK R N, SWAYZE G A. Evolution in imaging spectroscopy analysis and sensor signal-to-noise: an examination of how far we have come[C/OL]//The 6th Annual JPL Airborne Earth Science Work Shop, 1996. http://speclab.cr.usgs.gov/PAPERS.imspec. evol/aviris.evolution.html.