二維逐步正交匹配追蹤算法

邵 然，沈 軍

1.哈爾濱工程大學(xué) 信息與通信工程學(xué)院，哈爾濱150001

2.上海無線電設(shè)備研究所，上海200090

1 引言

2006年，Donoho、Candes等人[1-2]提出的壓縮感知理論，將信息的壓縮和采樣過程合為一體一并進行，開辟了計算成像新領(lǐng)域。作為最早也是最有名的應(yīng)用之一，萊斯大學(xué)的單像素相機[3]為壓縮感知成像指明了研究方向。但無論是何種系統(tǒng)，都會面臨圖像作為二維信號這一難關(guān)。直接的方法，比如在單像素相機中應(yīng)用DMD模塊，或者將輸入二維信號直接降維一行行掃描重構(gòu)成像，都面臨計算效率低、塊相應(yīng)明顯或者重構(gòu)圖像視覺效果差的問題[4]。而將采樣矩陣和稀疏表示對場景進行適配，比如Jiang等人[5]提出的基于玫瑰花瓣掃描的紅外成像系統(tǒng)，雖有較好的重構(gòu)效果，采樣矩陣卻受限于場景難以得到推廣。

在此基礎(chǔ)上，針對二維信號的采樣問題，Ghaffari等人[6]提出了二維觀測模型（2D Measurement Model，2DMM）。通過兩次采樣過程，使得任意類型的二維矩陣都可以變?yōu)樾辛邢∈璧姆叫尉仃嚕M而便利的壓縮感知處理。這種框架突破了原有一維信號壓縮感知的限制，在不進行各種變換的情況下保留了行間相關(guān)性信息，并進一步降低了采樣后數(shù)據(jù)量。在同一篇文章中，還提出了基于此框架的2D-SL0 重構(gòu)算法，該算法是基于一維信號的l0算法思想改進得到的。該文指出，根據(jù)傳統(tǒng)重構(gòu)算法設(shè)計思路得到二位重構(gòu)算法是可行的。隨即，F(xiàn)ang等人[7-8]提出了基于2DMM的OMP算法和2D-OMP算法，實現(xiàn)了較快速度下的高質(zhì)量重構(gòu)。在此基礎(chǔ)上，田文飚等提出了2DMM-SP 算法[9]，也就是二維觀測模型下的SP算法（為了對應(yīng)，下文稱其為2D-SP算法）。

本文將借鑒該設(shè)計思路，對現(xiàn)有算法進行分析與改進，提出二維逐步正交匹配追蹤重構(gòu)算法（2D Stagewise Orthogonal Matching Pursuit，2D-StOMP）。提出的算法能夠?qū)D像二維觀測得到的信號精準恢復(fù)，運算時長得以控制，信噪比和恢復(fù)效果則得到很大的提升。

2 二維圖像壓縮感知模型

2.1 壓縮感知理論與二維圖像應(yīng)用難題

經(jīng)典壓縮感知框架中，自然信號f 利用某種逆變換Ψ 來得到稀疏表示x，之后經(jīng)過測量矩陣Φ 采樣得到測量值y，即：

其中，感知矩陣A=ΦΨ 。若其滿足約束等距特性（Restricted Isometry Property，RIP 準則）[2，10-11]，則通過求解一個l1優(yōu)化問題便可以還原出原信號的稀疏表示，即稀疏表示的重構(gòu)信號：

通過巧妙的設(shè)計測量矩陣和變換矩陣，使應(yīng)用場景中信號的采樣率突破奈奎斯特定理的限制，是壓縮感知突出之處。但對于上述模型，若按照相關(guān)模型一貫的假設(shè)x ∈?N，y ∈?N，原二維信號為f ∈?M×N，即原信號為圖像這種二維信號，則Φ 與Ψ 則難以設(shè)計。

為了能夠?qū)嚎s感知應(yīng)用到圖像處理中，一般將圖像進行變換后，按行或者列進行采樣，轉(zhuǎn)換為?MN大小的一維信號進行分步處理。該方法存在很多問題：若有一行或者一列重構(gòu)失敗會導(dǎo)致黑色或白色條狀或塊狀紋路，影響圖像可讀性和視覺效果；圖像行間信息遭到破壞；由于一維信號進行小波變換等處理的要求，使得圖像行或列必須是2 的整數(shù)次冪，應(yīng)用前景被限制；且圖像行列間稀疏度往往有大幅度變化，簡單的循環(huán)無法妥善處理。雖然改進重構(gòu)算法，進行分塊處理，可以在一定程度上改善這種狀況，但在信息量被壓縮的情況下，行間信息仍很難得到保留，直到2DMM被提出。

2.2 2DMM

二維觀測模型（2DMM），由Ghaffari 等人[6]于2009年提出，不過類似的設(shè)計思想在關(guān)聯(lián)領(lǐng)域文獻里也早有論述[12]，比如，在雷達成像領(lǐng)域它有一個名字叫二維稀疏信號模型（2D Sparse Signal Model）[13]，其核心思想便是，經(jīng)過兩次變換和采樣來得到規(guī)整的采樣值，同時保留行間信息。本文將采用其中一例模型[8]。

首先，假設(shè)輸入二維信號f 在Ψ 變換域上是k 稀疏的(k ?N2)，其本身是一個N×N 大小的矩陣，稀疏表示出的矩陣大小不變的信號為x，即：

而對于相應(yīng)的測量值，則也需要經(jīng)過兩次測量：

感知矩陣A=ΦΨ 。進一步的，根據(jù)2DMM中的證明可以得到：

其中，字典（dictionary）矩陣D=A ?A，?表示兩個矩陣的克羅內(nèi)克積。此時各個矩陣大小如表1所示。

表1 矩陣大小

因此，該模型下的優(yōu)化問題則近似表示為求下式的l2范數(shù)最小值[7]：

下文中，通過分析現(xiàn)有重構(gòu)算法和改進型，可以得到上式的求解方法。

在2DMM下，對圖像的長寬等不再有限制，且采樣得到的數(shù)據(jù)量大幅度降低。但傳統(tǒng)的壓縮感知重構(gòu)算法不再適用于這一觀測模型，因此設(shè)計新的二維重構(gòu)算法便勢在必行。

3 改進的二維重構(gòu)算法

壓縮感知重構(gòu)算法始于凸優(yōu)化問題，而應(yīng)用最廣泛的便是以正交匹配追蹤算法（OMP）為代表的貪婪算法。在此基礎(chǔ)上，其每次選擇的原子（atoms）增加數(shù)量并進行正則化處理，便得到正則化正交匹配追蹤算法（ROMP）；而對每次選擇的原子的數(shù)量進行選取條件和復(fù)雜度計算后的調(diào)整，便可以得到相應(yīng)的子空間追蹤算法（SP）等；進一步的，將每一次循環(huán)得到的原子數(shù)量不規(guī)定為固定數(shù)，而是符合迭代閾值的若干數(shù)，為逐步正交匹配追蹤算法（StOMP）；最終，如果連稀疏度與每次原子數(shù)也不予規(guī)定，能夠直接應(yīng)用于各場景的，為稀疏度自適應(yīng)匹配追蹤（SAMP）。由此可見，每次迭代時的原子數(shù)量同閾值、算法效率都有密切關(guān)系[14-15]。因為每次迭代都應(yīng)該選擇與殘差相關(guān)性最佳的原子，單純的增加每次的原子數(shù)并不一定能夠增加效率，這種情況下殘差增加，既導(dǎo)致迭代次數(shù)增加，又會影響最終重構(gòu)精度。因此應(yīng)該選擇恰當閾值，使得每次迭代原子數(shù)量大致在稀疏度K 左右，此時其處于一個平衡點，使算法整體效率得到很大提升。

在二維壓縮感知重構(gòu)算法的設(shè)計上，同樣可以借鑒這種思路，增加2D-OMP 算法中每次循環(huán)的原子數(shù)量，在每次循環(huán)迭代求式（6）最小值的過程中，從構(gòu)造的字典每次選取符合閾值的、與殘差最相關(guān)的幾列原子，得到子集后一同進行處理，這便是2D-StOMP算法。該算法控制了復(fù)雜度，并提高了運算效率，其中，最關(guān)鍵的便是每次循環(huán)選取原子的閾值[15]。由文獻[16]可知，提出的StOMP 算法，其原子選取的閾值依據(jù)內(nèi)部噪聲指定。在一維信號處理過程中，經(jīng)分析其噪聲為。其中rs為重構(gòu)算法殘差，n 為重構(gòu)序列長度。在此之上再乘以2 到3 的增幅ts，便得出閾值。在二維重構(gòu)算法中，因為系統(tǒng)設(shè)計和觀測矩陣上高度相似性，每次循環(huán)重構(gòu)序列長度為M ，而增幅ts由經(jīng)驗[16]取值為2，因此新設(shè)計算法的閾值為。下面將仿照文獻[16]對此進行證明。

又因為：

故：

綜上：

其后，第S 次迭代按照公式（11）替代公式（7）中z，同理可證：

算法步驟如下。

（1）輸入：測量值y ∈?M×M；感知矩陣A ∈?M×N；稀疏度k ∈?+。

（2）初始化：輸入輸入量；殘差R=y；字典矩陣D=A ?A：

閾值增幅ts=2；原子計數(shù)δ=0。

主循環(huán) 循環(huán)數(shù)t 從1到k：

2.更新D=D?F；

4.選取C 中比Th 大的δ'個元素，δ'=δ+δ'；

5.更新( Iδ-δ'+1,Jδ-δ'+1)到( Iδ,Jδ)分別為上述元素的坐標；

6.更新標記矩陣F，上述( )

I,J 位置置零；

8.更新

4 理論分析與仿真

4.1 理論分析

2D-StOMP 算法的求解過程必定符合貪婪算法的內(nèi)在邏輯，對其算法進行的詳細分析既要包含有效性分析也要包含復(fù)雜度分析。在某種程度上，有效性影響算法重構(gòu)精度，復(fù)雜度影響算法運算時長。

4.1.1 有效性分析

當且僅當一個壓縮感知研究模型符合衡量其的三個定量指標，才可以說這個模型是嚴格壓縮感知意義上的完全可重構(gòu)模型。而這三個指標，即RIP準則[2，10]（RIC常數(shù)）、Spark 判別理論（Spark 常數(shù)）[17]和相關(guān)性判別理論（互相關(guān)系數(shù)）[18]，是有內(nèi)在聯(lián)系的。一般來講，一個模型中設(shè)計的觀測矩陣，只要符合其中一個指標的判定，就可以說這個觀測模型的設(shè)計是有效的，也可以直接調(diào)用相關(guān)的重構(gòu)算法進行信號恢復(fù)。本文將使用RIP準則來進行判定，將其作為引理1引入。

引理1（2k 階約束等距性）

當且僅當式（12）成立時，式（13）有解。其中，x ∈?N，為k 階稀疏信號；常數(shù)δ2k∈(0,1)。該公式證明參考文獻[2，10]。

引理2 如果A ∈?p×q，B ∈?r×s，則下式成立。

該公式證明參考文獻[19]定理3.7。引理2 提供了判定矩陣克羅內(nèi)克積的RIP參數(shù)計算方法。由引理1和引理2 可證明本文涉及到的2DMM 是有效的壓縮感知模型，如定理1所述。

定理1 對于二維矩陣信號x ∈?N×N，若其每行每列皆為k 稀疏，且感知矩陣A=ΦΨ 符合RIP 準則，則min‖ ‖?2有解。其中，字典矩陣D=A ?A。

4.1.2 復(fù)雜度分析

2D-StOMP 算法復(fù)雜度為每次循環(huán)運算量與循環(huán)次數(shù)的乘積，詳細分析如下。

（1）初始化：計算字典矩陣復(fù)雜度為M×N 矩陣的克羅內(nèi)克積運算的復(fù)雜度，即

（2）主循環(huán)：映射部分，由表1 中所列和文獻[8]中3.4.1 小節(jié)計算得，更新C 時運算量最大，為其余更新R 等運算量逐次增加，最大為

綜上所述，由于kδ ＜M2N2，所以主循環(huán)算法復(fù)雜度為。而采樣率不同，則δ 與M 的大小難以判斷，可將復(fù)雜度定為O(M2N2+kMN2) 。相比較而言，雖然2D-OMP文獻中復(fù)雜度為O(kMN2) ，但其未將初始化中克羅內(nèi)克積運算統(tǒng)計在內(nèi)。如果統(tǒng)一標準，兩種算法都移除克羅內(nèi)克積的運算過程，那么兩者復(fù)雜度相同，這也意味著，理論上兩種重構(gòu)算法的運行時長應(yīng)該相差無幾。另一方面，如果按照文獻[9]的標準，計算極限條件下（每次只選入一個原子）的算法復(fù)雜度，也不考慮克羅內(nèi)克積，則該算法復(fù)雜度比2D-SP 算法的要大。但極限條件下兩者實質(zhì)退化為類似于2D-OMP 的算法，并沒有實際意義；在非極限條件下兩者復(fù)雜度相近。

除此之外，新算法若與二維信號一維化后再進行OMP 重構(gòu)的傳統(tǒng)方法相比更大幅度地降低了復(fù)雜度。

4.2 仿真實驗

為了驗證新設(shè)計算法的優(yōu)越性，本文設(shè)計了仿真實驗。實驗采用二維隨機高斯矩陣為測量矩陣，二維DCT 變換矩陣為變換矩陣；重構(gòu)算法包括2D-StOMP、2D-SP、2D-OMP、2D-SL0 和1D-OMP；平臺電腦操作系統(tǒng)為64位Windows7，雙核2.10 GHz CPU，內(nèi)存6 GB。

重構(gòu)算法恢復(fù)質(zhì)量的評價標準有人眼視覺效果評估、峰值信噪比PSNR、運行時長等。實驗結(jié)果和分析如下。

4.2.1 視覺效果分析

實驗輸入圖像有四例：戰(zhàn)斗機的紅外圖像；圖像處理常用例子Lena 的灰度圖像；高光譜圖像測試集Pavia University Scene 的一幀部分高光譜圖像；實驗室自行采集的遠處客機的紅外圖像。圖片大小：前三者為256×256 像素，最后一張為128×128 像素。采用此四例能夠充分展示在不同圖像類型、背景、內(nèi)容下不同重構(gòu)算法恢復(fù)效果的細節(jié)。

視覺效果實驗采樣率為M/N=ε=0.55 。由于2DMM下實際上進行了橫縱兩次采樣，因而嚴格意義上的采樣率，即采樣信息量與輸入信息量的比值為ε2，即0.3。實驗添加了采樣率ε 為0.3情況下一維壓縮感知重構(gòu)算法作為對照組。實際上，由于壓縮感知理論限制，一維壓縮感知框架在采樣率0.3 以下就無法恢復(fù)出圖像，因此在0.3以下2DMM框架具有獨特的優(yōu)勢。

在采樣率相同，稀疏度同取64的情況下，各重構(gòu)算法恢復(fù)效果如圖1 所示。如圖所示，可見在低采樣率同等條件下，2D-StOMP的恢復(fù)效果和PSNR 優(yōu)于其他重構(gòu)算法：既沒有明顯的塊效應(yīng)，也沒有重構(gòu)失敗導(dǎo)致的黑色條紋，很大程度上保留了行間信息，降低了重構(gòu)中由于收斂不當引入的噪聲。特別的，與2D-SP 相比，雖然PSNR 大致接近，但在一些細節(jié)部分，諸如Lena 圖像中的帽子和臉部部分，前者視覺效果上要更加細膩（下文數(shù)據(jù)將會證明，隨著采樣率增加，兩者差距將進一步拉大）。

4.2.2 峰值信噪比與運行時長

仿真算法對同一圖像進行重構(gòu)，計算輸出圖像同輸入圖像的峰值信噪比（PSNR）和運行時長，對不同采樣率條件下，100次運算取得平均值。

MSE是指圖像均方根誤差。

由于采樣率在0.55（即一維壓縮感知采樣率為0.3）以下時，一維重構(gòu)算法很難恢復(fù)圖像，因此對比中未列入1D-OMP算法。運行結(jié)果如圖2所示。

由圖2 及表2 可知，2D-StOMP 算法在采樣率為0.3及以上時，PSNR要優(yōu)于其他算法。

圖1 視覺對比圖

圖2 重構(gòu)算法PSNR對比圖

圖3 重構(gòu)算法運行時長對比圖

表2 不同算法在稀疏度變化時恢復(fù)圖像所用時間和PSNR

由圖3 及表2 可知，2D-StOMP 算法運行時長要高于其他算法一到兩個數(shù)量級，但并未超過秒級。除非是某些對時間極其敏感的應(yīng)用場景，比如紅外制導(dǎo)、自動駕駛、武器控制等等，在其余場景（比如高光譜成像等），考慮到其相較其余二維重構(gòu)算法更為優(yōu)秀的恢復(fù)效果、相較一維重構(gòu)算法更加優(yōu)越的壓縮比，一定時間上的犧牲可以接受。并且，隨著電子元器件的進步和計算能力的不斷提高，上述算法運行時長間的差距會越發(fā)不明顯。

綜上所述，這里所設(shè)計的仿真實驗，有效地證明了2D-StOMP 算法相比較于其他算法的優(yōu)越性。這種優(yōu)勢體現(xiàn)在，同等采樣率下，該算法對同一圖像的重構(gòu)效果更佳，且重構(gòu)時長沒有過度增長；且經(jīng)過編程語言的更替與改進，重構(gòu)時長也不會有數(shù)量級上的差距。

歸結(jié)起來，這是二維觀測模型和靈活選取迭代原子數(shù)量疊加的雙重效果。前者保證其相對于一維壓縮感知的優(yōu)勢：即使經(jīng)過改進的稀疏矩陣，在壓縮感知步驟中的稀疏變換后對信號的采樣，亦會導(dǎo)致其丟失行間信息；相對應(yīng)的，由于二維觀測模型是針對二維矩陣信號，尤其是圖像而設(shè)計的，因而其更能適應(yīng)圖像行與行間、或者列與列間稀疏度大幅度變化的特殊情況，也能更好地適應(yīng)圖像壓縮感知應(yīng)用場景，與此同時，由于二維觀測模型中兩次變換的存在，觀測所需的信息量也被再次壓縮，促成更高效率的壓縮感知算法。后者則保證了其相對于現(xiàn)有二維重構(gòu)算法的優(yōu)勢：固然，對于2D-OMP等算法而言，單純的增加迭代次數(shù)及稀疏度能在一定程度上改善恢復(fù)效果；但迭代后的殘差卻沒有任何控制手段，由于殘差變化幅度大，選取原子整體上便容易“背離”原始數(shù)值，同時，由于同樣的原因，變化幅度大導(dǎo)致不易收斂到某一節(jié)點，平均計算單次迭代中計算量反而增加，直接導(dǎo)致算法運行時長大幅度增加（達到同等效果往往需要運算時長增加兩個數(shù)量級以上），這限制了通過增加迭代層次來改進算法的途徑，諸如2DSP 等固定單次迭代原子數(shù)量的算法也存在類似問題；而通過選取閾值靈活選擇原子，本文算法能夠控制每次迭代后殘差都保持在一個范圍內(nèi)，更為下一層迭代提供了更佳的閾值選取標準，因而理論上其在同等時長下理應(yīng)得到更高的恢復(fù)精度，且圖像行列間關(guān)聯(lián)性越高，稀疏度越大越有利于該算法，這在視覺仿真實驗中亦有所體現(xiàn)，人體面部細節(jié)上該算法明顯優(yōu)于其他算法。

結(jié)合前文分析，考慮到2D-StOMP 算法與2D-OMP等傳統(tǒng)算法在復(fù)雜度這一因素上相差不大，這意味著，若采用更加契合循環(huán)的程序語言對本文算法進行實現(xiàn)，其運算時長不會同其他算法產(chǎn)生數(shù)量級上的差距：眾所周知，Matlab 語言不善于處理多次循環(huán)，以及大小不確定數(shù)組的內(nèi)存空間分配。因此，相信經(jīng)過進一步的優(yōu)化和改進，屆時2D-StOMP算法的運行時長將會得到有效降低，算法的性能也能有所提高。

5 結(jié)束語

針對二維信號壓縮感知重構(gòu)算法存在的相關(guān)問題，本文結(jié)合一維重構(gòu)算法的設(shè)計思路，提出了2D-StOMP算法。提出的算法通過每次循環(huán)迭代時選擇契合閾值條件的原子來提高重構(gòu)準確性。理論分析和仿真實驗表明，本文算法解決了現(xiàn)有的塊效應(yīng)和部分行列重構(gòu)失敗問題，重構(gòu)效率得到很大提升，示例圖像的重構(gòu)質(zhì)量得到了明顯改善。