PCA聯合字典的稀疏系數NMF融合

孫小芳

(閩江學院地理科學系，福州 350121)

0 引言

稀疏表達又稱為稀疏分解，用冗余函數構造稀疏字典，在影像融合中對高、低空間分辨率影像分別計算稀疏系數，通過算法生成融合影像的融合稀疏系數，進而結合稀疏字典完成影像稀疏融合[1]。基于稀疏思想的融合中，稀疏字典與融合稀疏系數的生成是2個研究重點。

稀疏字典從用途來看，可分成待融合影像的稀疏字典與融合結果影像的稀疏字典。待融合影像是指用于融合的高、低空間分辨率影像，該類字典的來源有2種：①將標準正交基級聯得到超完備字典，通常包括傅里葉變換、小波變換、離散余弦變換、Gabor濾波、曲線波以及輪廓波等，例如將同一空間位置對應的同方向跨尺度小波基函數的線性組合作為新的基函數[2]；②通過待融合影像的樣本學習來構造過完備字典，參與字典學習的樣本類型包括隨機選擇樣本[3]、基于影像分割或分類所產生的區域影像建立樣本[4]、基于隨機共振和自適應的稀疏域選擇樣本[5]及選擇純像元建立字典學習樣本[6]等。融合結果影像的字典除了各種標準正交基函數外，還包括以下3類：①采用多光譜影像字典，例如在MODIS與ETM+影像融合中，由MODIS影像提供融合影像的字典[6]；②采用高空間分辨率影像字典，例如SPOT與TM影像融合中，融合字典由SPOT影像字典提供[7]；③基于多光譜字典與高空間分辨率影像字典利用規則生成融合字典，該種融合字典生成的方法包括正則項建立優化函數求融合字典[8-10]、對各種聚類子字典采用主成分分析(principal component analysis，PCA)方法構造融合字典[11]、利用融合字典與高空間分辨率影像字典存在的權重關系及融合字典與多光譜字典存在的模糊濾波關系構建融合字典[8]以及隨機選擇高空間分辨率影像與多光譜影像樣本建立聯合字典[3]。

目前字典的生成多數直接來源于影像樣本，但由于遙感影像的空間分辨率限制，在地物復雜地區大多數影像存在著混合像元問題，這使得字典的精度受到一定影響。2014年Huang等[6]利用MODIS影像中的純像元建立字典，但是MODIS影像中各種地類的純像元個數較少，影響了字典的應用效果。為了減少混合像元對字典建立的影響，本研究在像元線性分解影像中利用在線字典學習法建立字典，通過提取分解影像字典與全色影像字典的主成分分量建立PCA聯合稀疏字典，融合影像的稀疏系數采用非負矩陣分解(nonnegative matrix factor，NMF)融合多光譜影像與全色影像的稀疏系數生成。

1 基本理論

1.1 稀疏分解

利用盡量少的原子影像塊與非零值稀疏系數來完全或近似地表達原始影像的方式，就是影像稀疏表達。該種表達方式將影像投影到由稀疏字典組成的特征空間，影像的信息集中在較少的原子影像塊中，非零值稀疏系數表明了影像的內部結構及特征[12]。基于高度非線性逼近理論的稀疏表達公式中包含2個系數：一是根據信號的特點構造的原子庫即過完備字典D; 二是從字典中找到最佳線性組合所用到的稀疏系數α[13-14]。

y=D·α

(1)

由于實際處理的信號通常帶有噪聲，影像稀疏分解過程實際上是一種逼近過程，即

(2)

式中：yM為原始信號y的逼近信號；xr為殘差分量；dN為給定稀疏字典D∈RM×N中的一個原子；〈RNxrdN〉為信號在dN上的投影，構成信號在稀疏字典D上的稀疏系數α；RMx對應殘差分量xr。

為了求解最優α，且滿足殘差分量xr達到最小，將式(2)轉換為

(3)

式中：‖α‖0為α的l0范數；αN為給定稀疏字典中的一個稀疏系數；ΩΜ為稀疏字典大小；ε為誤差總和。式(3)表明在最小均方誤差約束條件下，求解稀疏系數的最少個數。由于l0范數的非凸性，為了求解這個NP-hard問題，學者通常用l1范數近似代替l0范數。

1.2 NMF融合

利用NMF算法將原始影像V分解成基向量W和權重系數矩陣H，當W的秩r與數據集的特征空間維數一致時，那么得到的W是對V最有效的體現，既反映了影像V中最基本的特征，同時有效抑制了影像的噪聲。利用NMF的這種特性，可以完成NMF的影像融合。將待融合的影像按行優先的方式存儲成V，即

(4)

式中：Ii為待融合的影像，i=1，2，…，m；M和N分別為每一幅待融合影像的行數和列數，并令n=M×N，則有V∈Rn×m。V可近似分解為非負矩陣Wn×r與非負矩陣Hr×m的乘積[15]。本文的稀疏系數融合是采用逐個多光譜波段與全色波段進行融合，所以選取r=1，即W的維數為MN×1。融合時采用標準梯度下降法迭代規則求解W和H，由于V=W×H+δ，在每次迭代后計算誤差δ，若δ小于設定的值，則迭代終止。

2 融合算法

本文融合算法的流程如圖1所示。

本研究改進了稀疏字典的生成方法，提出了利用不同地物的線性分解影像作為字典的來源影像。進行線性光譜分解時，所用到的純像元來源于純像元指數(pure pixel index，PPI)計算，再利用N維可視化器從PPI指數計算的像元中提取出感興趣區域(region of interest，ROI)純像元，作為線性光譜分解的端元。利用在線字典學習算法分別對各類地物分解影像進行處理，得到各類地物分解影像的稀疏字典，將每一幅分解影像稀疏字典與全色稀疏字典分別進行PCA處理并提取第一主成分分量作為PCA稀疏字典，并將8個PCA稀疏字典生成PCA聯合稀疏字典。利用PCA聯合稀疏字典與正交匹配追蹤(orthogonal matching pursuit，OMP)算法計算多光譜影像與全色影像的稀疏系數，將每一幅多光譜影像的稀疏系數與全色影像的稀疏系數進行NMF融合，得到該波段融合影像的稀疏系數，融合影像的稀疏系數與PCA聯合字典重構生成融合影像。

3 實驗與分析

采用2012年8月8日福建賽場的WorldView-2影像數據作為實驗數據。該數據擁有1個0.5 m空間分辨率的全色波段和8個2 m空間分辨率的多光譜波段，具體包括：海岸波段B1、紅光波段B2、藍光波段B3、紅邊波段B4、綠光波段B5、近紅外1波段B6、黃光波段B7和近紅外2波段B8。數據已完成輻射校正與幾何糾正，可以進行影像融合處理，實驗影像的大小為512像元×512像元。

在本研究中，稀疏字典的生成是基于地物分解影像。在利用PPI指數計算時，對所設定的參數進行調整，最終確定迭代次數為10 000次，閾值參數設為12，得到3 016個像元。PPI指數計算的結果給出了像元作為純像元的潛在性，但并沒有對純像元的類別做出判斷，這時就需要利用N維可視化器提取各類別的純像元。通過觀看各角度旋轉，找出各方向散點圖的尖角，選出作為各種類別的純像元。本研究最終確定8種類別地物：水體、裸土、農田、農作物、林地、不可滲透表面、白色屋頂建筑和藍色屋頂建筑。將8類地物的純像元樣本作為線性分解的端元輸入，得到8種地物的像元分解影像。

3.1 確定稀疏字典矩陣大小

進行稀疏字典計算時，先要確定稀疏字典的大小。本研究利用多光譜影像的稀疏分解與重構來探討適合的稀疏字典矩陣大小。利用在線字典學習算法分別對8類地物的分解影像進行處理，得到8幅影像的稀疏字典，并將其聯合起來，組成聯合稀疏字典。一個字典原子的大小為8像元×8像元，在考慮計算機運行能力與影像稀疏分解及重構精度問題的基礎上，探討合適的稀疏字典矩陣大小。設定各幅多光譜影像的字典個數分別為10，20，30，40，50和60，則對應的聯合字典的個數分別為80，160，240，320，400和480，利用該聯合字典與OMP算法分別計算各多光譜波段的稀疏系數，利用稀疏系數與聯合字典重構各多光譜波段，探討影像重構均方根誤差(root mean square error，RMSE)與字典個數的關系，從而確定本研究的字典個數。

表1給出6種聯合稀疏字典的參數及8個波段重構的平均RMSE，可以看出重構影像平均RMSE從80個字典的0.099下降到480個字典的0.054。說明字典個數對重構影像有著相關性影響，這從圖2中所表達的8個波段的重構RMSE中也得到反映。綜合考慮運算效率與影像的重構精度，本文選擇的聯合字典大小為480，即每個分解影像提供60個字典用于建立聯合字典。

表1 字典參數與重構RMSE

圖2 各波段稀疏重構RMSE

3.2 PCA聯合稀疏字典

利用在線字典學習法生成8個分解影像的字典與全色影像字典。利用PCA算法計算各分解影像稀疏字典與全色稀疏字典，取第一主成分分量作為PCA聯合稀疏字典。圖3分別為直接利用8個分解影像的字典組成的聯合字典以及依次對8個分解影像進行處理后得到的PCA聯合稀疏字典。

(a) 聯合字典

(b) PCA聯合稀疏字典

影像的方差可以反映影像的信息量，比較PCA聯合稀疏字典與聯合字典的方差(表2)可知，經過PCA處理得到的各影像字典及影像聯合字典的方差值均大于未經PCA處理的字典，說明經過PCA處理的字典融合進了全色影像字典的主要信息，從而保證PCA聯合稀疏字典能代表影像中主要特征。

表2 聯合字典與PCA聯合稀疏字典的方差比較

3.3 融合影像評價

利用PCA聯合稀疏字典完成多光譜波段(圖4(a))與全色波段(圖4(b))的OMP計算，分別得到各個多光譜波段與全色波段的稀疏系數α。利用NMF融合算法計算每個多光譜波段與全色波段的融合系數α，即各個波段融合影像的系數α，將該數值乘以PCA聯合稀疏字典，得到各個融合波段，如圖4(c)。在完成本文提出的融合算法后，利用聯合字典進行多光譜與全色稀疏分解，將得到的系數α利用NMF算法進行融合，得到聯合字典的NMF融合影像，如圖4(d)；同時為了對比分析本文方法與傳統方法的融合效果，分別完成小波融合與經典的PCA影像融合，如圖4(e)與圖4(f)。從目視效果上看，本文提出的融合算法與圖4(a)多光譜影像色調較一致，保證了多光譜信息。融合影像的紋理層次豐富與圖4(b)全色影像清晰度接近，信息量增加，同時具有多光譜與高空間分辨率影像的特征。對比其他3種融合算法，本文方法所融合的地物光譜表現與實際情況接近，田地、水體和居民區更清晰可辨，更有力地表現了細節。

(a) 多光譜影像 (b) 全色影像 (c) PCA聯合稀疏字典NMF融合影像

(d) 聯合字典NMF融合影像 (e) 小波融合影像 (f) PCA融合影像

通過計算4幅融合影像的5種定量評價指標來比較本文方法與其他經典方法的融合效果，具體如表3所示。

表3 融合評價指標

信息熵可以評價影像中信息量的多少，值越高，表示影像的信息越多。可以看出PCA聯合稀疏字典的NMF融合影像具有更多的信息，信息熵值在4種方法中最高。清晰度與空間頻率是通過相鄰像元的差別來反映影像的空間變化頻率，表現為細節的反映能力。相較于聯合字典NMF融合、小波融合與PCA融合，本文方法的這2個參數指標明顯更高，說明本文方法提高了融合影像的紋理細節信息，更加清晰地反映了地物的特征。扭曲程度與偏差指數反映融合影像與多光譜影像的差別程度，值越大說明與原始多光譜影像的差異越大。本文方法的這2個指標均小于其他3種融合方法，說明本文方法較好地保持了影像的光譜信息。

4 結論

1)本研究提出分解地物影像字典與全色影像字典的第一主成分分量構成PCA聯合稀疏字典。所生成的PCA聯合稀疏字典既能包含多光譜影像分解地物特征，也包含高空間分辨率影像特征。

2)利用重構影像與原始影像的RMSE探討合適的字典個數。隨著字典個數的增加，重構影像的誤差有明顯減小的趨勢。考慮到重構影像的效率與計算機運算的限制，最終確定稀疏字典個數為480。

3)利用NMF融合確定融合影像的稀疏系數。該稀疏系數最大限度地保留2幅原始影像的信息。采用5種定量評價指標分析，本文提出的融合算法包含更多的信息，提高了融合影像的紋理細節，并較好地保持了原始影像的多光譜特征。