宋毅梁

（四川大學(xué)計(jì)算機(jī)學(xué)院，成都610065）

0 引言

高光譜圖像是一種被廣泛運(yùn)用在植被勘探、地質(zhì)檢測(cè)、農(nóng)業(yè)、軍事、天氣等領(lǐng)域的成像手段。不同材質(zhì)的物體其分子和原子的排列決定了其本質(zhì)的特征，而當(dāng)電磁波攝入到不同的物質(zhì)上，其物質(zhì)的分子、原子排列特征、震動(dòng)等會(huì)讓物體呈現(xiàn)不同的光譜吸收與反射的特征。而高光譜成像便是運(yùn)用這一特點(diǎn)，運(yùn)用大量的窄波波段的電磁波（<10nm），得到一個(gè)帶有大量不同波段數(shù)據(jù)信息的三維數(shù)據(jù)塊。在空間圖像維的圖像與一般的圖像沒(méi)有分別，而對(duì)于光譜維來(lái)說(shuō)，圖像的每一個(gè)像元都有大量的不同波段的光譜信息。

在日常采集高光譜圖像信息的過(guò)程中，因?yàn)椴杉瘯r(shí)采集通道可能因?yàn)檫@個(gè)過(guò)程中產(chǎn)生的熱噪聲使得原本的圖像被噪聲污染，而高光譜圖像數(shù)據(jù)在有噪聲污染的情況下會(huì)給后續(xù)的數(shù)據(jù)分析帶來(lái)很大的障礙。因此去除這類加性噪聲會(huì)為之后的數(shù)據(jù)分析帶來(lái)便利。

傳統(tǒng)的噪聲抑制方法，把張量數(shù)據(jù)的二維圖像以光譜維展開(kāi)成二維的矩陣，之后運(yùn)用矩陣的低秩先驗(yàn)的特點(diǎn)，對(duì)圖像進(jìn)行PCA等降維去噪。而這種去噪方式并沒(méi)有考慮到圖像的光譜維之間的聯(lián)系，會(huì)損失掉不少的圖像結(jié)構(gòu)信息。因此基于張量代數(shù)的去噪方式是當(dāng)前高光譜圖像噪聲抑制的重點(diǎn)課題。

目前主要的基于張量代數(shù)的算法有Tucker分解和CP分解兩種方式?；赥ucker3的噪聲抑制方法包括低秩張量近似LRTA[1]、遺傳核張量分解[2]，以及多維維納濾波[3]。基于CP分解的噪聲抑制方法有并行因子分析（PARAFAC）[4]和秩-1張量分解[5]。除此之外基于高階張量圖像數(shù)據(jù)的小波分解噪聲抑制方法因其良好的細(xì)節(jié)保留特性也被大量提出[6-8]。這些噪聲處理方式對(duì)于受到高斯白噪聲影響的HSI圖像的信號(hào)空間還原是比較高效的。但是在遇到完全丟失像素信息的椒鹽噪聲、死像素等則不能有效處理。而近年以來(lái)提出的魯棒性主成分分析（RPCA）[9]，這種方法已經(jīng)被驗(yàn)證可以在傳統(tǒng)圖像上對(duì)稀疏噪聲進(jìn)行有效的噪聲抑制還原原始信號(hào)。而這種方法的高階推廣也被人提出并且用在稀疏噪聲的HSI圖像的噪聲抑制上。TRPCA[10]用 || ||Z*張量的核范數(shù)TNN被用于描述張量的秩。除此之外還有基于TV全變差等方法的提出。本文針對(duì)目前已經(jīng)出現(xiàn)的基于張量分解的噪聲抑制方法進(jìn)行了分析和效果對(duì)比。

1 張量基本概念與分解模型

1.1 張量的基本概念

（1）張量的mode-n展開(kāi)

張量的mode-n展開(kāi)是把張量的第n個(gè)mode的纖維（fibers）作為展開(kāi)矩陣的列向量然后把這些mode的纖維順序排列。一個(gè)張量X∈R I1xI2xI3的mode維展開(kāi)為x(n)∈R I n xM n，其中Mn=I1X….I n-1XI n+1….I N。

（2）張量的mode-n乘積

張量X∈R I1X…XI N與矩陣U∈R JXI n的mode-n的乘積被表示為XxnU，實(shí)質(zhì)上是張量的n-mode的纖維（fi?bers）與矩陣相乘。因此有公式：

（3）秩-1張量

秩-1張量表示N個(gè)矢量進(jìn)行外積得到的張量，而外積之后的張量階數(shù)取決于外積矢量的個(gè)數(shù)，其形式為：

此為N階的秩-1張量。

1.2 張量分解模型

（1）Tucker分解

Tucker分解是奇異值分解（SVD）的高階推廣，一個(gè)三階張量的Tucker分解會(huì)把張量X∈R I1×I2×I3會(huì)被分解為一個(gè)核張量G∈R M×N×J和三個(gè)分別對(duì)應(yīng)三個(gè)mode做乘積的系數(shù)矩陣A∈R I1×M、B∈R I2×N與C∈R I3×J。其分解的表示形式為：

此外對(duì)張量的mode積而言，對(duì)對(duì)應(yīng)的mode做乘積相當(dāng)于在這個(gè)mode上做投影。

（2）CP分解

CANDECOMP/PARAFAC分解（簡(jiǎn)稱CP分解）這種分解會(huì)把張量分解為秩-1張量的和。而張量的CP-秩被表示為分解的秩-1張量的個(gè)數(shù)。

2 張量噪聲抑制算法

2.1 基于Tucker分解的LRTA算法

LRTA算法運(yùn)用主成分分析（PCA）的思想擴(kuò)展到高階的HIS三階張量圖像數(shù)據(jù)上，通過(guò)估計(jì)信號(hào)空間的維度在降維后還原原始圖像。令加性噪聲在對(duì)每個(gè)mode的投影中與信號(hào)空間分離以此在達(dá)到收斂條件之后得到估計(jì)信號(hào)的結(jié)果。

一般的高斯噪聲模型認(rèn)為采集圖像所攜帶的噪聲是一種加性噪聲，而想要去除高斯白噪聲這樣的加性噪聲，重要的就是把采集到的數(shù)據(jù)投影在算法估計(jì)的信號(hào)空間上。在PCA里，這樣的投影工具是一個(gè)矩陣它把二維圖像數(shù)據(jù)降維投影得到最逼近原始信號(hào)的數(shù)據(jù)。而對(duì)于三階張量數(shù)據(jù)來(lái)說(shuō)，對(duì)各個(gè)mode進(jìn)行投影則是在各個(gè)mode上進(jìn)行mode-n乘積。因此問(wèn)題的關(guān)鍵在于得到每個(gè)mode對(duì)應(yīng)的投影矩陣p（i）。估計(jì)的信號(hào)空間數(shù)據(jù)被表示為：

而每個(gè)mode對(duì)應(yīng)的信號(hào)空間的投影需要先知道每個(gè)mode信號(hào)空間的維數(shù)K n，這些維數(shù)的確定依賴于張量Tucker-秩。(K1,…,K N)每個(gè)mode的秩被作為衡量Tucker分解的張量秩。而基于高階PCA思想的這種去噪方式所要處理的其實(shí)是如下的優(yōu)化問(wèn)題：

顯然求張量的秩是一個(gè)非凸的問(wèn)題，因此要變?yōu)榭汕蠼獾膬?yōu)化形式：

張量的核范數(shù)是張量沿每個(gè)mode展開(kāi)之后求得的奇異值之和。而最小化核張量SNN則是滿足低秩先驗(yàn)的條件。

LRTA算法描述：

（1）input：高光譜HSI張量數(shù)據(jù)R，以及被估計(jì)的投影空間(K1,…,K N)

2.2 基于CP分解的噪聲抑制算法

基于CP分解把張量分解為秩-1張量的想法，提出把張量分解為帶權(quán)值的秩-1張量。認(rèn)為張量的信號(hào)空間是一組高度相關(guān)的秩-1張量組成的，噪聲則在光譜域的相關(guān)性較小。

通過(guò)CP分解后得到的帶權(quán)值的秩-1張量：

其中M是采集信號(hào)R的CP分解所需要的秩-1張量的總個(gè)數(shù)，而αr則是每個(gè)秩-1張量所占信號(hào)空間表達(dá)的權(quán)值。要通過(guò)R求解未被噪聲污染的信號(hào)關(guān)鍵在于求解張量的秩。而與矩陣不同，張量的秩在目前還沒(méi)有明確的定義形式[11]。

不同于Tucker的各個(gè)mode上的秩表示張量的秩，而是分解后秩-1張量的個(gè)數(shù)被表示為了張量的秩序。又因?yàn)樵摲椒ㄕJ(rèn)為噪聲以加性信號(hào)為前提而高斯白噪聲在信號(hào)空間里是并不相關(guān)的因此存在的權(quán)值很小，而信號(hào)空間因?yàn)樵诠庾V空間上高度相關(guān)因此權(quán)值占得比重更大，得到的被估計(jì)的信號(hào)空間可以表示為：

K表示的是代表信號(hào)空間的秩-1張量的個(gè)數(shù)，也是權(quán)值最大的K個(gè)秩-1張量。

而W、V同理可以根據(jù)上邊的方式迭代跟新，最終當(dāng)這一次迭代減去上一次迭代的值滿足條件時(shí)就可以跳出得到估計(jì)的信號(hào)張量數(shù)據(jù)X?。

R1TD算法描述：

（1）input：高光譜HSI張量數(shù)據(jù)R，以及被估計(jì)的秩-1張量個(gè)數(shù)K，最大迭代數(shù)M

（2）循環(huán)t=1→M

2.3 Tucker-秩與CP-秩的估計(jì)方法

（1）AIC赤池信息準(zhǔn)則

AIC是一種優(yōu)秀的模擬統(tǒng)計(jì)信息的一種標(biāo)準(zhǔn)[12]，因?yàn)閺埩恐鹊亩x不像矩陣一樣明確，因此AIC用于估計(jì)擬合各個(gè)mode上的統(tǒng)計(jì)信息，即Tucker-秩。

（2）CP分解秩-1張量的估計(jì)

秩-1張量定義的秩不能通過(guò)計(jì)算每個(gè)mode的統(tǒng)計(jì)信號(hào)擬合來(lái)直接得出。需要利用信噪比得到：

Ki是由AIC得到的每個(gè)mode對(duì)應(yīng)的統(tǒng)計(jì)擬合SNR是張量的信噪比，而在輸入條件并不知道信號(hào)空間的時(shí)候，信噪比可以通過(guò)求統(tǒng)計(jì)算法的結(jié)果得到近似的評(píng)估。

（1）在張量中準(zhǔn)備一個(gè)移動(dòng)窗口，5×5、7×7等，計(jì)算部分圖像區(qū)域的方差，得到最小的方差，把這個(gè)方差作為噪聲的統(tǒng)計(jì)信息估計(jì)。

（2）求出信號(hào)統(tǒng)計(jì)信息的估計(jì)：

3 Tucker與CP的抑制效果對(duì)比

本節(jié)運(yùn)用的HIS圖像數(shù)據(jù)為Washington DC Mal,Indian Pine高光譜圖像數(shù)據(jù)集合，對(duì)圖像添加加性的高斯白噪聲在不同波段的圖像分別應(yīng)用LRTA和R1TD兩種基于張量分解的噪聲抑制方法進(jìn)行去噪聲，通過(guò)對(duì)圖片細(xì)節(jié)與峰值信噪比PSNR等圖像效果指標(biāo)進(jìn)行評(píng)估對(duì)比。

圖1

圖片從左到右依次為（1）原始圖像211光譜波段（2）加入高斯白噪聲之后的211光譜段圖像（3）LRTA處理之后的211波段圖像PSNR為29.603（4）R1TD處理之后的211波段圖像PSNR為33.011

圖像在進(jìn)行了張量的低秩分解處理之后噪聲得到抑制，不過(guò)兩種方法都基于低秩的去噪方法因此在細(xì)節(jié)上并沒(méi)有得到很好的保留。

4 結(jié)語(yǔ)

本文對(duì)主流的兩種張量分解的噪聲抑制方法進(jìn)行了闡述并分析了其工作原理，并且利用Washington DC Mal，Indian Pine高光譜數(shù)據(jù)集進(jìn)行了仿真。利用PSNR等方法對(duì)噪聲抑制之后的光譜波段進(jìn)行了對(duì)比。通過(guò)實(shí)驗(yàn)與對(duì)比之后發(fā)現(xiàn)在進(jìn)行這兩種方法之后，去噪過(guò)程實(shí)際上是在降低數(shù)據(jù)的維度，不可避免地會(huì)造成細(xì)節(jié)的丟失因此未來(lái)的工作可以如何對(duì)細(xì)節(jié)盡可能地保留等課題上。