壓縮感知中l(wèi)1問題的自適應(yīng)罰參數(shù)交替方向法

徐 麗 君,胡 瑞 戈,李 婷

(大連海事大學(xué) 理學(xué)院， 遼寧 大連 116026 )

0 引 言

Candes等[1]以及Donoho[2]提出的壓縮感知作為一個(gè)新的采樣理論，通過開發(fā)信號的稀疏特性，在遠(yuǎn)小于奈奎斯特采樣率的條件下，用隨機(jī)采樣獲取信號的離散樣本，通過重建算法完美地重建信號．壓縮感知理論一經(jīng)提出，就引起學(xué)術(shù)界和工業(yè)界的廣泛關(guān)注．它在信息論、圖像處理、地球科學(xué)、光學(xué)微波成像、模式識別、無線通信、生物醫(yī)學(xué)工程等領(lǐng)域受到高度關(guān)注．

假設(shè)x0∈Rn是希望得到的原始稀疏信號，不直接對x0進(jìn)行采樣，而是先得到一組線性觀測數(shù)據(jù)b=Ax0∈Rm，其中A∈Rm×n是編碼矩陣，m

這是一個(gè)NP難的組合優(yōu)化問題．因此，更好更可行的方案是求解凸優(yōu)化l1問題(基追蹤問題)：

因?yàn)樯鲜鰞蓚€(gè)問題在某些有利條件下具有相同的解[3]．

當(dāng)b包含噪聲或當(dāng)x0不完全稀疏而只是可壓縮時(shí)，對基追蹤問題進(jìn)行一定的松弛，使得基追蹤問題轉(zhuǎn)化為約束基追蹤去噪問題：

或非約束基追蹤去噪問題：

其中δ>0，μ>0是參數(shù)．從優(yōu)化理論角度出發(fā)，約束基追蹤去噪問題和非約束基追蹤去噪問題是等價(jià)的．而本文只對非約束基追蹤去噪問題進(jìn)行研究．

過去學(xué)者們構(gòu)造了一些迭代算法來求解上述稀疏優(yōu)化問題，如牛頓內(nèi)點(diǎn)法[4]、高效近塊坐標(biāo)同倫算法[5]、貪婪算法[6-7]、基于迭代收縮軟閾值理論的固定點(diǎn)連續(xù)算法(FPC)[8]、梯度投影算法[9-10]、快速迭代收縮軟閾值算法(FISTA)[11]、稀疏重構(gòu)算法(SpaRSA)[12]、基于梯度理論的分塊坐標(biāo)梯度下降法(CGD)[13]、光譜投影梯度法(SPGL1)[14]、Bregman迭代算法[15]等．Yang等[16]將交替方向法(alternative direction method，ADM)引到壓縮感知領(lǐng)域中解決l1問題，但在求解過程中，并沒有對增廣拉格朗日函數(shù)中的罰參數(shù)進(jìn)行深入研究．盡管理論上只要罰參數(shù)大于零就可保證ADM的全局收斂性，但在實(shí)際計(jì)算過程中罰參數(shù)的取值也會大大影響收斂速度．

針對上述問題，本文考慮根據(jù)迭代過程中某些量的數(shù)值關(guān)系，探討罰參數(shù)的動態(tài)調(diào)整準(zhǔn)則，設(shè)計(jì)一種自適應(yīng)的罰參數(shù)調(diào)整方法，以便能夠快速且穩(wěn)定地獲得最優(yōu)解．

1 l1問題的ADM

1.1 經(jīng)典的ADM

設(shè)函數(shù)f(x)：Rm→R和g(y)：Rn→R都是凸函數(shù)且相互獨(dú)立，x∈Rm，y∈Rn是優(yōu)化變量，等式約束中A∈Rp×m，B∈Rp×n，b∈Rp，凸優(yōu)化問題為

上述問題的增廣拉格朗日函數(shù)為

L(x，y，λ)=f(x)+g(y)-λT(Ax+By-b)+

其中λ是拉格朗日乘子，β>0是罰參數(shù)．為了獲得最優(yōu)解，通過增廣拉格朗日函數(shù)進(jìn)行迭代，其中第k步迭代如下：

對于上述迭代，需要同時(shí)計(jì)算x和y，這是比較困難的．而ADM把復(fù)雜的聯(lián)合最小化問題轉(zhuǎn)變?yōu)閮蓚€(gè)簡單子問題：第1步，已知(yk，λk)極小化函數(shù)L(x，y，λ)，求得xk+1；第2步，已知(xk+1，λk)極小化函數(shù)L(x，y，λ)，求得yk+1；第3步，已知(xk+1，yk+1)更新拉格朗日乘子λ．具體迭代如下：

1.2 應(yīng)用ADM于原始問題

基于ADM的思想和框架，引入輔助變量r，使非約束基追蹤去噪問題轉(zhuǎn)化為變量可分離問題：

于是，上式的增廣拉格朗日函數(shù)為

(1)

其中y∈Rm是拉格朗日乘子．對問題(1)使用交替極小化，具體步驟如下：

步驟1給定x=xk和y=yk，求問題(1)關(guān)于r的最小值rk+1為

rk+1=(μβ/(1+μβ))(yk/β-(Axk-b))

(2)

步驟2給定y=yk和r=rk+1，求問題(1)關(guān)于x的極小化問題，易知該問題可等價(jià)為下面的問題：

通過下面近似求解來替代對該問題精確求解：

其中τ>0是鄰近參數(shù)，且

gkAT(Axk+rk+1-b-yk/β)

最后得出解

xk+1=Shrink(xk-τgk，τ/β)=

max{|xk-τgk|-τ/β，0}°

((xk-τgk)/|xk-τgk|)

(3)

其中式(3)用的是軟閾值算法，“°”代表分量乘法．

步驟3給定x=xk+1和r=rk+1，更新乘子y：

yk+1=yk-γβ(Axk+1+rk+1-b)

(4)

其中γ>0是步長常數(shù)．

由以上推導(dǎo)，對非約束基追蹤去噪問題應(yīng)用ADM產(chǎn)生的迭代過程如下：

rk+1=(μβ/(1+μβ))(yk/β-(Axk-b))

xk+1=Shrink(xk-τgk，τ/β)

yk+1=yk-γβ(Axk+1+rk+1-b)

(5)

上述迭代過程實(shí)際上是一個(gè)不精確的ADM，因?yàn)樵谇蠼鈞k+1時(shí)用近似求解來代替精確求解．文獻(xiàn)[17]中只研究了γ=1時(shí)的收斂性，而文獻(xiàn)[16]對于γ>1的情況也給出了收斂結(jié)果并做出了證明．

1.3 應(yīng)用ADM于對偶問題

考慮非約束基追蹤去噪問題的對偶問題：

寫出增廣拉格朗日子問題：

其中x∈Rn是一個(gè)乘子(實(shí)際上是原變量)．這里假設(shè)矩陣A滿足AAT=I，對對偶問題用ADM，迭代如下：

(6)

值得注意的是，當(dāng)且僅當(dāng)r=μy時(shí)非約束基追蹤去噪問題的對偶問題中等式成立．因此，原始-對偶?xì)埐詈蛯ε奸g隙如下：

在計(jì)算中，迭代式(6)被下式終止：

其中ε>0．

本文主要針對對偶問題的ADM進(jìn)行改進(jìn)，即迭代式(6)．文獻(xiàn)[16]中已經(jīng)給出：對偶問題的ADM比原始問題的ADM更有效．

2 罰參數(shù)自適應(yīng)調(diào)整

罰參數(shù)β是ADM中的一個(gè)重要參數(shù)，雖然理論上只要β>0就可保證ADM的全局收斂性，但實(shí)際計(jì)算中β過大或過小都會影響收斂速度．以本文研究的對偶問題的ADM為例，即迭代式(6)，目前在罰參數(shù)β的選取上是固定的，即在迭代過程中不改變罰參數(shù)的取值，除非過大或過小．然而，在許多問題中，固定的β在迭代過程中，并非都是有效的．換言之，對目標(biāo)函數(shù)的尋優(yōu)和約束條件的懲罰不起作用．此時(shí)，依然使用該罰參數(shù)就會影響效率．

針對上述問題，本章首先分析罰參數(shù)在迭代過程中對目標(biāo)和約束的影響，然后根據(jù)其相對變化量，提出一種動態(tài)調(diào)整罰參數(shù)的策略，使得ADM在迭代過程中自動選取合適的β，使本文方法更有效．

2.1 罰參數(shù)分析

首先，將增廣拉格朗日子問題中給定的常值β改寫為與迭代步k相關(guān)的βk：

由于βk在目標(biāo)函數(shù)中相當(dāng)于對懲罰項(xiàng)(即對偶問題約束條件)的一種權(quán)重，其大小直接影響對偶問題目標(biāo)函數(shù)下降和約束條件的滿足情況．因此基于迭代過程中對偶問題目標(biāo)函數(shù)值的相對變化量與約束函數(shù)值的相對變化量，判斷βk是否應(yīng)該調(diào)整以及該如何調(diào)整．

(1)對偶問題目標(biāo)函數(shù)值的平均相對變化量：

(2)對偶問題約束函數(shù)值的平均相對變化量：

其中xk、rk表示經(jīng)過k次迭代后的值，xk+j、rk+j表示經(jīng)過k+j次迭代后的值．為了比較這兩個(gè)值的大小，采用計(jì)算兩個(gè)平均相對變化量之比

Δ=Δ1/Δ2

來進(jìn)行量化．

若βk取值恰當(dāng)，Δ1和Δ2的值是相近的，從而Δ會在閾值內(nèi)(Δ∈[δ1，δ2])，其中δ1<δ2是事先取定的值．若Δ不在此區(qū)間范圍，則說明當(dāng)前的罰參數(shù)βk在一定程度上不合適，需要進(jìn)一步分析調(diào)整．

顯然，當(dāng)βk過大時(shí)，Δ1<Δ2，即對偶問題目標(biāo)函數(shù)值下降速率相對較慢，通常會出現(xiàn)Δ<δ1；當(dāng)βk過小時(shí)，Δ1>Δ2，即對偶問題約束函數(shù)值下降速率相對較慢，通常會出現(xiàn)Δ>δ2．

根據(jù)上述分析，給出自適應(yīng)罰參數(shù)調(diào)整準(zhǔn)則：

當(dāng)Δ∈[δ1，δ2]時(shí)，βk+1=βk；

當(dāng)Δ<δ1時(shí)，調(diào)小罰參數(shù)：βk+1=aβk；

當(dāng)Δ>δ2時(shí)，調(diào)大罰參數(shù)：βk+1=cβk．

其中a<1，c>1是常數(shù)，稱其為調(diào)節(jié)參數(shù)．

2.2 自適應(yīng)罰參數(shù)ADM

將上述罰參數(shù)調(diào)整準(zhǔn)則應(yīng)用到原始的ADM，給出本文的自適應(yīng)罰參數(shù)ADM：

(1)輸入原始信號xs，觀測矩陣A，一維觀測量b，最大迭代次數(shù)，ε>0，調(diào)整頻率j>0，Δ閾值δ1和δ2，調(diào)節(jié)參數(shù)a<1，c>1．

(2)當(dāng)k小于等于最大迭代次數(shù)時(shí)，通過迭代式(6)繼續(xù)迭代．

(4)當(dāng)?shù)絢等于j的整數(shù)倍時(shí)，根據(jù)自適應(yīng)罰參數(shù)調(diào)整準(zhǔn)則更新罰參數(shù)；即每經(jīng)過j步，則更新一次罰參數(shù)．以更新的罰參數(shù)βk+1返回步驟(2)繼續(xù)迭代．

3 數(shù)值實(shí)驗(yàn)

本章通過實(shí)驗(yàn)分別對不同參數(shù)的影響做出分析，總結(jié)出一個(gè)快速魯棒的參數(shù)選取方案；再根據(jù)得出的方案，做出比較實(shí)驗(yàn)，通過設(shè)置不同的罰參數(shù)初始值觀察本文方法的魯棒性和有效性．

3.1 通過實(shí)驗(yàn)獲得參數(shù)方案

分別針對調(diào)整頻率j，Δ閾值，調(diào)節(jié)參數(shù)a、c進(jìn)行實(shí)驗(yàn)，研究不同的參數(shù)取值對自適應(yīng)罰參數(shù)調(diào)整準(zhǔn)則的影響．

由于罰參數(shù)過大或過小都會影響信號的恢復(fù)質(zhì)量，因此將罰參數(shù)設(shè)置為βk∈[5×10-5，50]．此外，取ε=5×10-4，最大迭代次數(shù)為9 999．

下面通過3組實(shí)驗(yàn)詳細(xì)研究準(zhǔn)則中的不同參數(shù)對本文方法效率的影響．

3.1.1 調(diào)整頻率對本文方法的影響 自適應(yīng)罰參數(shù)調(diào)整準(zhǔn)則中的頻率j表示根據(jù)每j步的信息計(jì)算準(zhǔn)則中的各個(gè)值．本次實(shí)驗(yàn)分別取j=5，10，20，30，實(shí)驗(yàn)結(jié)果見表1．

表1 不同調(diào)整頻率下的恢復(fù)效果指標(biāo)Tab.1 Recovery effect indicators with different adjustment frequencies

從表1可以看出，罰參數(shù)自適應(yīng)調(diào)整頻率的改變主要影響本文方法的運(yùn)行速度，隨著調(diào)整頻率的增大，達(dá)到終止條件需要的迭代次數(shù)變多，從而運(yùn)行時(shí)間變長．從表中很容易看出j=5時(shí)，本文方法所需的迭代次數(shù)和迭代時(shí)間最少，但是從所得解的質(zhì)量情況來看，并不是最好的．因此，綜合考慮，選取j=10為自適應(yīng)罰參數(shù)調(diào)整準(zhǔn)則中j的默認(rèn)取值．

3.1.2Δ閾值對本文方法的影響Δ閾值決定達(dá)到什么樣的條件才進(jìn)行罰參數(shù)的調(diào)整，間接影響罰參數(shù)的調(diào)整頻率，因此也被作為重要參數(shù)之一．對很多組不同的Δ閾值進(jìn)行了實(shí)驗(yàn)，本文僅列舉其中的4組實(shí)驗(yàn)，結(jié)果見表2．

表2 不同閾值下的恢復(fù)效果指標(biāo)Tab.2 Recovery effect indicators with different thresholds

根據(jù)本節(jié)實(shí)驗(yàn)結(jié)果，從所需的迭代次數(shù)和解的質(zhì)量等方面綜合考慮，選取Δ閾值[1/15，1/10]為自適應(yīng)罰參數(shù)調(diào)整準(zhǔn)則的缺省值．

3.1.3 調(diào)節(jié)參數(shù)對本文方法的影響 除了調(diào)整頻率、Δ閾值外，調(diào)節(jié)參數(shù)a<1，c>1決定了每次罰參數(shù)調(diào)整的大小，也是影響罰參數(shù)自適應(yīng)調(diào)整的一個(gè)重要因素．因此，在本實(shí)驗(yàn)中，設(shè)置不同的調(diào)節(jié)參數(shù)a、c，分析其對信號恢復(fù)效果的影響．部分實(shí)驗(yàn)結(jié)果見表3．

表3 不同調(diào)節(jié)參數(shù)下的恢復(fù)效果指標(biāo)Tab.3 Recovery effect indicators with different adjustment parameters

根據(jù)本節(jié)實(shí)驗(yàn)結(jié)果，可以看出本文方法對調(diào)節(jié)參數(shù)的大小不十分敏感．這也意味著本文所提出的自適應(yīng)罰參數(shù)調(diào)整準(zhǔn)則符合罰參數(shù)在本文方法迭代中的意義．選擇數(shù)值效果較好的a=0.9，c=2.00作為本文方法的默認(rèn)取值．

綜上所述，自適應(yīng)罰參數(shù)ADM選取的默認(rèn)參數(shù)方案：調(diào)節(jié)頻率j=10，Δ閾值為[1/15，1/10]，a=0.9，c=2.00．需要指出的是，本文方法的參數(shù)取值是經(jīng)過大量數(shù)值實(shí)驗(yàn)(本文只列舉了一組)得到的，一般不需要再次調(diào)試．

3.2 驗(yàn)證本文方法的魯棒性和有效性

本文方法是在程序包YALL1的基礎(chǔ)上改進(jìn)的，YALL1是基于l1對偶問題的ADM框架(即迭代式(6))編寫的Matlab程序，具體可參考文獻(xiàn)[16]．為了驗(yàn)證本文方法的魯棒性和有效性，在本例中設(shè)定不同量級的初始罰參數(shù)取值，對YALL1算法和本文方法所得的恢復(fù)效果指標(biāo)進(jìn)行對比，見表4．

表4 不同罰參數(shù)初始值下的恢復(fù)效果指標(biāo)Tab.4 Recovery effect indicators with different penalty parameter initial values

表4表明：(1)對于YALL1算法，當(dāng)罰參數(shù)初始值很小時(shí)，原始問題約束條件得到很好的滿足，但對偶問題的約束并未滿足(約束函數(shù)值達(dá)到102)，因此其恢復(fù)效果(e值)較差．當(dāng)罰參數(shù)初始值很大時(shí)，由于YALL1算法中對于罰參數(shù)較大的情形給出了一種減小罰參數(shù)的調(diào)整方案，因此，可以看到其恢復(fù)效果較好(相對誤差達(dá)到10-3)，但原始問題和對偶問題的約束條件滿足的精度仍然較低．(2)對于本文方法，無論初始罰參數(shù)取值過大或過小，其恢復(fù)效果都能夠達(dá)到10-3的誤差，并且原始問題和對偶問題的約束條件滿足的精度都較高，對偶間隙相對差都能達(dá)到10-4．(3)與YALL1算法較好的結(jié)果相比(第4、5組)，由于自適應(yīng)罰參數(shù)的動態(tài)調(diào)整，不僅所得解的性態(tài)更好，還使得本文方法的迭代次數(shù)更少，加快了運(yùn)行速度．

為了更加直觀地比較YALL1算法與本文方法，以初始罰參數(shù)β為例，畫出原始問題和對偶問題目標(biāo)函數(shù)值h隨迭代次數(shù)i的變化曲線，如圖1所示．

圖1 原始-對偶問題的目標(biāo)函數(shù)值隨迭代 次數(shù)的變化Fig.1 The value of the objective function of the primal- dual problem varing with the number of iterations

為了進(jìn)一步說明本文方法的有效性，將本文方法、YALL1算法、高效近塊坐標(biāo)同倫算法[5]進(jìn)行對比．本實(shí)驗(yàn)數(shù)據(jù)采用Matlab隨機(jī)生成的m×n矩陣A、加噪聲的b，以及長度n、稀疏度v的原始稀疏信號xs(利用KKT條件生成的真實(shí)信號)．高效近塊坐標(biāo)同倫算法[5]和YALL1算法以及本文方法均采用其默認(rèn)設(shè)置．對不同規(guī)模的數(shù)據(jù)進(jìn)行10次實(shí)驗(yàn)，結(jié)果見表5．(對所有實(shí)驗(yàn)結(jié)果取均值)

表5 不同算法下的恢復(fù)效果指標(biāo)Tab.5 Recovery effect indicators with different algorithms

從表5可以看出：(1)對于數(shù)據(jù)恢復(fù)的相對誤差，高效近塊坐標(biāo)同倫算法稍好一些，這是因?yàn)楦咝Ы鼔K坐標(biāo)同倫算法是基于KKT條件設(shè)計(jì)的，與生成的原始信號較吻合．另外，由于數(shù)據(jù)有噪聲，本文方法結(jié)果在容許范圍內(nèi)．(2)對于約束函數(shù)值，3種算法的結(jié)果都在精度范圍內(nèi)．(3)對于迭代步數(shù)，由于高效近塊坐標(biāo)同倫算法是將大規(guī)模問題分解成一系列小規(guī)模問題，期間濾除了大部分零分量，并且是基于KKT條件進(jìn)行更新迭代，因此大大減少了迭代步數(shù)．(4)對于迭代時(shí)間，雖然高效近塊坐標(biāo)同倫算法迭代步數(shù)很少，但每次迭代需要求解大量不同規(guī)模的線性方程組，因此用時(shí)較多．而本文方法在時(shí)間上優(yōu)勢明顯，從數(shù)值實(shí)驗(yàn)上驗(yàn)證了自適應(yīng)調(diào)整罰參數(shù)對YALL1算法的改進(jìn)成效．

4 結(jié) 語

本文研究了求解壓縮感知中的l1極小化問題的ADM，在經(jīng)典ADM框架基礎(chǔ)上，提出了一種自適應(yīng)罰參數(shù)交替方向法．通過平衡目標(biāo)函數(shù)下降速度和約束條件的滿足程度，對罰參數(shù)進(jìn)行動態(tài)調(diào)整，進(jìn)而使得信號恢復(fù)得快速、準(zhǔn)確．本文方法在一定程度上放松了對初始罰參數(shù)的范圍要求，加速了經(jīng)典ADM．下一步，將研究如何將本文方法推廣到帶有罰參數(shù)的實(shí)際優(yōu)化問題．