基于流形結(jié)構(gòu)化稀疏先驗(yàn)的高光譜圖像壓縮感知方法
【技術(shù)領(lǐng)域】
[0001] 本發(fā)明涉及一種高光譜圖像壓縮感知方法，特別是涉及一種基于流形結(jié)構(gòu)化稀疏 先驗(yàn)的高光譜圖像壓縮感知方法。
【背景技術(shù)】
[0002] 高光譜圖像包含成百上千的波段，每個(gè)像素包含一條連續(xù)的光譜。豐富的光譜信 息使得高光譜圖像在目標(biāo)檢測(cè)、識(shí)別等方面具有極大的優(yōu)勢(shì)，然而，高光譜圖像巨大的數(shù)據(jù) 量對(duì)圖像的采集、傳輸和處理提出了苛刻的要求，制約了其實(shí)際應(yīng)用。因此，高光譜圖像壓 縮是高光譜領(lǐng)域的熱點(diǎn)研究之一。壓縮感知成像理論證明僅需要采集少量的線性觀測(cè)值便 可以精確重建原始場(chǎng)景的圖像。相對(duì)于傳統(tǒng)的圖像壓縮算法，極大地減少了成像過(guò)程中的 資源消耗。
[0003]ChengBoLi等人在文南犬"Acompressivesensingandunmixingschemeforh yperspectraldataprocessing,IEEETransactionsonImageProcessing, 2012, 21 (3): 1200 - 1210"中公開(kāi)了一種高效的高光譜圖像壓縮感知算法。成像過(guò)程中，使用單像素相 機(jī)采集少量線性觀測(cè)值作為壓縮數(shù)據(jù)。重建過(guò)程中，基于線性混合模型，引入適量的端元光 譜，重建空間梯度稀疏的豐度值矩陣。最后，通過(guò)線性混合重建的豐度值矩陣和引入的端元 光譜重建高光譜圖像。然而，該算法僅考慮了空間的稀疏性，未能充分利用高光譜圖像的三 維結(jié)構(gòu)，重建精度受限；其次，該算法的稀疏性約束噪聲魯棒性差；此外，算法性能嚴(yán)重依 賴端元光譜的選擇，實(shí)用性受限。

【發(fā)明內(nèi)容】

[0004] 為了克服現(xiàn)有高光譜圖像壓縮感知方法精度低的不足，本發(fā)明提供一種基于流形 結(jié)構(gòu)化稀疏先驗(yàn)的高光譜圖像壓縮感知方法。該方法隨機(jī)采樣每個(gè)像素光譜的少量線性觀 測(cè)值作為壓縮數(shù)據(jù)，通過(guò)流形結(jié)構(gòu)化稀疏先驗(yàn)，同時(shí)刻畫高光譜圖像稀疏化后光譜維中的 稀疏性和空間維中的流形結(jié)構(gòu)；通過(guò)隱變量貝葉斯模型，將信號(hào)重建，稀疏先驗(yàn)學(xué)習(xí)以及噪 聲估計(jì)統(tǒng)一到一個(gè)正則化回歸模型進(jìn)行優(yōu)化求解。學(xué)習(xí)得到的稀疏先驗(yàn)既能充分地刻畫高 光譜圖像的三維結(jié)構(gòu)，又具有較強(qiáng)的噪聲魯棒性。利用該稀疏先驗(yàn)，實(shí)現(xiàn)了高光譜圖像的高 精度重建。在真實(shí)的高光譜遙感數(shù)據(jù)集Urban上的實(shí)驗(yàn)結(jié)果表明，當(dāng)在壓縮數(shù)據(jù)中加入高 斯白噪聲使得壓縮數(shù)據(jù)信噪比為15db，采樣率為0. 09時(shí)，獲得了 23db的峰值信噪比。
[0005] 本發(fā)明解決其技術(shù)問(wèn)題所采用的技術(shù)方案是：一種基于流形結(jié)構(gòu)化稀疏先驗(yàn)的高 光譜圖像壓縮感知方法，其特點(diǎn)是包括以下步驟：
[0006] 步驟一、針對(duì)包含nb個(gè)波段，每個(gè)波段包含η^行和η。列的高光譜圖像，將每一個(gè) 波段拉伸成為一個(gè)行向量，重新組成一個(gè)二維矩陣，尤(/%=?,.x?f)。其中，X的每 一列對(duì)應(yīng)每個(gè)像素的光譜；每一行對(duì)應(yīng)每個(gè)波段的所有像素值。將行和列分別稱為空間維 和光譜維。
[0007] 獲取壓縮數(shù)據(jù)過(guò)程中，利用列歸一化的高斯隨機(jī)觀測(cè)矩陣瓶￥~隨機(jī)采樣高光 譜圖像X的光譜維，獲得壓縮數(shù)據(jù)F總],mb為壓縮后波段長(zhǎng)度。
[0008] F=AX+N(1)
[0009] 其中，jVeKWiX~表示采樣中的噪聲。采樣率P定義為P=mb/nb。
[0010] 步驟二、利用Haar小波基對(duì)高光譜圖像的每個(gè)光譜進(jìn)行稀疏化，如X=ΨΥ，ψ為 小波基，Y為列稀疏的系數(shù)矩陣，模型（1)表示為F=AWY+N。假設(shè)采樣過(guò)程中噪聲N服從 ...ν?Λ^Ο,Σ,,,/)的矩陣正太分布，I為對(duì)應(yīng)大小的單位矩陣，模型（1)對(duì)應(yīng)的似然函數(shù)定義為
[0011]
⑵：
[0012]其中：
，Ση=diag(x)表示以λ為對(duì)角線元素的對(duì)角矩陣。
I示Q矩陣的加權(quán)跡范數(shù)。
[0013] 除過(guò)列稀疏性，高光譜圖像空間像素之間的相似性使得Υ中不同列的稀疏信號(hào)位 于一個(gè)結(jié)構(gòu)未知的流形結(jié)構(gòu)上。為了充分描述Υ的特性，假設(shè)Υ服從如下的矩陣正太分布
[0014]
(3)
[0015] 為描述Υ中列信號(hào)的稀疏性，令Σν=diag(Y)表示以γ為對(duì)角線元素的對(duì)角 矩陣，
[設(shè)γ服從如下的伽瑪分布
[0016]
(4)
[0017] 式子⑶中，Σα描述Υ中不同列信號(hào)之間的相關(guān)性，因此，式子⑶隱式地表示 Υ中不同的稀疏信號(hào)之間存在的流形結(jié)構(gòu)。為了更加靈活地學(xué)習(xí)Σα，進(jìn)一步假設(shè)Σα服從 如下的反威沙特分布
[0018]
CS)
[0019] 其中，1是給定的常量，表示自由度，Γηρ是多變量伽瑪函數(shù)，為參考協(xié)方 差矩陣。該先驗(yàn)通過(guò)最小化Σ。#Θ之間的布雷格曼散度，使得Σq趨近于Θ，從而減輕 了 學(xué)習(xí)過(guò)程中的過(guò)擬合問(wèn)題。
[0020] 步驟三、為使得流行結(jié)構(gòu)化稀疏先驗(yàn)?zāi)軌蚋玫仄ヅ鋱D像分布并具有較強(qiáng)的噪聲 魯棒性，通過(guò)隱變量貝葉斯模型對(duì)噪聲參數(shù)λ和先驗(yàn)參數(shù)γ，κ，2。￥和Θ進(jìn)行估計(jì)。令 f=vec(F)，y=vec(Y)，n=vec(N)和Φ= /?(/?Ψ)，vec(Q)表示將矩陣Q拉成列向量， 發(fā)表示克羅內(nèi)克積，則模型（2)等價(jià)于
[0021]
籩)
[0022] 同樣，模型⑶中關(guān)于Y的先驗(yàn)等價(jià)于
[0023]
(7)
[0024]根據(jù)公式（6)、公式（7)，所有的未知參數(shù)通過(guò)求解如下的優(yōu)化問(wèn)題得到
[0025]
[0026]通過(guò)積分，并引入_21og運(yùn)算，容易得知式子（8)，等價(jià)于最小化如下的式子
[0027]
[0028] 其中，tr(·)表示跡范對(duì)式子（9)的第一項(xiàng)做如下變形
[0029]?〇) J
[0030] 將式子（10)帶入到式子（9)中，得到如下等價(jià)于式子⑶的正則化回歸模型 12345678910 、L··!/· --· Νρ --· WI
(.--k'WII 2 該模型將信號(hào)重建、稀疏先驗(yàn)學(xué)習(xí)和噪聲估計(jì)統(tǒng)一到一個(gè)框架下。 3 步驟四、為提升算法效率，引入如下的近似關(guān)系， 4
[0034]
(12) 5 基于關(guān)系（12)，采用坐標(biāo)下降法將式子（11)分解為若干個(gè)子問(wèn)題進(jìn)行迭代求解， 每個(gè)子問(wèn)題中僅優(yōu)化一個(gè)變量而固定剩余的其他變量。具體步驟如下： 6 ①初始化λ，γ，κ為對(duì)應(yīng)長(zhǎng)度的全1向量，I，計(jì)數(shù)變量t= 0; 7 ②學(xué)習(xí)參考協(xié)方差矩陣Θ。定義關(guān)于觀測(cè)值矩陣F的權(quán)值矩陣Μ 8
(13) 9 My為Μ的i行j列的元素，Λ/7eir%xi2表示空間中以第i個(gè)像素為中心，大小為 k= 3的鄰域窗口中的所有光譜的觀測(cè)值。II·I|F表示弗羅貝尼烏斯范數(shù)，。=0. 7。參 考協(xié)方差矩陣Θ= (D-M) 4為對(duì)角陣， 10 ③固定λ和γ，根據(jù)式子（11)得到關(guān)于γ的子問(wèn)題，如下
[0041]
(14)
[0042]基于近似關(guān)系（12)，求解得到Y(jié)的更新規(guī)則如下，
[0043] Υ=ΣΓγΨτΑτ(Ση+ΑΨΣΓγΨτΑτ) (15)
[0044] ④固定Y，λ，κ和Σ"，利用近似關(guān)系（12)得到關(guān)于γ的子問(wèn)題，如下
[0045]
(1.6)
[0046] 其中，Yi.表示Υ的第i行，γΛγ的第i個(gè)元素，求解得到如下的更新形式：
[0047]
(17)
[0048] 其中，a=diag[2ry-2ryWTAT(2n+AW2ryWTAT)ΙΨΣ」，與之前不同，此處 diag( ·)表示取矩陣對(duì)角線元素組成向量，aα的第i個(gè)元素。
[0049] ⑤固定Y和Y，利用近似關(guān)系（12)得到5^的子問(wèn)題
[0050]
(18)
[0051] μ=mb+np+l+l，求解得到xcy的更新形式，如下：
[0052]
(19)
[0053] 為提升噪聲魯棒性，令# = 1丨儼+ + 。
[0054] ⑥固定Υ和丫，利用近似關(guān)系（12)得到關(guān)于λ的優(yōu)化子問(wèn)題，如下
[0055]
(20)
[0056] 求解得到如下的更新形式：
[0057]
(21)
[0058] 其中，λ;為λ的第i個(gè)元素，Q=ΑΨΥ-F，Q」表示Q的第i列，υ;為向量υ= diag[(Xn+AW2ryWTAT) 1的第i個(gè)元素，diag( ·)運(yùn)算和④步相同。
[0059] ⑦固定丫，得到關(guān)于κ的優(yōu)化子問(wèn)題，如下
[0060](22)
[0061] !^為κ的第
i個(gè)元素。求解得到如下的更新形式
[0062] κ1= 2/γ, (23)
[0063] ⑧假設(shè)上一次迭代重建得到的稀疏信號(hào)為Υ'，最新重建的稀疏信號(hào)為Υ，說(shuō)計(jì)算 更新前后的差異，n= |^ -y||f/|^ ||F，計(jì)數(shù)器t加1。如果計(jì)數(shù)器t彡2〇〇并且更 新差異η彡?ο4,則循環(huán)執(zhí)行步驟③至⑧；否則，退出循環(huán)。
[0064] 假設(shè)最終得到最優(yōu)估計(jì)的，則重建高光譜圖像，χ_=ΨΥ
[0065] 本發(fā)明的有益效果是：該方法隨機(jī)采樣每個(gè)像素光譜的少量線性觀測(cè)值作為壓縮 數(shù)據(jù)，通過(guò)流形結(jié)構(gòu)化稀疏先驗(yàn)，同時(shí)刻畫高光譜圖像稀疏化后光譜維中的稀疏性和空間 維中的流形結(jié)構(gòu)；通過(guò)隱變量貝葉斯模型，將信號(hào)重建，稀疏先驗(yàn)學(xué)習(xí)以及噪聲估計(jì)統(tǒng)一到 一個(gè)正則化回歸模型進(jìn)行優(yōu)化求解。學(xué)習(xí)得到的稀疏先驗(yàn)既能充分地刻畫高光譜圖像的三 維結(jié)構(gòu)，又具有較強(qiáng)的噪聲魯棒性。利用該稀疏先驗(yàn)，實(shí)現(xiàn)了高光譜圖像的高精度重建。在 真實(shí)的高光譜遙感數(shù)據(jù)集Urban上的實(shí)驗(yàn)結(jié)果表明，當(dāng)在壓縮數(shù)據(jù)中加入高斯白噪聲使得 壓縮數(shù)據(jù)信噪比為15db，采樣率為0. 09時(shí)，獲得了 23db的峰值信噪比。
[0066] 下面結(jié)合【具體實(shí)施方式】對(duì)本發(fā)明作詳細(xì)說(shuō)明。
【具體實(shí)施方式】
[0067] 本發(fā)明基于流形結(jié)構(gòu)化稀疏先驗(yàn)的高光譜圖像壓縮感知方法具體包括以下步 驟：
[0068] 針對(duì)包含nb個(gè)波段，每個(gè)波段包含η^行和η。列的高光譜圖像，將每一個(gè)波段拉伸 成為一個(gè)行向量，重新組成一個(gè)二維矩陣，工e \X' )。其中，X的每一列對(duì)應(yīng) 每個(gè)像素的光譜；每一行對(duì)應(yīng)每個(gè)波段的所有像素值。將行和列分別稱為空間維和光譜維。 本發(fā)明主要包含以下四個(gè)步驟：
[0069] 1、獲取壓縮數(shù)據(jù)。
[0070] 壓縮過(guò)程中，利用列歸一化的高斯隨機(jī)觀測(cè)矩陣J 隨機(jī)采樣高光譜圖像X 的光譜維，獲得壓縮數(shù)據(jù)Fe ~，mb為壓縮后波段長(zhǎng)度。
[0071]F=AX+N(1)其中，iVelT，表示采樣中的噪聲。采樣率P定義為P=mb/ nb〇
[0072] 2、建立基于流形結(jié)構(gòu)化稀疏先驗(yàn)的壓縮感知模型。
[0073] 利用Haar小波基對(duì)高光譜圖像的每個(gè)光譜進(jìn)行稀疏化，如X=ΨΥ，Ψ為小波基， Y為列稀疏的系數(shù)矩陣。因此，模型（1)可表示為F=AWY+N。假設(shè)采樣過(guò)程中噪聲N服
從Λ心\qo,Σs，i)的矩陣正太分布，I為對(duì)應(yīng)大小的單位矩陣（下同）。因此，模型⑴對(duì)應(yīng) 的似然涵撒πτw宙々為 12345 ⑵ V"}j 2 其中：
，Ση=diag(A)表示以λ為對(duì)角線元素的對(duì)角矩陣。
I;示Q矩陣的