蔡榮榮，王 斌

(1. 復(fù)旦大學(xué) 電磁波信息科學(xué)教育部重點(diǎn)實(shí)驗(yàn)室，上海 200433；2. 復(fù)旦大學(xué) 信息科學(xué)與工程學(xué)院 智慧網(wǎng)絡(luò)與系統(tǒng)研究中心，上海 200433)

一種基于低秩表示的高光譜圖像恢復(fù)方法

蔡榮榮1,2，王 斌1,2

(1. 復(fù)旦大學(xué) 電磁波信息科學(xué)教育部重點(diǎn)實(shí)驗(yàn)室，上海 200433；2. 復(fù)旦大學(xué) 信息科學(xué)與工程學(xué)院 智慧網(wǎng)絡(luò)與系統(tǒng)研究中心，上海 200433)

高光譜圖像(HSI)在獲取過程中不可避免地受到各種噪聲的干擾，如高斯白噪聲、沖擊噪聲、壞死的線條等．為了確保后續(xù)應(yīng)用能夠順利進(jìn)行，高光譜圖像的恢復(fù)是一項(xiàng)重要的預(yù)處理過程．文中提出一種新的高光譜圖像恢復(fù)方法，主要有以下兩點(diǎn)貢獻(xiàn)： 其一，將低秩表示模型引入到高光譜圖像恢復(fù)中；其二，利用高光譜圖像的相似圖像塊進(jìn)行聯(lián)合低秩表示．模擬和實(shí)際高光譜圖像數(shù)據(jù)的實(shí)驗(yàn)結(jié)果表明，所提出的方法能有效地去除各種噪聲干擾，同時較好地恢復(fù)圖像細(xì)節(jié)，因此可作為一種實(shí)現(xiàn)高光譜圖像恢復(fù)的有效手段．

高光譜圖像；低秩表示；相似圖像塊；聯(lián)合；恢復(fù)

近年來高光譜圖像(Hyperspectral Imagery, HSI)由于其豐富的光譜信息，有著廣泛的應(yīng)用，特別是在精準(zhǔn)農(nóng)業(yè)、植被生態(tài)、地質(zhì)礦產(chǎn)、大氣環(huán)境監(jiān)測等領(lǐng)域中有較大的應(yīng)用價值．但是高光譜圖像在獲取過程中不可避免地受到各種干擾，例如高斯白噪聲、稀疏的沖擊噪聲、壞死的線條等．這不僅使圖像的目視效果變差，更對高光譜圖像的上述應(yīng)用造成很大的影響．因此高光譜圖像的恢復(fù)是一項(xiàng)很重要的預(yù)處理過程，其圖像恢復(fù)的效果將直接影響高光譜圖像后續(xù)應(yīng)用的精度．

高光譜圖像的恢復(fù)研究一直是高光譜領(lǐng)域比較活躍的研究課題．迄今為止，國內(nèi)外學(xué)者已經(jīng)提出了多種不同的高光譜圖像恢復(fù)方法．傳統(tǒng)的高光譜圖像的恢復(fù)方法將每一個波段視為獨(dú)立的二維圖像，對各波段逐一應(yīng)用比較成熟的二維圖像的去噪方法，比如基于小波的方法[1]、基于全變分(Total-Variation, TV)的方法[2]、三維塊匹配(Block-Matching 3-D, BM3D)算法[3]等．這樣的逐波段恢復(fù)方法只利用了各個波段自身的空間信息，忽略了不同波段之間的相關(guān)性，對于某些受到強(qiáng)烈噪聲影響的波段恢復(fù)效果較差．

不同于逐波段對高光譜圖像進(jìn)行處理，有學(xué)者將高光譜圖像看作一個三維圖像，提出了多維度的高光譜圖像恢復(fù)方法．例如，Letexier和Bourennane[4]將高光譜圖像看作是一個三維張量，使用Tucker張量分解模型，引入了廣義的多維維納濾波對其進(jìn)行恢復(fù)；而文獻(xiàn)[5]使用PARAFAC(Parallel Factor Analysis)張量分解模型對高光譜圖像進(jìn)行恢復(fù)．更進(jìn)一步，有學(xué)者將BM3D方法拓展到三維圖像上，提出了視頻三維塊匹配(Video Block-Matching 3-D, VBM3D)算法[6]．VBM3D對于各個波段統(tǒng)計特性相同的噪聲恢復(fù)效果較好，但是在高光譜圖像中，各種類型噪聲的存在導(dǎo)致各個波段的噪聲強(qiáng)度不一，VBM3D對于這種情況下的恢復(fù)效果就顯得不夠了．

近幾年來，國內(nèi)外學(xué)者將低秩矩陣的概念引入到高光譜圖像恢復(fù)領(lǐng)域中[7-9]，取得較好的效果．文獻(xiàn)[7]引入低秩矩陣的概念去除高光譜圖像的壞死線噪聲，考慮到局部的幾何結(jié)構(gòu)并添加了圖的正則化．文獻(xiàn)[8]基于低秩矩陣恢復(fù)(Low-Rank Matrix Recovery, LRMR)較為成功地去除了高光譜圖像的混合噪聲，其正是得益于無噪的高光譜圖像數(shù)據(jù)低秩的特性．文獻(xiàn)[9]在LRMR的框架下添加了TV的正則化，進(jìn)一步提升了圖像恢復(fù)效果．但是，LRMR框架是基于魯棒主成分分析(Robust Principal Component Analysis, RPCA)模型．在RPCA的基礎(chǔ)上，Liu等將RPCA發(fā)展到低秩表示(Low-Rank Representation, LRR)模型[10]，認(rèn)為低秩的數(shù)據(jù)可以源自于多個子空間．相較RPCA而言，LRR模型更為符合高光譜圖像的線性混合模型．這兩年來LRR模型已經(jīng)被引入到高光譜圖像處理領(lǐng)域，用作估計子空間的數(shù)目[11]、異常探測[12]，由于其模型更加符合高光譜圖像的特點(diǎn)，所以這些方法都較為成功．另外，文獻(xiàn)[8-9]的恢復(fù)方法只利用了高光譜圖像的每一小塊圖像的局部信息．如果高光譜圖像在某些區(qū)域受到強(qiáng)烈噪聲影響導(dǎo)致某些塊的噪聲強(qiáng)度很大以及局部信息受到巨大影響，不引入新的信息將很難得到較好的恢復(fù)結(jié)果．

針對典型低秩分解方法只利用單一圖像塊的問題，本文提出一種基于LRR的高光譜圖像恢復(fù)方法．這里，我們將LRR模型引入到高光譜圖像恢復(fù)領(lǐng)域；更進(jìn)一步，利用高光譜圖像的非局部相似塊的信息，以彌補(bǔ)單個塊局部信息不足、樣本容量不夠的缺陷．在所提出的方法中，首先，通過最小二乘估計得到初步無噪的高光譜圖像，由于其能夠比較好表示數(shù)據(jù)子空間，從中選取相應(yīng)圖像塊的估計結(jié)果作為LRR字典可以更好地凸顯LRR模型的優(yōu)勢；再根據(jù)各個波段的噪聲水平加權(quán)融合出單幅圖像，對單幅圖像進(jìn)行圖像塊分類，將分類結(jié)果對應(yīng)到三維高光譜圖像中，由于融合圖像基本不含噪聲，分類結(jié)果較為準(zhǔn)確；最后對于每一高光譜圖像塊，在同一類中找出多個與之相似的塊，進(jìn)行聯(lián)合LRR，因?yàn)橄嗨茍D像塊的具有相似的低秩結(jié)構(gòu)，聯(lián)合LRR可以加大低秩樣本容量，利于恢復(fù)出圖像塊的低秩結(jié)構(gòu)．低秩恢復(fù)得到的重構(gòu)矩陣不是單一的圖像塊，而是多個相似的圖像塊．一方面使用了更符合高光譜圖像數(shù)據(jù)的LRR模型，另一方面使用非局部相似塊的信息．仿真和實(shí)際高光譜圖像的實(shí)驗(yàn)結(jié)果都表明，所提出方法的性能表現(xiàn)優(yōu)于現(xiàn)有的典型高光譜圖像恢復(fù)方法，驗(yàn)證了所提出方法的有效性．另外，對于實(shí)驗(yàn)參數(shù)的選擇，本文也做了詳細(xì)的分析和討論．

1 高光譜圖像的低秩模型

1.1無噪高光譜圖像的低秩特性

基于光譜混合模型，無噪的高光譜圖像通常遵循線性混合模型[13]，表示如下：

X=AS，

(1)

其中：X∈n×b為無噪高光譜圖像數(shù)據(jù)，n為像素個數(shù)(如果高光譜圖像的空間維數(shù)為h和w，那么n=h×w)，b為波段個數(shù)；S∈P×b為端元矩陣，表示為S=[s1,s2,…,sP]T，si∈b,P為端元個數(shù)；A∈n×P為端元矩陣，它的元素ai,j表示在線性混合模型中第j個端元在第i個像素中所占的比例系數(shù)．在大部分情況下，高光譜圖像中端元的個數(shù)遠(yuǎn)小于波段數(shù)或者像素個數(shù)[13]，即P·min(n,b)．又由于X是A和S的乘積，所以X的秩滿足：

(2)

顯然min(rank(A),rank(S))≤P，所以得到如下結(jié)論：

(3)

這意味著矩陣X是低秩的，或者說矩陣X的列所張成的光譜空間的維數(shù)小于或等于端元個數(shù)P．這就證明了無噪高光譜圖像的低秩特性．

1.2低秩表示模型

當(dāng)高光譜圖像受到高斯白噪聲的微小擾動時，含噪的高光譜圖像數(shù)據(jù)Y∈n×b可以分解成一個低秩的矩陣X和一個微小擾動的矩陣N如下：

Y=X+N.

(4)

主成分分析(Principle Component Analysis, PCA)通過奇異值分解(Singular Value Decomposition, SVD)來完成數(shù)據(jù)降維，可以有效地從Y中恢復(fù)出X．然而當(dāng)數(shù)據(jù)中存在比較大的噪聲或者異常時，N就用E代替，E的元素為任意大小，此時PCA的效果就不夠準(zhǔn)確．針對這種情況，研究學(xué)者提出了RPCA[14]，把問題歸結(jié)為以下凸優(yōu)化問題：

(5)

其中：λ>0為正則化參數(shù)；‖·‖*表示矩陣的核范數(shù)，為矩陣的奇異值之和；‖·‖1表示矩陣的l1范數(shù)，為矩陣元素絕對值之和．

min‖Z‖*+λ‖E‖2，1
s.t.Y=DZ+E，

(6)

其中：λ>0為正則化參數(shù)；D為能夠生成數(shù)據(jù)子空間的字典；‖·‖2,1是l2,1范數(shù)，表示矩陣每一列的l2范數(shù)之和．

解決問題(6)中的優(yōu)化問題可以得出最優(yōu)的Z*和E*．Z*可以看成數(shù)據(jù)Y對于字典D的最低秩表示．由于rank(DZ*)≤rank(Z*)，DZ*可以看成是原始數(shù)據(jù)Y低秩恢復(fù)的結(jié)果．噪聲由于不能被數(shù)據(jù)子空間表示，被分離出來，用矩陣E*表示．令D=I，優(yōu)化問題(6)就退化成(5)．所以LRR可以看成是RPCA模型的更一般的形式，換言之RPCA可以看成是使用了標(biāo)準(zhǔn)基作為字典的LRR模型．

使用增廣的拉格朗日乘子法(Augmented Lagrange Multipliers, ALM)求解上述優(yōu)化問題，先將其轉(zhuǎn)化為以下等價問題：

(7)

通過最小化以下增廣拉格朗日函數(shù)L可得出最優(yōu)的Z*和E*：

(8)

其中：Y1和Y2為拉格朗日乘子；μ>0為懲罰參數(shù)．更詳細(xì)的求解過程和算法收斂性證明可參考文獻(xiàn)[10]．

2 基于低秩表示的恢復(fù)方法

2.1空間自適應(yīng)的相似塊的搜尋

圖1 含噪高光譜圖像中的相似塊Fig.1 Similar blocks in the noisy HSI

如前所述，典型的低秩恢復(fù)方法都是對單一圖像塊分別處理的，忽略了高光譜圖像中存在非局部的相似塊．這些相似圖像塊對應(yīng)的地物類型相似且具有很類似的結(jié)構(gòu)特征，導(dǎo)致這些圖像塊具有相似的空間光譜信息．因此，這些三維圖像塊在轉(zhuǎn)換為二維圖像矩陣后具有近似的低秩特征，也可以理解為由這些二維圖像矩陣組成的大圖像矩陣的秩與單個圖像矩陣的秩是很近似的．如果將這些二維圖像放在一起進(jìn)行聯(lián)合LRR，無疑可以增加圖像矩陣的樣本容量，這對于優(yōu)化問題(6)中從Y中恢復(fù)出低秩的部分是有利的．另一方面，由于相似圖像塊多次參與LRR表示，因此這些圖像塊有多個低秩恢復(fù)的結(jié)果，對這些結(jié)果加權(quán)平均得到最終恢復(fù)結(jié)果．例如，圖1為高光譜圖像中一些相似塊，在恢復(fù)紅框里的高光譜圖像塊時，引入綠框和藍(lán)框里的高光譜圖像塊來進(jìn)行聯(lián)合LRR，這對于恢復(fù)紅框里無噪圖像塊的幫助是顯而易見的．另一方面，在恢復(fù)綠框和藍(lán)框里的圖像塊時，紅框里的圖像塊被引入到聯(lián)合LRR，這樣紅框里的圖像塊就得到了3次LRR的結(jié)果．

本節(jié)提出了一種空間自適應(yīng)的相似塊的搜尋方法，分為以下兩步： 首先需要將高光譜圖像根據(jù)各個波段的噪聲情況融合成單幅圖像．然后將融合后的圖像分成相互重疊的圖像塊，提取每一圖像塊的局部自適應(yīng)回歸核(Locally Adaptive Regression Kernels, LARKs)[15-16]特征，再將這些圖像塊聚成K類．

2.2.1 估計各個波段的噪聲，將高光譜圖像融合成單幅圖像

噪聲估計在高光譜圖像中一直是一個重要問題，也有許多方法可以解決，文中采用基于多重回歸理論的方法[17]來進(jìn)行噪聲估計．

用Y=[y1,y2,…，yb]表示含噪高光譜圖像，其中Y的每一列yi是高光譜圖像第i個波段拉伸成的矢量．用Y?i=[y1,y2,…，yi-1,yi+1,yi+2,…,yb]表示不含第i個波段的含噪高光譜圖像．假設(shè)yi可以用剩余的b-1個波段的線性組合表示如下：

yi=Y?iβi+ξi，

(9)

其中：Y?i∈n×(b-1)；βi∈(b-1)×1為回歸向量；ξi∈n×1為噪聲向量．對于每一個i∈{1,2,…，b}，回歸向量的最小二乘估計由下式給出：

(10)

然后可得到第i個波段的信號和噪聲的估計如下：

(11)

(12)

最后各個波段的噪聲方差可由下式得出：

(13)

為了盡可能保留高光譜圖像每一個波段的信息，采用加權(quán)平均的思想得到融合圖像．加權(quán)的策略為噪聲強(qiáng)度大的波段權(quán)值較小，噪聲強(qiáng)度小的波段權(quán)值較大，具體如下：

(14)

圖2 高光譜圖像的融合Fig.2 Fusing of HSI

其中權(quán)值ω(yi)是通過高斯核函數(shù)定義：

(15)

其中：h=max(σi)．上式定義的權(quán)值使得那些受噪聲影響小的高質(zhì)量的波段對于融合圖像有顯著的貢獻(xiàn)，使得融合圖像的質(zhì)量得到提升，為之后的圖像塊分類做好準(zhǔn)備．圖2為高光譜圖像融合成單幅圖像的例子，從圖中可以看出融合以后的圖像噪聲強(qiáng)度大大減小，對之后的相似圖像塊的搜尋大有幫助，也驗(yàn)證了提出的融合策略是正確的．

本小節(jié)的目的有兩個： 其一，產(chǎn)生高光譜圖像的初步無噪估計，為之后LRR模型字典D的選取提供基礎(chǔ)；其二，把高光譜圖像通過加權(quán)融合為單張圖像，為了方便后續(xù)相似圖像快的搜尋．

2.1.2 相似塊的搜尋

圖3 相似的圖像塊的聚類過程Fig.3 Clustering of the similar blocks

在2.1.1中已經(jīng)得到了一幅融合后的高光譜圖像，利用滑動窗將融合后的圖像分成相互重疊的圖像塊，每個像素點(diǎn)對應(yīng)以其為中心的一個正方形圖像塊．這樣，一幅圖像分成的圖像塊個數(shù)等于這幅圖像的像素點(diǎn)個數(shù)．對于每一個圖像塊，均提取它們的LARK特征，使用K-means[18]方法將這些圖像塊分成K類，每一類都包含幾何特征很相似的圖像塊．圖3所示為K=9時的融合高光譜圖像塊的聚類結(jié)果，著色方案為： 如果以像素i為中心的圖像塊Yi屬于第k類(k=1,2,…,K)，那么像素i就著以第k類所對應(yīng)的顏色．從圖中可以看出即便圖像存在噪聲，使用LARK特征對圖像塊進(jìn)行聚類得到的效果也是較準(zhǔn)確的．使用LARK特征的聚類結(jié)果比較依賴參數(shù)K，一方面K比較小時會導(dǎo)致結(jié)構(gòu)不相似的圖像塊被聚在同一類中，對后續(xù)處理產(chǎn)生影響；另一方面，K比較大的時候會導(dǎo)致每一類的圖像塊數(shù)量太少，使得后續(xù)處理的圖像塊樣本量不夠．參數(shù)K對于實(shí)驗(yàn)結(jié)果的影響，將在實(shí)驗(yàn)中進(jìn)行進(jìn)一步分析．

對融合后的圖像完成相似圖像塊的聚類之后，將此結(jié)果對應(yīng)回三維高光譜圖像中，就可以得到三維高光譜圖像塊的聚類結(jié)果了．對于每一高光譜圖像塊，之后相似塊的選取均在和它同類的塊中進(jìn)行尋找．

2.2基于聯(lián)合低秩表示的恢復(fù)方法

提出的基于聯(lián)合LRR的高光譜圖像恢復(fù)算法可以細(xì)分為以下3步： 圖像塊的組合、聯(lián)合矩陣LRR的求解和高光譜圖像的恢復(fù)．

(16)

(17)

最后，通過低秩恢復(fù)得到的圖像矩陣，重構(gòu)出無噪高光譜圖像．

(18)

(19)

其中：τ是一個調(diào)節(jié)圖像塊相似性對于權(quán)值影響的參數(shù)；di,ij表示圖像塊Yi和Yij的塊相似度．

3 實(shí)驗(yàn)驗(yàn)證

在本節(jié)中，分別使用模擬圖像數(shù)據(jù)和實(shí)際圖像數(shù)據(jù)來證明提出方法的圖像恢復(fù)效果，并對比了VBM3D[6]和LRMR[8]方法．進(jìn)一步地，為了證明LRR模型的有效性，對每一高光譜圖像塊，使用LRR模型恢復(fù)出其低秩結(jié)構(gòu)，將此方法稱為LRR方法；另外，為了證明引入相似圖像塊的作用，對于每一高光譜圖像塊，找出和它最相似的一些塊，進(jìn)行聯(lián)合LRMR，將此方法稱為聯(lián)合LRMR方法，即聯(lián)合LRMR和聯(lián)合LRR方法的唯一區(qū)別在于進(jìn)行低秩分解時一個采用LRMR模型，另一個采用LRR模型．我們將LRR和聯(lián)合LRMR方法也作為對比方法．在實(shí)驗(yàn)開始之前，將高光譜圖像數(shù)據(jù)的每一個波段的像素值都?xì)w一化為[0,1]．VBM3D方法中噪聲的方差由2.1.1中估計得到；LRMR方法的參數(shù)采用文獻(xiàn)[8]的選擇；聯(lián)合LRMR方法中圖像塊大小與LRMR一致．

3.1模擬圖像實(shí)驗(yàn)及評價指標(biāo)

模擬圖像實(shí)驗(yàn)選用了2個高光譜圖像數(shù)據(jù)集： 帕維亞大學(xué)(Pavia University)數(shù)據(jù)集和華盛頓特區(qū)中心(Washington D.C. Mall)數(shù)據(jù)集．Pavia University數(shù)據(jù)集的空間大小為610×340，包含103個波段．實(shí)驗(yàn)1中選取了大小為100×100×103的子圖像為實(shí)驗(yàn)對象．Washington D.C. Mall數(shù)據(jù)集的空間大小為1208×307，包含191個波段．實(shí)驗(yàn)2中選取了大小為100×100×103的子圖像為實(shí)驗(yàn)對象．

在模擬實(shí)驗(yàn)中，對2個數(shù)據(jù)集均添加以下3種噪聲： 每個波段都添加σ=0.1的高斯白噪聲，模擬實(shí)際中的熱噪聲；在波段20，40，41，60，61，80，81隨機(jī)選擇5%的像素點(diǎn)加入σ=0.7稀疏高斯白噪聲，模擬實(shí)際中的像元壞死；在波段21，31，41，51，61，71，81波段隨機(jī)添加4條壞死的線條，模擬實(shí)際中的傳感器故障．

本文采用2個評價指標(biāo)，峰值信噪比(Peak Signal-to-Noise Ratio, PSNR)和結(jié)構(gòu)相似性(Structural Similarity, SSIM)指標(biāo)[19]定量地對恢復(fù)結(jié)果進(jìn)行評價．對無噪圖像和恢復(fù)圖像逐波段計算PSNR和SSIM值，然后取平均得到平均PSNR(Mean PSNR,mPSNR)和平均SSIM(Mean SSIM,mSSIM)值．各個指標(biāo)的定義如下：

(20)

(21)

(22)

(23)

3.1.1 實(shí)驗(yàn)1

實(shí)驗(yàn)1采用Pavia University數(shù)據(jù)集，提出方法的實(shí)驗(yàn)參數(shù)選取如下： LRR的參數(shù)λ=0.35、滑動窗即圖像塊的大小Wsize=11、圖像塊聚類的類數(shù)K=10、相似圖像塊的數(shù)量p=5、式(19)的參數(shù)τ=10．Pavia University數(shù)據(jù)集的第61波段無噪圖像如圖4(a)所示，圖4(b)是添加了3種噪聲的含噪圖像，圖4(c)～(g)分別是使用VBM3D、LRMR、LRR、聯(lián)合LRMR和聯(lián)合LRR方法恢復(fù)后的圖像．從圖中可以看出，VBM3D方法對于強(qiáng)烈的稀疏噪聲和壞死的線條去除效果較差，因?yàn)樗褂玫哪Ｐ椭荒苋コ咚拱自肼暎愀獾氖?，由于VBM3D使用聯(lián)合維納濾波對圖像進(jìn)行處理，導(dǎo)致它對圖像有過平滑的趨勢，例如圖4(c)右上角的模糊說明了這一點(diǎn)．LRMR、LRR、聯(lián)合LRMR和聯(lián)合LRR方法由于使用低秩模型，不僅可以去除高斯白噪聲，也可以去除強(qiáng)烈的稀疏噪聲和壞死的線條，恢復(fù)出比較好的圖像．不過，從圖中可以看出，LRMR方法不能完全去除所有噪聲，LRR方法相較LRMR方法的恢復(fù)效果要好一些，說明LRR模型比LRMR模型更適合于高光譜圖像的恢復(fù)．但是圖4(e)紅色圈出的區(qū)域受到噪聲影響，LRR方法的恢復(fù)效果并沒有很好，存在一些涂抹痕跡．提出的方法將相似圖像塊組合到一起進(jìn)行聯(lián)合LRR，引入了相似圖像塊的信息增加了低秩矩陣的樣本容量，更好的幫助恢復(fù)出不含噪聲的高光譜圖像．從圖4(g)看出所提出方法不僅可以去除所有類型的噪聲，同時還比較好地恢復(fù)出了圖像的空間結(jié)構(gòu)，如紅色圈出的區(qū)域的恢復(fù)效果就好于LRR方法．另外，聯(lián)合LRMR相比LRMR方法恢復(fù)效果的提升也能說明引入相似圖像塊的作用．

圖4 Pavia University無噪圖、含噪圖和不同方法恢復(fù)結(jié)果的第61波段Fig.4 Band 61 of Pavia University clean data, noisy data and restoration results by different methods

為了進(jìn)一步定量地評價提出方法的恢復(fù)效果，逐波段計算了不同方法恢復(fù)后的高光譜圖像的PSNR和SSIM值，結(jié)果如圖5所示．在波段20，21，31，40，41，51，60，61，71，80，81上，所有方法的PSNR和SSIM值均有不同程度的凹陷，原因?yàn)樵谶@些波段上存在多種類型噪聲，導(dǎo)致圖像質(zhì)量急劇下降．從圖中可以看出，提出方法的各波段的PSNR和SSIM值均好于其他方法，特別是在上述波段上，提出方法取得了更大的改善．VBM3D方法只能去除高斯白噪聲，所以在具有別的噪聲的波段上PSNR和SSIM值均比較低．其余4種方法雖然在不只有高斯白噪聲的波段上PSNR和SSIM值有所降低，但是總體效果還是比較好的．尤其是提出方法，不僅在含有3種噪聲類型的波段上保持較高的PSNR和SSIM值，在只含高斯白噪聲的波段的PSNR和SSIM值也比其余方法要高，這說明引入相似圖像塊，對于各種噪聲的去除效果均具有提升．

圖5 不同方法恢復(fù)的圖像的PSNR和SSIMFig.5 PSNR and SSIM of restoration results by different methods

表1含噪圖像和不同方法恢復(fù)結(jié)果的mPSNR和mSSIM

注： 加黑表示性能最優(yōu)．

表1為各種方法恢復(fù)圖像的mPSNR和mSSIM值．LRR方法的mPSNR和mSSIM均比LRMR高，說明了LRR模型的優(yōu)勢．聯(lián)合LRMR方法的mPSNR和mSSIM也比LRMR高，說明引入相似圖像塊對于低秩矩陣恢復(fù)也是有效果的．

3.1.2 實(shí)驗(yàn)1的參數(shù)分析

實(shí)驗(yàn)過程涉及到了5個參數(shù)的選?。?式(6)中的正則化參數(shù)λ，滑動窗的大小Wsize，圖像塊聚類的類數(shù)K，聯(lián)合LRR的相似圖像塊的數(shù)量p和式(19)中的參數(shù)τ．為了分析討論這5個參數(shù)對于實(shí)驗(yàn)結(jié)果的影響，對其進(jìn)行了更詳細(xì)的實(shí)驗(yàn)，采用mPSNR和mSSIM2個評價指標(biāo)來衡量實(shí)驗(yàn)結(jié)果．

1) 參數(shù)λ

式(6)所述的LRR模型表現(xiàn)的好壞很大程度上取決于參數(shù)λ，尤其在實(shí)際環(huán)境中，由于缺乏數(shù)據(jù)的先驗(yàn)知識，如何選取合適的λ是一個值得探討的問題．文獻(xiàn)[20]指出在沒有先驗(yàn)知識的情況下，λ按照以下式子選擇，LRR模型可以得到比較的結(jié)果，其中n為高光譜圖像的波段數(shù)．

(24)

圖6 參數(shù)λ對于恢復(fù)結(jié)果的影響Fig.6 Impact of parameter λ to restoration results

2) 滑動窗的大小Wsize

圖7給出了恢復(fù)結(jié)果和Wsize的關(guān)系，適中大小的Wsize對于恢復(fù)的結(jié)果比較好．當(dāng)Wsize比較小的時候，圖像塊的像素個數(shù)太少，樣本量不夠，導(dǎo)致恢復(fù)效果較差．當(dāng)Wsize比較大的時候，圖像塊過大時，它的紋理特征比較復(fù)雜，導(dǎo)致應(yīng)用LARK特征對圖像塊聚類時準(zhǔn)確度下降，進(jìn)而恢復(fù)效果變差．

圖7 滑動窗窗邊的大小Wsize對于恢復(fù)結(jié)果的影響Fig.7 Impact of edge size Wsize of sliding window to restoration results

3) 圖像塊聚類的類數(shù)K

K對于恢復(fù)結(jié)果的影響見表2．從表中可以看出恢復(fù)結(jié)果對于K的選取并不敏感，原因是只要選取適當(dāng)?shù)腒，使每一類的圖像塊個數(shù)大于等于p．那么對于每一個圖像塊都能在它所屬的類中找出和它最相似的p-1個塊，聯(lián)合進(jìn)行LRR，得到好的恢復(fù)效果．

4) 相似圖像塊的數(shù)量p

表3描述了p對于恢復(fù)效果的影響．當(dāng)p=1時，聯(lián)合LRR就退化成LRR了．從表3中，可以看出當(dāng)p從1跳躍到2時，恢復(fù)效果有了明顯提升，說明引入相似圖像塊的信息對于圖像恢復(fù)是有效果的．當(dāng)p大于5時，mPSNR和mSSIM有些許下降，說明聯(lián)合LRR的相似塊不是越多越好．圖像塊取多了，比較不相似的塊就有可能被引入到聯(lián)合LRR中，這會影響恢復(fù)效果．

表2 圖像塊聚類的類數(shù)K對于恢復(fù)結(jié)果的影響

表3 相似圖像塊的數(shù)量p對于恢復(fù)結(jié)果的影響

注： 加黑表示性能最優(yōu).

5) 參數(shù)τ

參數(shù)τ對于恢復(fù)結(jié)果的影響較小，從表4中可以看出τ的選取對mPSNR的影響不到0.1，對mSSIM的影響不到0.001．因?yàn)閳D像塊的權(quán)值是經(jīng)過歸一化的，說明τ對于各圖像塊的融合結(jié)果幾乎沒有影響．

表4 參數(shù)τ對于恢復(fù)結(jié)果的影響

注： 加黑表示性能最優(yōu).

總之，從上述的分析可以得出以下幾點(diǎn)結(jié)論： 參數(shù)λ可以在求出α后，通過式(24)得到；滑動窗的大小Wsize適中選取(圖7)；相似圖像塊的數(shù)量p在一定范圍內(nèi)選取(表3)；圖像塊聚類的類數(shù)K和參數(shù)τ，提出的方法對它們不敏感(表2和表4)．

3.1.3 實(shí)驗(yàn)2

實(shí)驗(yàn)2采用Washington D.C. Mall數(shù)據(jù)集，提出方法的實(shí)驗(yàn)參數(shù)選取如下： LRR的參數(shù)λ按照式(24)得出、滑動窗即圖像塊的大小Wsize=11、圖像塊聚類的類數(shù)K=10、相似圖像塊的數(shù)量p=5、式(19)的參數(shù)τ=10．

Washington D.C. Mall數(shù)據(jù)集的第81波段無噪圖像如圖8(a)所示，圖8(b)是添加了3種噪聲的含噪圖像，圖8(c)～(g)分別是使用不同方法恢復(fù)后的圖像．和實(shí)驗(yàn)1類似，VBM3D方法去除不了強(qiáng)烈的稀疏噪聲和壞死的線條，對圖像有過平滑的趨勢，存在涂抹痕跡．從圖中可以看出，LRR、聯(lián)合LRMR和聯(lián)合LRR方法對圖像的恢復(fù)效果均比LRMR好．圖8用紅色框出的區(qū)域受到嚴(yán)重噪聲影響，幾乎失去了空間結(jié)構(gòu)．VBM3D和LRMR方法均沒有重構(gòu)出這個結(jié)構(gòu)，LRR方法稍許恢復(fù)出了結(jié)構(gòu)但是比較模糊，而引入了相似塊的聯(lián)合LRMR和聯(lián)合LRR方法均恢復(fù)出了這個結(jié)構(gòu)，進(jìn)一步說明了LRR模型的優(yōu)勢和引入相似圖像塊的作用．

圖8 Washington D.C. Mall無噪圖、含噪圖和不同方法恢復(fù)結(jié)果的第81波段Fig.8 Band 81 of Washington D.C. Mall clean data, noisy data and restoration results by different methods

圖9為逐波段計算不同方法恢復(fù)后的圖像的PSNR和SSIM值的結(jié)果．在部分含多種噪聲的波段上，各種方法的PSNR和SSIM值均有不同程度的凹陷，原因和實(shí)驗(yàn)一相同．從圖中同樣可以看出，提出方法的各波段的PSNR和SSIM值均好于其他方法．在不只有高斯白噪聲的波段上，所有方法的PSNR和SSIM值均有所降低，VBM3D方法下降最嚴(yán)重．在只含有高斯白噪聲的波段上，VBM3D方法的PSNR和SSIM值和其余方法比較接近．LRMR和聯(lián)合LRMR方法在20～60波段表現(xiàn)不如LRR和聯(lián)合LRR方法穩(wěn)定，說明LRR模型比LRMR模型適合于高光譜數(shù)據(jù)．

圖9 不同方法恢復(fù)結(jié)果的PSNR和SSIMFig.9 PSNR and SSIM of restoration results by different methods

表5為各種恢復(fù)方法的mPSNR和mSSIM值．提出方法的mPSNR較VBM3D和LRMR方法分別提高4.16和2.90dB，mSSIM較VBM3D和LRMR方法分別提高0.07和0.03．說明本文提出的方法對于各種噪聲去除的效果都很穩(wěn)定．

3.2實(shí)際圖像實(shí)驗(yàn)

實(shí)際圖像實(shí)驗(yàn)選用了1992年6月拍攝的印第安納(Indiana)州西北部農(nóng)業(yè)區(qū)220個波段的AVIRIS高光譜航空影像，光譜分辨率和空間分辨率分別為10nm和20m，包含了145×145個像素，去除水吸收波段104～108，150～163，220后剩余200個波段，其中好幾個波段嚴(yán)重受到噪聲干擾．該數(shù)據(jù)集通常被用于高光譜圖像的分類實(shí)驗(yàn)中．這里，選取了9類地物，各類地物名稱和樣本數(shù)如表6所示．

噪聲的存在會對高光譜圖像分類的精度造成影響，所以可以用分類精度來評價對實(shí)際圖像的圖像恢復(fù)效果．本文使用支撐向量機(jī)(Support Vector Machine, SVM)作為分類器，分類總體精度(Overall Accuracy, OA)作為評價指標(biāo)(見表7)．每一類地物隨機(jī)選取10%的樣本作為訓(xùn)練樣本，其余作為測試樣本．Indiana數(shù)據(jù)集的地物真實(shí)與使用原始數(shù)據(jù)和不同方法處理后的圖像的分類結(jié)果見圖10．

表5 含噪圖像和不同方法恢復(fù)結(jié)果的mPSNR和mSSIM

注： 加黑表示性能最優(yōu).

表6 Indiana數(shù)據(jù)中地物各類別樣本數(shù)量

注： 加黑表示性能最優(yōu).

表7使用原始數(shù)據(jù)和不同方法恢復(fù)結(jié)果的分類總體精度

注： 加黑表示性能最優(yōu).

提出方法的實(shí)驗(yàn)參數(shù)選取如下： LRR的參數(shù)λ按照式(24)得出、滑動窗即圖像塊的大小Wsize=11、圖像塊聚類的類數(shù)K=10、相似圖像塊的數(shù)量p=5、式(19)的參數(shù)K=10．從OA值和分類結(jié)果圖可以知道，原始數(shù)據(jù)由于其受到噪聲影響嚴(yán)重，導(dǎo)致分類精度較差．在經(jīng)過不同方法處理后，分類精度均有所上升．由于實(shí)際高光譜圖像數(shù)據(jù)受到各種噪聲的影響，只能去除白噪聲的VBM3D方法對于提升分類精度的能力有限．而之后的4種方法能夠去除的噪聲類型比較多，所以能夠較大地提升高光譜圖像的分類精度．提出方法的分類精度較LRMR、VBM3D方法分別高了0.97%和7.64%，LRR和聯(lián)合LRMR方法的分類精度均比LRMR方法要高，說明LRR模型的確也可以較好地分離噪聲和信號，也說明了引入相似圖像塊在一定程度上是有助于實(shí)際圖像恢復(fù)的．而提出的方法綜合了LRR模型的優(yōu)勢和引入相似圖像塊的作用，對于實(shí)際高光譜圖像的恢復(fù)是有效的．

圖10 Indiana數(shù)據(jù)集的地物真實(shí)與使用原始數(shù)據(jù)和不同方法的恢復(fù)結(jié)果的分類結(jié)果Fig.10 Ground truth of Indiana dataset and results of classification using raw data and restoration results by different methods

4 結(jié) 論

本文提出了一種新的高光譜圖像的恢復(fù)方法．該方法先利用最小二乘估計出各個波段的噪聲水平，融合出單幅圖像．再對單幅圖像提取每一圖像塊的LARK特征，將圖像塊聚類．最后對于每一高光譜圖像塊，在類內(nèi)尋找相似圖像塊，組合在一起進(jìn)行聯(lián)合LRR．其中由于最小二乘估計的結(jié)果已經(jīng)能比較好地表示數(shù)據(jù)子空間，所以從中選擇相應(yīng)塊的最小二乘估計作為LRR的字典．更進(jìn)一步，該方法不僅引入了更符合高光譜圖像特性的LRR模型，還突破了典型低秩恢復(fù)方法只利用單一圖像塊信息的缺點(diǎn)，利用相似圖像塊更精確地恢復(fù)出高光譜圖像的低秩結(jié)構(gòu)．實(shí)驗(yàn)結(jié)果表明，所提出方法能夠在保留圖像細(xì)節(jié)的前提下去除多種類型的噪聲，在目視效果和客觀評價指標(biāo)上均優(yōu)于目前的高光譜圖像恢復(fù)方法．這對于處理受到各種噪聲干擾的高光譜圖像來說具有重要的實(shí)際意義．

[1] KOPSINIS Y, MCLAUGHLIN S. Development of EMD-based denoising methods inspired by wavelet thresholding [J].IEEETransSignalProcess, 2009,57(4)： 1351-1362.

[2] WEN Y W, NG M K, HUANG Y M. Efficient total variation minimization methods for color image restoration [J].IEEETransImageProcess, 2008,17(11)： 2081-2088.

[3] DABOV K, FOI A, KATKOVNIK V,etal. Image denoising by sparse 3-D transform-domain collaborative filtering [J].IEEETransImageProcess, 2007,16(8)： 2080-2095.

[4] LETEXIER D, BOURENNANE S. Noise removal from hyperspectral images by multidimensional filtering [J].IEEETransGeosciRemoteSens, 2008,46(7)： 2061-2069.

[5] LIU X, BOURENNANE S, FOSSATI C. Nonwhite noise reduction in hyperspectral images [J].IEEEGeosciRemoteSensLett, 2012,9(3)： 368-372.

[6] DABOV K, FOI A, EGIAZARIAN K. Video denoising by sparse 3D transform-domain collaborative filtering [C]∥IEEE Conf Signal Process, European, 2007： 145-149.

[7] LU X, WANG Y, YUAN Y. Graph-regularized low-rank representation for destriping of hyperspectral images [J].IEEETransGeosciRemoteSens, 2013,51(7)： 4009-4018.

[8] ZHANG H, HE W, ZHANG L,etal. Hyperspectral image restoration using low-rank matrix recovery [J].IEEETransGeosciRemoteSens, 2014,52(8)： 4729-4743.

[9] HE W, ZHANG H, ZHANG L,etal. Total-variation-regularized low-rank matrix factorization for hyperspectral image restoration [J].IEEETransGeosciRemoteSens, 2016,54(1)： 178-188.

[10] LIU G, LIN Z, YAN S,etal. Robust recovery of subspace structures by low-rank representation [J].IEEETransPatternAnalMachIntell, 2013,35(1)： 171-184.

[11] SUMARSONO A, DU Q. Low-rank subspace representation for estimating the number of signal subspaces in hyperspectral imagery [J].IEEETransGeosciRemoteSens, 2015,53(11)： 6286-6292.

[12] XU Y, WU Z, LI J,etal. Anomaly detection in hyperspectral images based on low-rank and sparse representation [J].IEEETransGeosciRemoteSens, 2016,54(4)： 1990-2000.

[13] IORDACHE M D, BIOUCAS-DIAS J M, PLAZA A. Sparse unmixing of hyperspectral data [J].IEEETransGeosciRemoteSens, 2011,49(6)： 2014-2039.

[15] TAKEDA H, FARSIU S, MILANFAR P. Kernel regression for image processing and reconstruction. [J].IEEETransImageProcess, 2007,16(2)： 349-366.

[16] CHATTERJEE P, MILANFAR P. Clustering-based denoising with locally learned dictionaries [J].IEEETransImageProcess, 2009,18(7)： 1438-1451.

[17] GAO L, DU Q, ZHANG B,etal. A comparative study on linear regression-based noise estimation for hyperspectral imagery [J].IEEEJSelTopicsApplEarthObservRemoteSens, 2013,6(2)： 488-498.

[18] LLOYD S P. Least squares quantization in PCM [J].IEEETransonInfTheory, 1982,28(2)： 129-137.

[19] WANG Z, BOVIK A C, SHEIKH H R,etal. Image quality assessment： From error visibility to structural similarity [J].IEEETransImageProcess, 2004,13(4)： 600-612.

[20] LIU G, XU H, TANG J,etal. A deterministic analysis for LRR [J].IEEETransPatternAnalMachIntell, 2016,38(3)： 417-430.

Abstract： Hyperspectral imagery(HSI) is unavoidably corrupted by various noises, such as Gaussian white noise, impulse noise, dead lines and so on, in the acquisition process. In order to ensure the performance of subsequent applications, HSI restoration is an essential preprocessing process. A new method for HSI restoration is proposed in this paper with two main contributions： 1) introducing low-rank representation model into the HSI restoration; 2) utilizing similar image blocks for joint low-rank representation. Experimental results on synthetic and real hyperspectral datasets demonstrate that the proposed method can effectively remove various noises and reconstruct the details of image, thus, it should also be used as an effective technique for HSI restoration.

Keywords： hyperspectral imagery; low-rank representation; similar image blocks; joint; restoration

rank(X)≤min(rank(A),rank(S))，

rank(X)≤min(rank(A),rank(S))≤P·min(n,b)，

min‖X‖*+λ‖E‖1
s.t.Y=X+E，

HyperspectralImageryRestorationBasedonLow-RankRepresentation

CAI Rongrong1, 2, WANG Bin1, 2

(1.KeyLaboratoryforInformationScienceofElectromagneticWaves(MoE),FudanUniversity,Shanghai200433,China; 2.ResearchCenterofSmartNetworksandSystems,SchoolofInformationScienceandTechnology,FudanUniversity,Shanghai200433,China)

A

0427-7104(2017)03-0346-13

2016-06-07

國家自然科學(xué)基金(61572133)；北京師范大學(xué)地表過程與資源生態(tài)國家重點(diǎn)實(shí)驗(yàn)室開放基金(2017-KF-19)

蔡榮榮(1991—)，男，碩士研究生；王 斌，男，教授，通信聯(lián)系人，E-mail： wangbin@fudan.edu.cn.