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

蔡榮榮，王 斌

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

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

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

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

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

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

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

1 高光譜圖像的低秩模型

1.1無(wú)噪高光譜圖像的低秩特性

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

X=AS，

(1)

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

(2)

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

(3)

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

1.2低秩表示模型

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

Y=X+N.

(4)

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

(5)

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

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

(6)

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

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

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

(7)

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

(8)

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

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

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

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

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

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

2.2.1 估計(jì)各個(gè)波段的噪聲，將高光譜圖像融合成單幅圖像

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

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

yi=Y?iβi+ξi，

(9)

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

(10)

然后可得到第i個(gè)波段的信號(hào)和噪聲的估計(jì)如下：

(11)

(12)

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

(13)

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

(14)

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

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

(15)

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

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

2.1.2 相似塊的搜尋

圖3 相似的圖像塊的聚類(lèi)過(guò)程Fig.3 Clustering of the similar blocks

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

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

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

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

(16)

(17)

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

(18)

(19)

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

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

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

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

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

在模擬實(shí)驗(yàn)中，對(duì)2個(gè)數(shù)據(jù)集均添加以下3種噪聲： 每個(gè)波段都添加σ=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條壞死的線(xiàn)條，模擬實(shí)際中的傳感器故障．

本文采用2個(gè)評(píng)價(jià)指標(biāo)，峰值信噪比(Peak Signal-to-Noise Ratio, PSNR)和結(jié)構(gòu)相似性(Structural Similarity, SSIM)指標(biāo)[19]定量地對(duì)恢復(fù)結(jié)果進(jìn)行評(píng)價(jià)．對(duì)無(wú)噪圖像和恢復(fù)圖像逐波段計(jì)算PSNR和SSIM值，然后取平均得到平均PSNR(Mean PSNR,mPSNR)和平均SSIM(Mean SSIM,mSSIM)值．各個(gè)指標(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、滑動(dòng)窗即圖像塊的大小Wsize=11、圖像塊聚類(lèi)的類(lèi)數(shù)K=10、相似圖像塊的數(shù)量p=5、式(19)的參數(shù)τ=10．Pavia University數(shù)據(jù)集的第61波段無(wú)噪圖像如圖4(a)所示，圖4(b)是添加了3種噪聲的含噪圖像，圖4(c)～(g)分別是使用VBM3D、LRMR、LRR、聯(lián)合LRMR和聯(lián)合LRR方法恢復(fù)后的圖像．從圖中可以看出，VBM3D方法對(duì)于強(qiáng)烈的稀疏噪聲和壞死的線(xiàn)條去除效果較差，因?yàn)樗褂玫哪Ｐ椭荒苋コ咚拱自肼暎愀獾氖?，由于VBM3D使用聯(lián)合維納濾波對(duì)圖像進(jìn)行處理，導(dǎo)致它對(duì)圖像有過(guò)平滑的趨勢(shì)，例如圖4(c)右上角的模糊說(shuō)明了這一點(diǎn)．LRMR、LRR、聯(lián)合LRMR和聯(lián)合LRR方法由于使用低秩模型，不僅可以去除高斯白噪聲，也可以去除強(qiáng)烈的稀疏噪聲和壞死的線(xiàn)條，恢復(fù)出比較好的圖像．不過(guò)，從圖中可以看出，LRMR方法不能完全去除所有噪聲，LRR方法相較LRMR方法的恢復(fù)效果要好一些，說(shuō)明LRR模型比LRMR模型更適合于高光譜圖像的恢復(fù)．但是圖4(e)紅色圈出的區(qū)域受到噪聲影響，LRR方法的恢復(fù)效果并沒(méi)有很好，存在一些涂抹痕跡．提出的方法將相似圖像塊組合到一起進(jìn)行聯(lián)合LRR，引入了相似圖像塊的信息增加了低秩矩陣的樣本容量，更好的幫助恢復(fù)出不含噪聲的高光譜圖像．從圖4(g)看出所提出方法不僅可以去除所有類(lèi)型的噪聲，同時(shí)還比較好地恢復(fù)出了圖像的空間結(jié)構(gòu)，如紅色圈出的區(qū)域的恢復(fù)效果就好于LRR方法．另外，聯(lián)合LRMR相比LRMR方法恢復(fù)效果的提升也能說(shuō)明引入相似圖像塊的作用．

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

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

圖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高，說(shuō)明了LRR模型的優(yōu)勢(shì)．聯(lián)合LRMR方法的mPSNR和mSSIM也比LRMR高，說(shuō)明引入相似圖像塊對(duì)于低秩矩陣恢復(fù)也是有效果的．

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

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

1) 參數(shù)λ

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

(24)

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

2) 滑動(dòng)窗的大小Wsize

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

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

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

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

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

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

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

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

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

5) 參數(shù)τ

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

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

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

總之，從上述的分析可以得出以下幾點(diǎn)結(jié)論： 參數(shù)λ可以在求出α后，通過(guò)式(24)得到；滑動(dòng)窗的大小Wsize適中選取(圖7)；相似圖像塊的數(shù)量p在一定范圍內(nèi)選取(表3)；圖像塊聚類(lèi)的類(lèi)數(shù)K和參數(shù)τ，提出的方法對(duì)它們不敏感(表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)得出、滑動(dòng)窗即圖像塊的大小Wsize=11、圖像塊聚類(lèi)的類(lèi)數(shù)K=10、相似圖像塊的數(shù)量p=5、式(19)的參數(shù)τ=10．

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

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

圖9為逐波段計(jì)算不同方法恢復(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)定，說(shuō)明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．說(shuō)明本文提出的方法對(duì)于各種噪聲去除的效果都很穩(wěn)定．

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

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

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

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

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

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

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

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

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

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

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

4 結(jié) 論

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

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

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