基于冗余緊框架的?2/?1極小化塊稀疏壓縮感知

張楓,王建軍

(西南大學數學與統計學院,重慶 400715)

1 引言

壓縮感知[1-3](compressed sensing)作為一種新的采樣理論,它利用信號的稀疏特性,在遠小于Nyquist采樣率的條件下,用隨機采樣獲取信號的離散樣本,通過非線性重建算法完美地恢復信號.壓縮感知理論一經提出,就引起了學者廣泛關注.目前已在壓縮成像[4],醫學成像[5],模式識別[6],圖像處理[7]等領域得到了廣泛應用.

在壓縮感知中,主要考慮以下模型:y=Ax+w,其中,A∈Rm×n是測量矩陣,y∈Rm是已知的線性測量,x∈Rn是待重構的未知信號,w∈Rn是噪聲(‖w‖2≤ε).壓縮感知的核心思想是依賴于信號是否是稀疏的或者近似稀疏的,即信號x的非零元素的個數是否遠小于x的長度.然而,在現實中常見的自然信號不一定都具有稀疏特性,甚至這類信號在某些正交基上都不能夠進行稀疏表示.自然地,上述信號重構過程不能直接應用于自然信號的重構.研究表明,一些自然信號在某些緊框架D∈Rn×N(n≤N,DD?=In)[8-9]上能夠稀疏表示,即x=Dα,其中α ∈RN是(近似)稀疏的.從數學的角度講,上述問題可以寫成以下形式

其中D?表示矩陣D的共軛轉置,‖D?x‖0表示向量D?x中非零元素的個數,若‖D?x‖0≤k,那么就稱向量D?x為k-稀疏信號.由于問題(1)是一個 NP-hard問題,即在多項式時間內,計算機無法有效求解.所以一種更為實際的和易于處理的凸優化方法被提出

定義 1.1對于任意k-稀疏信號v(‖v‖0≤k),若存在 0<δk<1,使得

那么稱矩陣A滿足k階D-限制性等容性質(D-RIP),最小的δk稱為D-限制性等容常數(D-RIC).

進一步,他們指出當測量矩陣A滿足D-限制性等容性質且δ2k<0.08時通過求解有約束優化問題(2)可以實現信號的魯棒重構.之后,文獻[11]將上述條件作了進一步改善,得到δ2k<0.4931.已知問題(2)可以轉換為以下無約束優化問題

文獻[8]提出,若A滿足D-限制性等容性質(D-RIP),D是一個緊框架,當D-RIP常數δ2k<0.1907并且時,問題(3)的解滿足

其中‖D?x?(D?x)k‖1是最佳k-項?1逼近誤差.C1,C2是兩個數值型常數并且C1由 D-RIP 常數和‖D?D‖1,1:=sup{‖D?Dv‖1:v∈RN,‖v‖1≤1} 共同量化,C2僅取決于D-RIP常數.

然而現實世界中的自然信號其結構千變萬化,一種常見的結構方式是自然信號在某些冗余字典下是塊稀疏的,即其非零元素以塊的形式出現,例如彩色圖像處理[7]和 DNA-微列陣[12]等.從數學的角度看,給定分塊τ={τ1,τ2,τ3,···,τd},對于任意向量α∈RN都可以被描述為

其中α[i]表示向量α的第i個子塊,而αi則表示向量α的第i個分量元素.如果向量α最多有k個非零塊,即‖α‖2,0≤k,則稱向量α為塊k-稀疏信號.特別地,當d=1時,等價于傳統的帶字典的壓縮感知問題.相應地,測量矩陣A∈Rm×n和冗余字典D∈Rn×N也可以分別被描述為

其中A[i],D[i]分別表示矩陣A和D的第i個子塊(矩陣),而Ai,Di則分別表示矩陣A和D的第i個列向量.

然而使用原始的?1極小化方法來恢復此類塊稀疏信號不能充分利用信號的結構性特征,即非零元素是以塊的形式出現的這一特性.為此,一些學者對傳統的壓縮感知方法進行了針對性的改進.文獻[13]提出并研究了如下的?2/?1問題

定義 1.2對于任意塊k-稀疏信號v(‖v‖2,0≤k),若存在 0<δk|τ<1,使得

那么稱矩陣A滿足k階 Block D-限制性等容性質 (Block D-RIP),最小的δk|τ稱為Block D-限制性等容常數(Block D-RIC).

文獻 [14]指出當測量矩陣A滿足 Block D-限制性等容性質 (Block D-RIP)且δ2k|τ<0.4142時,塊稀疏信號能夠通過求解有約束優化問題(4)進行魯棒重構.其后,文獻[15]將上述條件作了進一步改善,得到了δ2k|τ<0.4931.類似地,問題(4)可以轉化為以下塊無約束優化問題

假設A滿足Block D-限制性等容性質(Block D-RIP),D是一個緊框架,本文研究表明,當 Block D-RIP 常數δ2k|τ<0.2 并且時,問題(5)的解滿足

其中C1,C2是兩個數值型常數并且C1由D-RIP常數和

共同量化,C2僅取決于D-RIP常數.

2 預備知識

為了方便介紹后文,首先給出以下記號.

?給定正整數d,記索引集T?{1,2,···,d}且|T|=k,Tc表示T在{1,2,···,d}中的補集.

?DT∈Rn×|T|表示從D中取出索引集T對應的列所組成的矩陣,記

?記T1為包含的k個最大2范數塊的索引集,T2為包含的k個最大2范數塊的索引集,等.

接下來,為了證明主要定理,需要以下引理.

引理 2.1

證明由Tj的構造,有

上式兩邊取2范數,得

所以,有

引理2.2令測量y=Ax+w,h=?x,若,則問題(5)的解滿足

證明由實值凸的低階半連續函數的次微分的定義和是問題(5)的解可知,滿足問題(5)的子梯度最優化條件,即

其中v∈RN,若,則,否則‖vi‖2≤1.因此存在v∈RN使得‖v‖2,∞≤1,進一步,有

引理得證.

引理 2.3令D∈Rn×N為滿足DD?=In的矩陣,A∈Rm×n為滿足Block DRIP條件的矩陣.令索引集T?{1,2,···,d}恰有k個塊,測量

證明令T0為包含D?h的k個最大2范數塊的索引集,因為,所以取T=T0就能充分證明該引理.

因為A滿足 Block D-RIP條件,不失一般性,假設存在

使得‖u‖2=‖c‖2=1,因此有

注意到D是一個緊框架,即DD?=In,故

因此

結合(6)式,(9)式-(10)式,有

引理得證.

引理2.4令D∈Rn×N為滿足DD?=In的矩陣,測量y=Ax+w,h=，若,則問題(5)的解滿足

將h=?x,y=Ax+w代入上式,由D是一個緊框架可知

由三角不等式,易知

整理后得(11)式,引理得證.

3 主要結果

定理 3.1令D∈Rn×N為滿足DD?=In的矩陣,A∈Rm×n為滿足Block DRIP(0<δ2k|τ<0.2)條件的矩陣,(D?x)k為由D?x的k個最大 2范數塊組成的向量,測量y=Ax+w,h=?x,若,則問題(5)的解滿足

證明由引理2.1和引理2.4,有

由于 1?3β2>0(0<δ2k|τ<0.2),故

由 (8)式和 (14)式,有

由(13)式-(15)式知,(12)式成立.

4 數值實驗

為了驗證理論結果,本文分別做了兩組實驗:(a)設計滿足定理3.1條件的算法;(b)理論誤差上界對比實驗.實驗在CPU為Inter(R)Core(TM)i3,內存為2GB的臺式電腦上進行,運行軟件為MATLAB(R2014a).實驗中,測量矩陣A∈R128×256服從標準高斯分布,字典D∈R256×1024由傅里葉變換矩陣和單位矩陣合并而成,并且滿足DD?=In,測量誤差w∈R128服從正態分布,取正則化參數λ=1e?3,從而滿足,待重構的信號為x∈R256,在字典D下的塊稀疏信號α∈R1024中的非零塊位置隨機產生.為了克服實驗結果的偶然性,所有實驗將獨立重復地進行100次.在本文中,將冗余字典和Block-IRLS算法[16-18]相結合提出D-Block-IRLS算法:輸入:分塊τ={τ1,τ2,···,τd},測量矩陣A,字典D,觀測信號y,塊稀疏度估計k.

輸出:重構信號x.

(a)選擇適當的懲罰參數λ(0<λ<1).

(b)初始化迭代向量α(0),使其滿足ADα(0)=y.設置?0=1.

(c)開始迭代j=0.

(d)通過α(j)解決下面的線性問題

(e)當α(j)滿足停機條件,將Dα(j)作為輸出賦值給x,同時結束算法,否則執行下一步.

其中,r(α)表示將向量α的分塊取?2范數后,再由大到小依次排列形成的向量.而r(α)k+1表示向量r(α)的第k+1個分量值.

(g)j=j+1,并返回到第4步繼續執行.

首先針對D-Block-IRLS算法,采用兩種分塊形式,均分256塊和非均分256塊.圖1分別研究了均勻分塊和非均勻分塊的情況下通過D-Block-IRLS算法得到的誤差‖?x‖2以及理論誤差上界與塊稀疏度k的關系.由圖可知,無論是均勻分塊還是非均勻分塊通過D-Block-IRLS算法得到的誤差都遠小于理論誤差上界,換言之,利用D-Block-IRLS算法來重構塊稀疏信號可以滿足實驗設計的需要并且從側面印證了本文理論分析的正確性.

圖1 D-Block-IRLS算法誤差與理論誤差上界對比

由文獻 [8]知,當D∈Rn×N為滿足DD?=In的矩陣,A∈Rm×n為滿足 DRIP(0<δ2k<0.1907)條件的矩陣,(D?x)k為由D?x的k個最大元素組成的向量,測量

現在分別用塊和非塊的方式來處理塊稀疏信號α,取

其他參數保持一致[9].圖2表明,無論是均勻分塊還是非均勻分塊,所獲理論誤差上界都明顯優于(16)式的誤差上限.

5 小結

本文采用?2/?1極小化方法研究了基于冗余緊框架下的塊稀疏信號的恢復,獲得了該方法魯棒重構原始信號的充分條件和誤差上界估計,所獲結果表明,誤差上界可以通過正則化參數λ,k-項逼近和塊稀疏度來控制.仿真實驗證明了理論結果的準確性,該結果對于推動壓縮感知的進一步發展具有一定的理論價值和借鑒作用.