改進(jìn)的Haar 子帶變換雙濾波器自適應(yīng)算法

文昊翔，洪遠(yuǎn)泉，羅 歡

（韶關(guān)學(xué)院智能工程學(xué)院，韶關(guān)，512000）

引 言

自適應(yīng)算法[1?2]在系統(tǒng)辨識(shí)領(lǐng)域有廣泛應(yīng)用。而在某些系統(tǒng)辨識(shí)應(yīng)用(如回聲消除、有源噪聲控制)，因系統(tǒng)延時(shí)不確定，導(dǎo)致目標(biāo)系統(tǒng)的沖激響應(yīng)(Impulse response, IR)序列通常具有稀疏性[3]。即序列雖然包含大量系數(shù)，但其中大部分系數(shù)為零值系數(shù)或極小值以模擬延時(shí)；真正產(chǎn)生能量的大幅值系數(shù)在時(shí)域聚集，且數(shù)量只占系數(shù)總量極小一部分。

傳統(tǒng)的最小均方(Least mean square, LMS) 算法與歸一化LMS(Normalized LMS, NLMS)算法不考慮目標(biāo)系統(tǒng)的稀疏性，用一個(gè)長(zhǎng)的濾波器辨識(shí)IR 序列。過(guò)長(zhǎng)的濾波器導(dǎo)致自適應(yīng)系統(tǒng)收斂速度下降，計(jì)算復(fù)雜度增加[4]。

利用稀疏性可有效解決或部分緩解濾波器過(guò)長(zhǎng)導(dǎo)致的各種缺陷[5]。雙濾波器結(jié)構(gòu)(Dual filter struc?ture, DFS)[6?7]是其中一種有效的解決方案。它用一級(jí)濾波器W1(k)以低分辨率辨識(shí)目標(biāo)系統(tǒng)以定位活躍系數(shù)位置；然后一個(gè)短的二級(jí)濾波器W2(k)在時(shí)域覆蓋所有活躍系數(shù)，并精確辨識(shí)活躍系數(shù)。W2(k)通過(guò)排除大量非活躍系數(shù)參與自適應(yīng)算法以降低濾波器長(zhǎng)度。

目前對(duì)DFS 的研究主要針對(duì)W1(k)，需解決的問(wèn)題是：在不降低定位精度的條件下，如何降低W1(k)的長(zhǎng)度與計(jì)算復(fù)雜度。傳統(tǒng)的解決辦法是對(duì)W1(k)的輸入信號(hào)進(jìn)行降采樣以降低W1(k)長(zhǎng)度[8]。但降采樣將導(dǎo)致信號(hào)頻譜混疊，高頻分量丟失，不利于W1(k)定位活躍系數(shù)。

傳統(tǒng)的DCT、FFT 變換等將頻域局部化而時(shí)域信息完全丟失，因此無(wú)法應(yīng)用于時(shí)域定位。小波變換[9]可以同時(shí)兼顧信號(hào)的時(shí)頻特性。Haar 小波作為最基本的小波，具有支撐區(qū)間不重疊，計(jì)算復(fù)雜度低等優(yōu)良特性。通過(guò)求取Haar 小波與輸入信號(hào)的內(nèi)積，可以同時(shí)提取信號(hào)的時(shí)、頻信息。Haar 變換為延時(shí)估計(jì)、自適應(yīng)辨識(shí)稀疏系統(tǒng)提供了新思路[10]。

針對(duì)回聲消除應(yīng)用，Ho 等[11]較早將Haar 應(yīng)用于辨識(shí)稀疏系統(tǒng)。受此啟發(fā)，Bershad 等結(jié)合DFS，提出Haar 子帶變換 (Partial Haar transform, PHT)[12]算法。PHT 先在某個(gè)尺度對(duì)W1(k)的輸入信號(hào)進(jìn)行PHT 變換以壓縮信號(hào)長(zhǎng)度，然后再進(jìn)行自適應(yīng)濾波，并且Bershad 等從理論上證明W1(k)將最終收斂于IR 序列在該尺度下的Haar 變換系數(shù)。

以PHT 為基礎(chǔ)，Kechichian 等[13]引入模糊邏輯判定活躍系數(shù)位置；Ribas 等[14]則嘗試進(jìn)一步降低算法的計(jì)算復(fù)雜度。在眾多學(xué)者的努力下，PHT 已經(jīng)成為一類自適應(yīng)辨識(shí)稀疏系統(tǒng)的有效算法。

本文提出改進(jìn)PHT (Improved?PHT, I?PHT)算法。新算法以PHT 為基礎(chǔ)，通過(guò)為W1(k)引入時(shí)變收斂步長(zhǎng)并優(yōu)化該步長(zhǎng)使后驗(yàn)誤差最小化以加快W1(k)收斂速度；然后通過(guò)分時(shí)保存、計(jì)算W1(k) 的歸一化因子以降低算法的計(jì)算復(fù)雜度；最后以回聲消除為應(yīng)用背景對(duì)新算法進(jìn)行實(shí)驗(yàn)仿真以驗(yàn)證新算法的有效性。

1 問(wèn)題描述與相關(guān)工作

1.1 自適應(yīng)系統(tǒng)辨識(shí)

目標(biāo)系統(tǒng)IR 序列記為H=[h1,h2,…,hN]T，其中N為序列長(zhǎng)度。記k時(shí)刻輸入信號(hào)向量為X(k)=[x(k-1),x(k-2),…,x(k-N)]T，向量X(k)包含最近N個(gè)時(shí)刻的輸入信號(hào)。X(k)與H卷積，卷積結(jié)果加上背景噪聲v(k)形成系統(tǒng)輸出d(k)。

系統(tǒng)辨識(shí)需解決的問(wèn)題為：如何根據(jù)X(k)與d(k)估計(jì)H。根據(jù)自適應(yīng)理論，可以用一個(gè)FIR 濾波器W(k)=[w1(k),w2(k),…,wN(k)]以自適應(yīng)算法更新其系數(shù)以逼近、辨識(shí)目標(biāo)系統(tǒng)的IR 序列。其中LMS 算法系數(shù)更新方程為

式中：μ為收斂因子，用于平衡算法的收斂速度與辨識(shí)精度。濾波器W(k)按自適應(yīng)算法更新其系數(shù)，經(jīng)過(guò)若干次迭代，W(k)將收斂于目標(biāo)系統(tǒng)的IR 序列。若記E[]為取數(shù)學(xué)期望，則有

式(4)為基本的自適應(yīng)理論，與之相關(guān)的詳細(xì)內(nèi)容可參考文獻(xiàn)[1]。

圖1 典型稀疏系統(tǒng)的IR 序列Fig.1 Impulse response of typical sparse sys?tem

1.2 稀疏性與DFS

為保證充分建模，W(k)的長(zhǎng)度必須大于等于H的長(zhǎng)度。但當(dāng)目標(biāo)系統(tǒng)具備稀疏性時(shí)，若能獲得系統(tǒng)的延時(shí)估計(jì)，只需精確辨識(shí)活躍系數(shù)即可完成對(duì)整個(gè)目標(biāo)系統(tǒng)的自適應(yīng)辨識(shí)。一個(gè)典型稀疏系統(tǒng)的IR 序列如圖1 所示[6]。它是電子回聲路徑的IR 序列，采樣頻率8 kHz。

DFS 可利用稀疏性提高算法效率，其原理如圖2[6]所示。W1(k)以低分辨率辨識(shí)目標(biāo)系統(tǒng)并進(jìn)行延時(shí)估計(jì)；然后一個(gè)短的W2(k)在時(shí)域覆蓋所有活躍系數(shù)，并精確辨識(shí)活躍系數(shù)。

1.3 Haar 子帶變換

離散Haar 小波母函數(shù)定義為

離散Haar 小波定義為

圖2 雙濾波器結(jié)構(gòu)Fig.2 Dual filter structure

式中：t僅取整數(shù)值；m為尺度因子，取值范圍1≤m≤log2N，N為輸入信號(hào)長(zhǎng)度，m越大，對(duì)應(yīng)的Haar 時(shí)間窗越寬，以利于提取信號(hào)低頻分量，m越小，則時(shí)間窗越窄以提取高頻分量；n為時(shí)移因子，通過(guò)調(diào)整n平移小波母函數(shù)以覆蓋不同的時(shí)域范圍。

為便于闡述，以輸入信號(hào)長(zhǎng)度N=8 為例說(shuō)明。Haar 變換矩陣定義如下

式中：第1 列對(duì)應(yīng)m=3，與輸入信號(hào)求內(nèi)積后提取信號(hào)低頻分量；第2、3 列對(duì)應(yīng)m=2，分別提取信號(hào)時(shí)域前半、后半部分頻率信息；第4～7 列對(duì)應(yīng)m=1，分別提取信號(hào)對(duì)應(yīng)4 部分的高頻分量。

由式(7)可得Haar 子帶變換矩陣定義如下

利用Haar 子帶變換，Bershad 等[12]提出PHT 算法，W1(k)在Haar 域?qū)π盘?hào)進(jìn)行自適應(yīng)濾波，以降低濾波器長(zhǎng)度。系數(shù)更新過(guò)程為

同樣以N=8 為例，Haar 尺度因子M定義為M=log2N=3。當(dāng)取m=3 時(shí)，按照式(11) PHT 變換的定義，式(8)的變換矩陣將與輸入信號(hào)向量X(k)進(jìn)行矩陣相乘，W1(k)長(zhǎng)度N1壓縮為N1=1；若m=2，則N1=2；若m=1 時(shí)，則N1=4。可得，式(11)的 PHT 變換可大幅壓縮W1(k)的長(zhǎng)度。

經(jīng) Bershad 等證明[12]，W1(k)按式(11—13)進(jìn)行自適應(yīng)算法，收斂后有

因此PHT 算法既可壓縮W1(k)的長(zhǎng)度，又可以獲得與目標(biāo)系統(tǒng)H相關(guān)的信息。以W1(k)的信息為基礎(chǔ)，再利用峰值系數(shù)位置[13]、模糊邏輯[14]或移動(dòng)窗積分法[15]即可進(jìn)行延時(shí)估計(jì)以定位活躍系數(shù)。

式(5—14)為PHT 算法的基本原理，與之相關(guān)的詳細(xì)理論推導(dǎo)可參考文獻(xiàn)[12]。

2 改進(jìn)PHT 算法

以PHT(式(11—13))為基礎(chǔ)，提出I?PHT 算法提高W1(k)收斂速度。新算法先為PHT 引入時(shí)變收斂步長(zhǎng)β(k)，然后以后驗(yàn)誤差最小化為目的優(yōu)化β(k)。

2.1 推導(dǎo)過(guò)程

以式(13)為基礎(chǔ)引入時(shí)變收斂步長(zhǎng)β(k)可得

自適應(yīng)算法的先驗(yàn)誤差e1(k)定義如式(12)。后驗(yàn)誤差ε1(k)定義為

新算法目的為求取β(k)最優(yōu)值，使式（17）極大化。

其物理意義為在先驗(yàn)誤差e1(k)確定的條件下，使后驗(yàn)誤差ε1(k)能量最小化。以下為化簡(jiǎn)f(β(k))過(guò)程。由式(12)可得

同理，由式(16)可得

將式(15)代入式(19)，化簡(jiǎn)、整理后得

式(18)與式(20)兩式相減可得

由式(21)可得，f(β(k))是關(guān)于變量β(k)的二次函數(shù)。對(duì)函數(shù)f(β(k))關(guān)于β(k)求導(dǎo)，當(dāng)導(dǎo)數(shù)為零時(shí)，f(β(k))取得極大值，此時(shí)β(k)即為最優(yōu)解。即

解得

將式(23)代入式(15)可得

再次引入收斂步長(zhǎng)μi平衡算法收斂速度與穩(wěn)態(tài)失調(diào)。最后，新算法的系數(shù)更新方程為

2.2 計(jì)算復(fù)雜度優(yōu)化

I?PHT(式(24))與 PHT(式(13))相比，主要是式(24)的分母 [ZTm(k)Zm(k)](即歸一化因子)，引入額外的計(jì)算復(fù)雜度。

在NLMS 算法中，因X(k)與X(k-1)是簡(jiǎn)單的線性移位關(guān)系，可以利用式(25)快速計(jì)算歸一化因子。

但在I?PHT 中，因Zm(k)與Zm(k-1)不再滿足線性移位關(guān)系，因此[(k)Zm(k)]不能直接套用式(25)快速計(jì)算獲得，需另外尋找新的高效計(jì)算方法。

為方便討論，記Zm(k)=[z1(k),z2(k),…,zN1(k)]。定義壓縮比為r=N/N1，式中N為X(k)長(zhǎng)度，N1為Zm(k)的長(zhǎng)度。通過(guò)PHT 變換(式(11))，X(k)中相鄰的r個(gè)信號(hào)被壓縮為Zm(k)中的1 個(gè)信號(hào)。因此雖Zm(k)與Zm(k-1)不再保持線性移位關(guān)系，但Zm(k)與Zm(k-r)仍保持線性移位關(guān)系。針對(duì)該關(guān)系，建立向量P= [p1(k),p2(k),…,pr(k)]。記i=1+(kmodr)，式中mod 為求余運(yùn)算。將(k)Zm(k)值保存在pi(k)中，有pi(k)=(k)Zm(k)。即把算法歸一化因子按不同時(shí)刻分別保存在向量P中，并以r為循環(huán)。此時(shí)有

至此，新算法的實(shí)現(xiàn)可歸納為按步驟執(zhí)行式(11,12)與式(27—29)。

與 PHT 算法相比，I?PHT 每次迭代僅額外引入 3 次加/減法、2 次乘法。

3 實(shí)驗(yàn)仿真

以回聲消除為應(yīng)用背景進(jìn)行實(shí)驗(yàn)仿真以檢驗(yàn)新算法的性能。目標(biāo)系統(tǒng)H采用圖1 所示的實(shí)際回聲路徑，長(zhǎng)度為N=512?；钴S系數(shù)定位的判斷邏輯采用移動(dòng)窗積分法[15]。

壓縮比r決定W1(k)的長(zhǎng)度N1。若r取值過(guò)小，則不能有效縮短W1(k)的長(zhǎng)度；若r取值過(guò)大，則分辨率過(guò)低，無(wú)法有效定位活躍系數(shù)。本實(shí)驗(yàn)按文獻(xiàn)[13]取r=4，即N1=128。

實(shí)驗(yàn)仿真分兩部分進(jìn)行，實(shí)驗(yàn)1 檢驗(yàn)新算法對(duì)W1(k)性能的影響，實(shí)驗(yàn)2 檢驗(yàn)新算法對(duì)整個(gè)DFS 的影響。

3.1 實(shí)驗(yàn)1

以方差為1 的高斯白噪聲作為輸入信號(hào)X(k)。X(k)與H卷積，再以另一高斯白噪聲作為背景噪聲v(k)，XT(k)H與v(k)二者信噪比為 20 dB。

實(shí)驗(yàn) 1 主要比較 PHT 與 I?PHT 兩種算法。對(duì) I?PHT，μi=0.5；對(duì) PHT，μp=μi/E[(k)Zm(k)]。其中E[(k)Zm(k)]根據(jù)輸入信號(hào)的先驗(yàn)知識(shí)統(tǒng)計(jì)獲得。對(duì)W1(k)的性能評(píng)價(jià)采用兩個(gè)指標(biāo)：失調(diào)與有效定位活躍系數(shù)所需迭代次數(shù)。

(1)失調(diào)，單位為dB，計(jì)算公式為

式(30)中H1的計(jì)算公式為H1=ψM,m H。

兩算法的W1(k)失調(diào)曲線如圖3 所示。由圖3 可得，I?PHT 與PHT 相比，穩(wěn)態(tài)失調(diào)相似，但收斂速度有明顯提高。

（2)定位活躍系數(shù)所需迭代次數(shù)。當(dāng)W2(k)開(kāi)始覆蓋目標(biāo)系統(tǒng)的峰值系數(shù)時(shí)，可認(rèn)為W1(k)已準(zhǔn)確定位活躍系數(shù)。由圖4 所展示的實(shí)驗(yàn)得，為準(zhǔn)確定位活躍系數(shù)，PHT 約需408 次迭代，而 I?PHT 僅需約 101 次迭代。

圖3 一級(jí)濾波器失調(diào)曲線對(duì)比Fig.3 Misalignment curve comparison of the first-order filter

3.2 實(shí)驗(yàn)2

因W1(k)僅用于低精度辨識(shí)目標(biāo)系統(tǒng)以定位活躍系數(shù)，精確辨識(shí)目標(biāo)系統(tǒng)還需W2(k)完成，所以W2(k)的性能改善才是整個(gè)自適應(yīng)系統(tǒng)的最終目標(biāo)，對(duì)W2(k)性能的評(píng)估與W1(k)相比，更顯重要。

在計(jì)算W2(k)算法失調(diào)時(shí)，DFS 需先將短的W2(k)在首尾補(bǔ)零擴(kuò)展至長(zhǎng)度為N=512，記擴(kuò)展后的序列為W(k)，再按式（31）計(jì)算失調(diào)。

重復(fù)實(shí)驗(yàn) 1 的參數(shù)設(shè)置，PHT 與 I?PHT 的W2(k)均采用LMS 算法更新系數(shù)，W2(k)的長(zhǎng)度N2=64。μ=0.077。兩算法的W2(k)失調(diào)曲線如圖4 所示。

結(jié)合圖4 與實(shí)驗(yàn) 1 可得，因 I?PHT 可使W1(k)獲得更快的收斂速度，在更短的時(shí)間內(nèi)準(zhǔn)確定位活躍系數(shù)，所以最終導(dǎo)致W2(k)的收斂速度亦有大幅提高。因PHT 的W1(k)約需400 次迭代才可獲得穩(wěn)定的延時(shí)估計(jì)，因此其W2(k)只有在400 次迭代之后才進(jìn)入快速收斂階段。而I?PHT 的W2(k)僅需約100 次迭代即進(jìn)入快速收斂階段，獲得比PHT 更快的整體收斂速度。

為更好地模擬回聲消除系統(tǒng)工作環(huán)境，以實(shí)際語(yǔ)音信號(hào)激勵(lì)目標(biāo)系統(tǒng)，并再次進(jìn)行實(shí)驗(yàn)仿真。系統(tǒng)信噪比為20 dB。為檢驗(yàn)算法跟蹤能力，在20 s 時(shí)，目標(biāo)系統(tǒng)IR 序列整體左移 200 個(gè)單位。針對(duì)W1(k)的收斂步長(zhǎng)，I?PHT 設(shè)為μi=0.1；PHT 設(shè)為μp=μi/E[(k)Zm(k)]。W2(k)的收斂步長(zhǎng)，I?PHT 與 PHT 均統(tǒng)一設(shè)為μ=0.03。單濾波器 LMS算法作為經(jīng)典的自適應(yīng)算法亦納入對(duì)比，濾波器長(zhǎng)度設(shè)為N=512，步長(zhǎng)設(shè)置為μ=0.03/8。其他設(shè)置與實(shí)驗(yàn)1 一致。3 種算法的W2(k)失調(diào)如圖5 所示。

圖4 二級(jí)濾波器失調(diào)曲線對(duì)比Fig.4 Misalignment curve comparison of two-stage filter

圖5 3 種自適應(yīng)算法失調(diào)曲線對(duì)比Fig.5 Misalignment curve comparison of three adaptive algorithms

由于DFS 能有效降低濾波器長(zhǎng)度，其W2(k)長(zhǎng)度僅為L(zhǎng)MS 算法W(k)長(zhǎng)度的1/8。隨著濾波器長(zhǎng)度的下降，算法的收斂速度迅速提高。由圖5 得，I?PHT、PHT 的收斂速度顯著高于LMS 算法。

而同為 DFS，針對(duì)W1(k)而言，I?PHT 優(yōu)于 PHT 的地方主要有以下 3 點(diǎn)：

(1) 因I?PHT 的W1(k)能更快地定位活躍系數(shù)，因此其整體收斂速度明顯優(yōu)于PHT。圖3—5 可證明該結(jié)論。

(2) 實(shí)驗(yàn)發(fā)現(xiàn)當(dāng)輸入信號(hào)是真實(shí)語(yǔ)音信號(hào)時(shí)，I?PHT 的W1(k)收斂步長(zhǎng)μi取值范圍大于PHT 的μp。取μi<0.8 仍能保證W1(k)收斂，但 I?PHT 必須滿足μp<0.1/E[ZTm(k)Zm(k)]才能保證W1(k)收斂。因此 I?PHT 的穩(wěn)定性高于 PHT 算法。

(3) I?PHT 的步長(zhǎng)μi取值范圍有較清晰的指導(dǎo)值，其理論取值范圍為0<μi<1。當(dāng)輸入信號(hào)為高斯噪聲時(shí)，該理論范圍可保證W1(k)收斂；當(dāng)輸入為強(qiáng)相關(guān)、非平穩(wěn)的語(yǔ)音信號(hào)時(shí)，各種自適應(yīng)算法的穩(wěn)定性均會(huì)減弱，此時(shí)I?PHT 的μi實(shí)驗(yàn)取值范圍約為0<μi<0.8。但PHT 算法需先獲得E[ZTm(k)Zm(k)]的先驗(yàn)知識(shí)，才能確定μp的取值范圍。因此，與PHT 算法相比，I?PHT 更易于實(shí)現(xiàn)。

4 結(jié)束語(yǔ)

針對(duì)稀疏系統(tǒng)辨識(shí)應(yīng)用，DFS 是一類有效的解決方案。它用W1(k)定位活躍系數(shù)位置，用一個(gè)短的W2(k)精確辨識(shí)活躍系數(shù)。它通過(guò)降低濾波器有效長(zhǎng)度，以達(dá)到提高收斂速度，降低計(jì)算復(fù)雜度的目的。

本文主要針對(duì)W1(k)進(jìn)行討論，提出I?PHT 算法。新算法先對(duì)W1(k) 的輸入信號(hào)進(jìn)行PHT 變換以壓縮信號(hào)長(zhǎng)度；然后引入時(shí)變步長(zhǎng)，并以后驗(yàn)誤差最小化為目標(biāo)函數(shù)優(yōu)化時(shí)變步長(zhǎng)；最后將新算法的歸一化因子分時(shí)、循環(huán)保存到一個(gè)向量中，并通過(guò)自回歸方式維護(hù)歸一化因子以降低計(jì)算復(fù)雜度。

以回聲消除為應(yīng)用背景對(duì)I?PHT 進(jìn)行實(shí)驗(yàn)仿真，仿真結(jié)果表明，DFS 能有效降低濾波器長(zhǎng)度以提高收斂速度。而與PHT 相比，I?PHT 算法收斂速度與穩(wěn)定性均有明顯提高。