黃宏偉,謝正光,蔣小燕,蔡 旭
(南通大學(xué) 電子信息學(xué)院,江蘇 南通226019)
近幾年一種高效的信號(hào)采集技術(shù)壓縮感知 (compressed sensing,CS)[1]逐漸被人們學(xué)習(xí)了解。它的模式是邊采集邊壓縮,即在數(shù)據(jù)采集過程中就進(jìn)行適當(dāng)?shù)膲嚎s,這樣信號(hào)采集速度能低于奈奎斯特采樣速度,節(jié)省大量資源[2,3]。
壓縮感知測(cè)量過程可以用方程描述為y =Φx,其中y是長(zhǎng)度為M 的觀測(cè)信號(hào),Φ =[,…]是測(cè)量矩陣,大小為M×N ,x為K 階稀疏信號(hào),長(zhǎng)度為N,K <M <N,這個(gè)過程相當(dāng)于把原始信號(hào)投影到[,…]張成的空間上[4]。
壓縮感知中重構(gòu)算法是核心內(nèi)容,是從測(cè)量矩陣Φ和觀測(cè)信號(hào)y中恢復(fù)原始稀疏信號(hào)x。這個(gè)過程可以當(dāng)作l0最小化問題來求解:x=argmin x0subject to y=Φx,理論上,沒有噪聲的前提下,K 階稀疏信號(hào)僅需M ≥2 K 次測(cè)量就可以恢復(fù)。但是l0最小化問題是NP-h(huán)ard問題。Donoho和Candes指出不必求助于l0解稀疏解,而是通過一個(gè)更簡(jiǎn)單的l1最優(yōu)化方法:x =argmax x1subject to y =Φx,但是要求Φ滿足RIP[5]。雖然線性規(guī)劃方法重構(gòu)效果很好,但是它的高復(fù)雜度限制了其實(shí)際應(yīng)用,這就要求我們?cè)O(shè)計(jì)一種快速的求解方法。
貪婪算法因?yàn)槠漭^低的復(fù)雜度而受到越來越多人的關(guān)注。包括匹配追蹤 (matching pursuit,MP),OMP[6-8]、正則化正交匹配追蹤 (regularized OMP,ROMP)[9]、分段式正交匹配追蹤 (stagewise OMP,StOMP)、壓縮采樣匹配追蹤算法 (compressive sampling matching pursuit,Co-SaMP)[10]、SP[11]等等。這些方法的基本思想是通過相關(guān)性從不相關(guān)的原子中尋找正確的支撐集,每一個(gè)迭代過程都是從原子序列中找出認(rèn)為足夠相關(guān)的原子索引并添加到支撐集中。MP算法每次尋求的是局部的最優(yōu)解,雖然執(zhí)行速度快,但是最終不一定收斂到全局最優(yōu),所以重構(gòu)精度不是很高。帶有回溯思想的SP算法有兩個(gè)特征,一是當(dāng)信號(hào)非常稀疏時(shí),SP 算法比OMP 算法計(jì)算復(fù)雜度低;二是SP算法的重構(gòu)精確度與LP 方法相當(dāng)[11]。但是,SP 算法沒有考慮代理信號(hào)和殘差關(guān)系,每次迭代都選擇K 個(gè)原子,這可能導(dǎo)致迭代過程中選擇錯(cuò)誤支撐集數(shù)目增多。經(jīng)典的OMP算法有兩個(gè)缺點(diǎn):第一是每次選擇的最大相關(guān)原子不一定是正確原子;第二是每次迭代添加的支撐集原子索引都是不會(huì)被移除的,即無(wú)自主糾錯(cuò)能力。針對(duì)OMP這兩個(gè)缺點(diǎn),我們提出RIPTMP 算法的原子選擇方案。RIPTMP首先根據(jù)RIP和殘差條件一次選出多個(gè)原子,其次算法能夠根據(jù)信號(hào)特點(diǎn)剔除可能錯(cuò)誤的原子。在該方案下,算法提高了尋找支撐集的精度。
首先簡(jiǎn)單說明本文涉及的符號(hào)。supp(x)表示向量x 中非零元素的位置,即向量x 的支撐集。K 表示稀疏度。ΦT表示Φ的轉(zhuǎn)置,T0定義為初始的支撐集,表示初始的殘差,Tl表示第l次迭代得出的支撐集,表示Tl的補(bǔ)集。表示第l次迭代得出的信號(hào)殘差,表示x 對(duì)應(yīng)Tl位置上的元素集合。ΦI表示I 中所有索引對(duì)應(yīng)原子的集合,如果ΦI可逆,那么yp表示y 在ΦI上的投影yp=<y,,其中,代表ΦI的偽逆。δK+1為K+1階RIP常數(shù),表示第l次迭代的代理信號(hào)表示代理信號(hào)中的元素,其中<·>表示內(nèi)積運(yùn)算,ei為單位矩陣的第i列。
RIPTMP算法在原子選擇上帶有回溯思想,經(jīng)典的SP算法也具有同樣思想。SP算法需要待重構(gòu)稀疏信號(hào)的稀疏度K。算法在每次迭代過程中,原子添加的和刪減的數(shù)量都為K。具體算法見表1。
表1 子空間追蹤算法
SP有兩個(gè)缺點(diǎn),一方面在原子選擇數(shù)量上沒有考慮殘差變化,每次選取的數(shù)量為固定值K,這就可能導(dǎo)致當(dāng)前選擇多個(gè)錯(cuò)誤的原子,從而引起重構(gòu)精度下降;另一方面算法需要待重構(gòu)稀疏信號(hào)的稀疏度。
RIPTMP算法除了需要測(cè)量矩陣Φ和觀測(cè)信號(hào)y,還需要算法停止參數(shù)ε,支撐集添加閾值E,支撐集保持閾值P,和最大迭代次數(shù)imax。算法具體流程見表2。
表2 基于RIP閾值機(jī)制的匹配追蹤算法
RIPTMP算法中有兩個(gè)參數(shù),E 和P。下面我們給出這兩個(gè)參數(shù)選擇依據(jù):
1.3.1 預(yù)備知識(shí)
(1)有限等距性質(zhì) (RIP)定義[11]
如果矩陣Φ ∈RM×N滿足參數(shù)為 (K,δ)的有限等距性質(zhì),K ≤M ,0≤δ<1。對(duì)于任意的k 階稀疏信號(hào)x ∈
(2)最大相關(guān)選擇原理
假設(shè)Φ的所有列向量都是完全正交的,即Φ為正交矩陣,那么所求代理信號(hào)系數(shù)<,y>中非零值對(duì)應(yīng)的位置就是原信號(hào)的支撐集。但是大部分Φ中列向量都不是完全正交的,(1)說明滿足RIP條件的Φ近似為正交矩陣,這時(shí)候最大相關(guān)系數(shù)對(duì)應(yīng)的位置索引不一定在原信號(hào)支撐集中。但是,正確原子和觀測(cè)信號(hào)的相關(guān)系數(shù)和非正確原子和觀測(cè)信號(hào)的相關(guān)系數(shù)是有界的。因此我們可以根據(jù)相關(guān)系數(shù)的界尋找支撐集。
(3)引理
引理1[9]對(duì)于任意的兩個(gè)正整數(shù)m≤n,有δm≤δn。
引理2[5]設(shè)supp(x),supp(y){1,2,…N},|supp(x)|≤m,|supp(y)|≤n,supp(x)∩supp(y)=,那么|<Φx,
引理3 設(shè)向量x∈RN,|supp(x)|≤m,那么
1.3.2 E,P參數(shù)估計(jì)
第一次迭代:
假設(shè)矩陣Φ∈RM×N滿足參數(shù)為 (K +1,δ)的有限等距性質(zhì),根據(jù)引理2 有:isupp(x),|h1(i)|=|<
同理對(duì)于第l次迭代
結(jié)合第一次迭代和第l次迭代,我們把i∈supp(x)的|hl(i)|和jsupp(x)的|hl(j)|的界進(jìn)行估計(jì)。即對(duì)于任意jsupp(x),有
通過上面推導(dǎo)過程可以發(fā)現(xiàn)E 和P 的邊界條件經(jīng)過了放縮,E估計(jì)的界是在實(shí)際的界基礎(chǔ)上進(jìn)行了放大,P 估計(jì)的界是在實(shí)際界的基礎(chǔ)上進(jìn)行了縮小,因此在實(shí)驗(yàn)仿真階段可以在估計(jì)出來的界基礎(chǔ)上適當(dāng)?shù)姆趴s。
雖然RIPTMP算法是鑒于SP算法兩個(gè)缺點(diǎn)提出來的,但是它的理論基礎(chǔ)卻來源于OMP 算法。OMP 算法每次迭代只選擇最相關(guān)的原子,同時(shí)選擇的原子不會(huì)被移除,而實(shí)際上,最相關(guān)原子不一定是正確原子,當(dāng)算法引入了錯(cuò)誤原子時(shí)應(yīng)該通過回溯過程進(jìn)行剔除。因此,我們先把所有可能正確的原子都選出來,然后再剔除可能錯(cuò)誤的原子。通過理論分析,添加和刪減原子的數(shù)量和殘差,RIP 常數(shù)密切相關(guān),所以我們通過殘差和RIP 常數(shù)條件求出相應(yīng)閾值來篩選原子。
為了全面測(cè)試算法重構(gòu)性能,對(duì)一維信號(hào)和二維信號(hào)都進(jìn)行了相關(guān)實(shí)驗(yàn)。在一維信號(hào)部分,選擇了3種不同類型的稀疏信號(hào)分別為高斯信號(hào), “0-1”信號(hào)和均勻分布信號(hào)。二維信號(hào)則是標(biāo)準(zhǔn)的Lena,Boat測(cè)試圖像。
令測(cè)量矩陣為標(biāo)準(zhǔn)正態(tài)分布矩陣,測(cè)量次數(shù)M=128,信號(hào)長(zhǎng)度N=256。OMP 和RIPTMP 迭代停止閾值相同,均取ε=10-6,RIPTMP 的最大迭代次數(shù)imax=0.5 M,BP算法重構(gòu)通過cvx工具箱實(shí)現(xiàn)。FBP 算法給前長(zhǎng)步長(zhǎng)α=0.3 M,后退步長(zhǎng)β=α-1,停止閾值ε=10-6,最大迭代次數(shù)imax=0.5 M。實(shí)驗(yàn)中RIPTMP 算法有3種不同的設(shè)置,分別是(E=0.2,P=0.7)、(E=0.1,P=0.6)、(E=0.1,P=0.4)。重構(gòu)成功條件:均方誤差,實(shí)驗(yàn)次數(shù)times取500。
實(shí)驗(yàn)結(jié)果如圖1所示。
為了衡量算法對(duì)稀疏圖像信號(hào)重構(gòu)性能,選擇標(biāo)準(zhǔn)的256×256的Lena,Boat作為測(cè)試圖像。隨機(jī)產(chǎn)生大小為M×N 的測(cè)量矩陣Φ,Φ為歸一化的準(zhǔn)高斯矩陣,對(duì)圖像信號(hào)進(jìn)行小波變換選擇sym6小波基。保留小波變換后系數(shù)矩陣每列最大的q個(gè)小波系數(shù),其余置零。峰值信噪比計(jì)算
圖1 一維信號(hào)仿真實(shí)驗(yàn)
分別用SP,BP,RIPTMP,F(xiàn)BP這4種重構(gòu)算法重構(gòu)圖像信號(hào),其中SP稀疏度設(shè)置為q,RIPTMP參數(shù)設(shè)置為(E=0.1,P=0.6),F(xiàn)BP停止參數(shù)和最大迭代次數(shù)相同。重構(gòu)結(jié)果如圖2所示。
圖2 二維信號(hào)仿真實(shí)驗(yàn)
(1)一維信號(hào)重構(gòu)實(shí)驗(yàn)
在圖1 (a)高斯信號(hào)仿真實(shí)驗(yàn)中,E=0.2,P=0.7的RIPTMP算法的重構(gòu)誤差比OMP,SP 低,成功重構(gòu)概率明顯高于OMP,SP算法。SP成功重構(gòu)概率從K=40開始下降,而RIPTMP從K=43才開始下降。當(dāng)我們把前面推導(dǎo)并估計(jì)出的界E 適當(dāng)收縮到0.1,同時(shí)P 放大到0.6,RIPTMP效果明顯提升。一方面重構(gòu)均方誤差低于包括BP在內(nèi)的其它所有算法,當(dāng)K=50才由0逐步上升,而且在稀疏度取值區(qū)間[1,2,3...56]內(nèi)都保持較低的值;另一方面成功重構(gòu)概率在相同稀疏度條件下比其它算法都要高,當(dāng)稀疏度達(dá)到52 以上時(shí),才出現(xiàn)重構(gòu)失敗情況。參數(shù)為(E=0.1,P=0.4)的RIPTMP 重構(gòu)情況和參數(shù)為 (E=0.1,P=0.6)RIPTMP重構(gòu)情況相當(dāng)。
在圖1 (b)“0-1”信號(hào)實(shí)驗(yàn)中,OMP效果最差,從仿真數(shù)據(jù)可以看出 (為了方便顯示,曲線前面部分經(jīng)過適當(dāng)省略),當(dāng)稀疏度到10左右誤差就開始從0開始上升,而其它算法依然保持0的水平。當(dāng)K 增加到25左右SP算法和參數(shù)為 (E=0.2,P=0.7)的RIPTMP 和StOMP 重構(gòu)誤差開始從0增加,到K 增加到28左右,SP和FBP算法誤差開始從0上升。當(dāng)K 增加到31左右參數(shù)為 (E=0.1,P=0.4)和 (E=0.1,P=0.6)的RIPTMP算法重構(gòu)誤差開始從0上升。對(duì)應(yīng)的重構(gòu)百分,OMP從K=8處重構(gòu)開始出現(xiàn)失敗。參數(shù)為 (E=0.2,P=0.7)和StOMP 算法從K=28處重構(gòu)開始出現(xiàn)失敗。當(dāng)K 增加到29 時(shí),F(xiàn)BP和SP出現(xiàn)重構(gòu)失敗情況。當(dāng)K 達(dá)32左右,參數(shù)為 (E=0.1,P=0.4)和 (E=0.1,P=0.6)的RIPTMP 算法出現(xiàn)重構(gòu)失敗。BP算法重構(gòu)結(jié)果最好。
在圖1 (c)均勻信號(hào)實(shí)驗(yàn)中,RIPTMP 算法重構(gòu)誤差明顯小于OMP,接近BP 算法。參數(shù)為 (E=0.1,P=0.4)和 (E=0.1,P=0.6)的RIPTMP 算法在稀疏度[1,2,3…52]區(qū)間內(nèi)重構(gòu)的均方誤差和BP 算法相差不大,緊接著是SP算法,參數(shù)為 (E=0.2,P=0.7)RIPTMP算法等等。重構(gòu)百分比曲線上可以得出參數(shù)為 (E=0.1,P=0.6)和 (E=0.1,P=0.4)的RIPTMP 算法在稀疏度達(dá)到46 才開始下降,緊接著是FBP 算法、參數(shù)為(E=0.2,P=0.7)RIPTMP算法、SP算法等等。
實(shí)驗(yàn)說明,RIPTMP 算法在一定條件下,相比OMP,StOMP,SP,F(xiàn)BP算法有更高的重構(gòu)精度和重構(gòu)成功概率。雖然E,P參數(shù)推導(dǎo)過程得出相同的閾值條件,但是E=P≈0.4的閾值條件并不是準(zhǔn)確的。E推導(dǎo)的界比真實(shí)的界大,P推導(dǎo)的界比真實(shí)的界小,因此適當(dāng)改變E,P可以得出精度更高的支撐集估計(jì)效果,實(shí)驗(yàn)部分正好證明了這一點(diǎn)。
(2)二維信號(hào)重構(gòu)實(shí)驗(yàn)
圖2中,我們列出了SP,BP,F(xiàn)BP 和RIPTMP 算法重構(gòu)效果。Lena圖像重構(gòu)實(shí)驗(yàn)中,SP算法重構(gòu)信號(hào)的峰值信噪比為32.3dB,BP為33.4dB,F(xiàn)BP為33.2dB,參數(shù)為 (E=0.1,P=0.6)的RIPTMP 為33.4dB。Boat圖像重構(gòu)實(shí)驗(yàn)中,SP 算法重構(gòu)信號(hào)的峰值信噪比為31.7dB,BP為32.0dB,F(xiàn)BP 為31.8dB,參數(shù)為 (E=0.1,P=0.6)為32.0dB。可以看出RIPTMP算法的峰值信噪比高于SP算法,當(dāng)E=0.1,P=0.6時(shí),達(dá)到了BP算法水平。FBP 算法峰值信噪比也接近BP 算法,略低于RIPTMP算法。
二維圖像信號(hào)重構(gòu)中,雖然用小波變化對(duì)圖像進(jìn)行稀疏化處理,但是得到的并不是嚴(yán)格的稀疏矩陣。因此我們對(duì)小波變化后的矩陣進(jìn)行了處理,即保留小波變換后系數(shù)矩陣每列最大的q 個(gè)小波系數(shù),其余置零。因?yàn)闇y(cè)量矩陣的列向量不是完全正交的,所以可能導(dǎo)致數(shù)值較小的小波系數(shù)所在列的索引被當(dāng)成重要的原子來重構(gòu),而真正大的小波系數(shù)被忽略,這就會(huì)影響重構(gòu)的精度。在信號(hào)重構(gòu)部分,經(jīng)典的SP算法需要稀疏度條件,而圖像稀疏化過程中如果不進(jìn)行特定的處理不能獲取稀疏度的條件,因而不能確定每次迭代選取原子個(gè)數(shù)。RIPTMP 算法是通過最大迭代次數(shù)控制迭代過程,在稀疏度較小的時(shí)候,經(jīng)過很少的迭代就可以找出所有的原子,因此我們并不需要待重構(gòu)稀疏信號(hào)的稀疏度。所以,二維圖像信號(hào)重構(gòu)上RIPTMP 算法比SP算法更加合理。FBP算法雖然不需要稀疏度條件,但是每次更新原子的步長(zhǎng)是不變的,前長(zhǎng)步長(zhǎng)α=0.3 M ,后退步長(zhǎng)β=α-1。這種選擇原子同樣沒有考慮信號(hào)殘差的變化因素,且前向步長(zhǎng)的確定條件并沒有給出理論依據(jù),而RIPTMP 給出原子選擇的依據(jù),而且仿真效果好于FBP,因此更合理。
本文基于SP算法重構(gòu)精度低,需要已知稀疏度的缺點(diǎn)提出RIPTMP算法。首先算法在尋找支撐集上和SP不同,RIPTMP算法是根據(jù)RIP條件估計(jì)正確支撐集代理信號(hào)的幅值范圍,選出所有滿足的原子,而不是僅僅只選擇代理信號(hào)的最大值對(duì)應(yīng)的原子。其次,算法不需要稀疏度作為算法輸入條件。實(shí)驗(yàn)結(jié)果表明,在一定條件下RIPTMP 算法重構(gòu)精度高于OMP,SP,F(xiàn)BP 等經(jīng)典算法,達(dá)到了BP算法水平,而且算法不需要待重構(gòu)信號(hào)稀疏度。
[1]Candès E J,Wakin M B.An introduction to compressive sampling [J].Signal Processing Magazine,IEEE,2008,25(2):21-30.
[2]SHI Guangming,LIU Danhua,GAO Dahua,et al.Advances in theory and application of compressed sensing [J].Acta Electronica Sinica,2009,37 (5):1070-1081 (in Chinese). [石光明,劉丹華,高大化,等.壓縮感知理論及其研究進(jìn)展[J].電子學(xué)報(bào),2009,37 (5):1070-1081.]
[3]DAI Qionghai,F(xiàn)U Changjun,JI Xiangyang.Research on compressed sensing [J].Journal of Computers,2011,34(3):425-434 (in Chinese).[戴瓊海,付長(zhǎng)軍,季向陽(yáng).壓縮感知研究 [J].計(jì)算機(jī)學(xué)報(bào),2011,34 (3):425-434.]
[4]FANG Hong,YANG Hairong.Greedy algorithms and compressed sensing [J].Acta Automatica Sinica,2011,37 (12):1413-1421 (in Chinese).[方紅,楊海蓉.貪婪算法與壓縮感知理論 [J].自動(dòng)化學(xué)報(bào),2011,37 (12):1413-1421.]
[5]Candes E J.The restricted isometry property and its implications for compressed sensing [J].Comptes Rendus Mathematique,2008,346 (9):589-592.
[6]Wu R,Huang W,Chen D R.The exact support recovery of sparse signals with noise via orthogonal matching pursuit[J].Signal Processing Letters,IEEE,2013,20 (4):403-406.
[7]SHEN Y,LI S.Sparse signals recovery from noisy measurements by orthogonal matching pursuit [J].arXiv Preprint arXiv:1105.6177,2011.
[8]Ding J,Chen L,Gu Y.Perturbation analysis of orthogonal matching pursuit [J].IEEE Transactions on Signal Processing,2013,61 (2):398-410.
[9]Needell D,Vershynin R.Uniform uncertainty principle and signal recovery via regularized orthogonal matching pursuit[J].Foundations of Computational Mathematics,2009,9 (3):317-334.
[10]Needell D,Tropp J A.CoSaMP:Iterative signal recovery from incomplete and inaccurate samples [J].Applied and Computational Harmonic Analysis,2009,26 (3):301-321.
[11]Dai W,Milenkovic O.Subspace pursuit for compressive sensing signal reconstruction [J].IEEE Transactions on Information Theory,2009,55 (5):2230-2249.
[12]MO Q,SHEN Y.A remark on the restricted isometry property in orthogonal matching pursuit[J].IEEE Transactions on Information Theory,2012,58 (6):3654-3656.