PCD算法對(duì)音頻信號(hào)降噪的參數(shù)選擇

何選森， 徐 麗， 樊躍平

(1.廣州商學(xué)院 信息技術(shù)與工程學(xué)院, 廣東 廣州 511363; 2.湖南大學(xué) 信息科學(xué)與工程學(xué)院, 湖南 長(zhǎng)沙 410082)

降噪是信號(hào)處理最重要的任務(wù)之一，而對(duì)信號(hào)的分解有助于降噪。利用字典可以將信號(hào)分解成一系列原子的線性組合[1]。以超完備字典[2]為基礎(chǔ)，冗余(稀疏)表示提供了一種強(qiáng)有力的信號(hào)模型[3]，它將信號(hào)近似為超完備字典中最少數(shù)量非零原子的線性組合[4]。信號(hào)的稀疏表示具有潛在的降噪能力[5]，以冗余表示為基礎(chǔ)的收縮[6]成為具有吸引力的降噪技術(shù)。小波收縮[7]使用正交變換實(shí)現(xiàn)降噪，而新的趨勢(shì)則是采用超完備變換。在眾多信號(hào)分解方法中，基本追蹤(basis pursuit, BP)[8]法是將信號(hào)分解為字典原子的最佳疊加，并使表示系數(shù)的l1范數(shù)(或l0范數(shù))最小化[9]。對(duì)于消除Gauss白噪聲，基本追蹤降噪(BP denoising, BPDN)[10]是最佳的。基于超完備字典和冗余表示，迭代收縮(iterative shrinkage, IS)算法[11]成為一種高效的數(shù)值技術(shù)。文獻(xiàn)[12]提出了四種IS算法，其中并行坐標(biāo)下降(parallel coordinate descent, PCD)法具有最快的搜索速度。文獻(xiàn)[13]對(duì)PCD進(jìn)行了收斂性分析，并引入正則化函數(shù)的形式，使得對(duì)坐標(biāo)優(yōu)化提供解析的解。PCD在圖像處理[14]和圖像降噪[15]方面獲得了廣泛的應(yīng)用。

盡管PCD被證明是一種優(yōu)秀的圖像降噪算法，但它在音頻降噪領(lǐng)域并沒有得到相應(yīng)的運(yùn)用。這是因?yàn)樵谝纛l或語(yǔ)音信號(hào)降噪時(shí)，需要事先對(duì)PCD算法的參數(shù)進(jìn)行設(shè)定，而目前還沒有文獻(xiàn)給出相關(guān)參數(shù)的設(shè)置標(biāo)準(zhǔn)。正是基于這種考慮，本文利用仿真實(shí)驗(yàn)，研究PCD算法相關(guān)參數(shù)的選擇依據(jù)。主要通過(guò)比較音頻信號(hào)外加Gauss噪聲的強(qiáng)度、軟收縮函數(shù)的閾值、算法的迭代次數(shù)、分解信號(hào)的超完備字典與PCD降噪性能之間的關(guān)系，給出以上參數(shù)選擇的基本原則，以擴(kuò)展PCD算法的應(yīng)用范圍。

1 去噪背景

設(shè)x∈RN為N維觀測(cè)信號(hào)矢量，s∈RN為信源，T∈RN×N為正交變換，則觀測(cè)信號(hào)表示為x=Ts。對(duì)于稀疏信源s，在同一時(shí)刻，只有很少分量si∈s是非零值，其它多數(shù)分量都為零值。在實(shí)際應(yīng)用中，信源是不可觀測(cè)的，只能通過(guò)傳感器采集來(lái)獲得觀測(cè)信號(hào)x。由于采集信道噪聲的存在，x中不可避免地包含了噪聲。對(duì)于稀疏信號(hào)，小波分析中的收縮降噪方法[16]是假定絕對(duì)值很小的分量為噪聲并將其設(shè)置為零，而僅保留較大絕對(duì)值的分量以實(shí)現(xiàn)降噪目的。

若信號(hào)x∈RN可以由幾個(gè)(單位能量的)原子dk的線性組合來(lái)近似，則由所有原子組成的集合{dk：k=1, 2, …,K}就是一個(gè)字典D∈RN×K。設(shè)矢量z∈RK是稀疏的，則x可表示為：

x=Dz

(1)

式中：z是信號(hào)x的重構(gòu)系數(shù)。

當(dāng)K=N時(shí)，字典D是正交的，則式(1)的表示是唯一的。當(dāng)K>N時(shí)，字典D是超完備的，這種情況下式(1)存在有無(wú)限多個(gè)可能的系數(shù)集。這種非唯一性表示提供了極大的選擇性，即從眾多表示中選擇最適合的表示。尋找具有最少量非零系數(shù)的解是一種非常具有吸引力的優(yōu)化技術(shù)[17]。當(dāng)使用超完備字典分解信號(hào)時(shí),式(1)給出了欠定的線性方程組，它的最稀疏解可表示為如下的優(yōu)化問(wèn)題[17]：

(2)

或者

(3)

式中：ε是誤差容限；‖?‖p為lp范數(shù)；l0范數(shù)是對(duì)一個(gè)向量非零元素的計(jì)數(shù)。

式(2)～(3)中的優(yōu)化問(wèn)題可以通過(guò)正交匹配追蹤法[12]求解。基本追蹤(BP)是尋找系數(shù)具有最小l1范數(shù)的信號(hào)表示[9]，若在式(2)中用l1范數(shù)代替l0范數(shù)，則有：

(4)

式中：l1范數(shù)是求向量各元素絕對(duì)值的和。

因l1范數(shù)是l0范數(shù)的凸化[17]，也就是說(shuō)式(4)是式(2)的凸化。為了更接近于實(shí)際應(yīng)用的環(huán)境，考慮一個(gè)含噪聲的信號(hào)矢量

y=x+v

(5)

式中：y∈RN是信號(hào)x的含噪聲觀測(cè)量；v∈RN是Gauss零均值白噪聲向量，并且限制它的能量為有限值。對(duì)信號(hào)y降噪的優(yōu)化問(wèn)題表示為[18-19]：

(6)

此問(wèn)題可以由Lagrange形式替代[9]：

(7)

2 超完備字典

在PCD算法中，超完備字典起著非常重要的作用，選擇合適的超完備字典是PCD實(shí)現(xiàn)收縮降噪的前提條件。

在音頻信號(hào)處理領(lǐng)域，有兩種典型的超完備字典，分別是離散余弦變換(discrete cosine transform, DCT)字典[20]和Gabor字典[20]。

一個(gè)窗口式DCT字典表示為DDCT=[d1,d2, …,dK]，它的第k個(gè)原子定義為[20]

(8)

式中：1≤k≤K；0≤t≤N-1；N是每個(gè)音頻幀的長(zhǎng)度；K是DCT字典的大小(尺寸)，也就是離散頻率的數(shù)量；wd是加權(quán)窗口。

DCT原子通過(guò)引入初始相位φ可推廣為Gabor原子，這樣做是為了增強(qiáng)信號(hào)結(jié)構(gòu)的潛在適應(yīng)性和信號(hào)的稀疏性。

窗口式Gabor字典表示為DG=[d(k,φ)](k,φ)∈Ψ，原子以連續(xù)參數(shù)集Ψ=[1,K]×[0, 2π]為索引，定義為[20]：

(9)

式中：1≤k≤K；0≤t≤N-1；N是幀的長(zhǎng)度；K是Gabor字典的尺寸大小；φ是它的相位；wd是其加權(quán)窗口。在無(wú)任何近似的情況下，Gabor字典連續(xù)索引的原子d(k,φ)的分解可以用離散字典中的原子對(duì)來(lái)表示。一對(duì)的原子可以是相同頻率且相位為零的余弦和正弦對(duì)：

(10)

(11)

一般地，對(duì)DCT和Gabor詞典的列向量(音頻幀)都進(jìn)行歸一化處理以簡(jiǎn)化其計(jì)算。例如式(10)～(11)所示Gabor字典的正弦和余弦原子，采用它們對(duì)應(yīng)的單位范數(shù)版本。

3 PCD算法

實(shí)際上，式(7)就是PCD算法的目標(biāo)函數(shù)[12]。對(duì)于含噪聲信號(hào)y，去除其中的Gauss噪聲可歸納為對(duì)以下函數(shù)的最小化[12]：

(12)

當(dāng)采用最大后驗(yàn)概率估計(jì)時(shí)，函數(shù)f(x)可從Bayes的觀點(diǎn)來(lái)解釋。在式(12)中，等號(hào)右邊第一項(xiàng)通常被稱為對(duì)數(shù)似然，它描述干凈信號(hào)x與其含噪聲觀測(cè)量y之間的關(guān)系，而第二項(xiàng)ρ(x)表示未知信號(hào)x的先驗(yàn)值。對(duì)不同的函數(shù)ρ(x)，式(12)中f(x)具有不同的展開式。如果字典D是超完備的，則式(12)的目標(biāo)函數(shù)將變成以下的優(yōu)化問(wèn)題[12]：

(13)

將優(yōu)化式(13)與式(7)比較可知，在PCD算法中，函數(shù)ρ(z)采用的是l1范數(shù)。在文獻(xiàn)[12]中詳細(xì)地介紹了PCD算法的推導(dǎo)過(guò)程，這里僅給出最終的更新迭代式：

zk+1=zk+μ{S[v(D,y,zk)]-zk}=zk+μhk

(14)

式中：hk是可計(jì)算的向量；zk是z的第k次解，zk+1是其更新值；μ是步長(zhǎng)參數(shù)，若μ=1，則采用完全收縮，若μ<1，則效果為松弛步長(zhǎng)；S是軟收縮(閾值)函數(shù)，其定義為：

(15)

式中：θ為閾值。

在式(14)中，向量v(D,y,zk)定義為：

v(D,y,zk)=diag-1(DTD)DT(y-Dzk)+zk

(16)

由式(14)可知收縮函數(shù)S以閾值λ?diag-1(DTD)·1K對(duì)向量v(D,y,zk)逐項(xiàng)操作，其中1K表示元素為1的K×1矩陣(列向量)。由于字典D的列是歸一化的，則diag-1(DTD)=I，其中I為單位矩陣。于是，式(16)可簡(jiǎn)化為：

v(D,y,zk)=DT(y-Dzk)+zk

(17)

PCD算法的搜索速度可以通過(guò)線搜索[11]或序列子空間優(yōu)化[19]來(lái)進(jìn)一步加快。線搜索即一維搜索，是指單變量函數(shù)的優(yōu)化過(guò)程。在PCD算法中，線搜索用于尋找最佳的步長(zhǎng)μ，其目的是求解以下的優(yōu)化問(wèn)題：

(18)

而式(18)的解是通過(guò)求下面的方程獲得：

(19)

式中：sign為符號(hào)函數(shù)。式(15)也具有收縮的結(jié)構(gòu)。

4 仿真結(jié)果與分析

4.1 仿真環(huán)境與性能指標(biāo)

為了驗(yàn)證PCD算法對(duì)音頻信號(hào)的降噪性能，通過(guò)一系列仿真實(shí)驗(yàn)，歸納出算法相關(guān)參數(shù)的選擇原則。仿真PC機(jī)配置為：Intel(R) Celeron(R) 1007U-1.5 GHz的CPU，4 GB內(nèi)存，Windows 10操作系統(tǒng)。所有仿真都是在MATLAB 9 (R2016a)上運(yùn)行。用于測(cè)試的音頻源信號(hào)為兩通道的音樂信號(hào)，其類型為wav，信號(hào)的采樣率為44.1MHz。每個(gè)信源樣本數(shù)為L(zhǎng)=216=65 536。五個(gè)音頻信源X=[x1,x2,x3,x4,x5]T的時(shí)域波形見圖1。

圖1 五路音樂音頻源信號(hào)的波形Fig.1 The time domain waveform of 5 source signals

為了測(cè)度PCD算法的降噪性能，采用兩種經(jīng)典的技術(shù)指標(biāo)，分別是信噪比(signal to noise ratio, SNR)和均方根誤差(root mean square error, RMSE)：

(20)

(21)

信號(hào)降噪的仿真過(guò)程如下。在時(shí)域中音頻信源與Gauss噪聲相加生成含噪聲的觀測(cè)信號(hào)。每個(gè)觀測(cè)信號(hào)乘以長(zhǎng)度為64的Hamming窗以生成音頻(語(yǔ)音)幀，使連續(xù)兩幀之間的重疊樣本數(shù)為32，即每一幀的前后各有16個(gè)樣本與相鄰幀的樣本重疊。音頻幀的分割過(guò)程是通過(guò)使用Voicebox中的MATLAB函數(shù)enframe來(lái)實(shí)現(xiàn)。將PCD算法應(yīng)用于分幀后的觀測(cè)信號(hào)進(jìn)行降噪，每次使用PCD是對(duì)一幀進(jìn)行處理。對(duì)每個(gè)降噪的幀乘以逆Hamming窗，且在相鄰幀重疊的樣本中，分別在每一幀兩邊各去除50%的樣本，僅保留每一幀中間部分的樣本，最后將經(jīng)過(guò)處理的各個(gè)片段簡(jiǎn)單串聯(lián)起來(lái)即可得到重新合成的降噪后的信號(hào)。

PCD算法的迭代過(guò)程如下。設(shè)定最大迭代次數(shù)M。將迭代次數(shù)初始化為k=1，仿真開始。在PCD算法每次迭代后計(jì)算并記錄相應(yīng)的RMSE和SNR值。然后，將迭代次數(shù)更新為k←k+1，以重復(fù)執(zhí)行PCD算法。最后，判斷仿真過(guò)程是否結(jié)束，若k

4.2 噪聲強(qiáng)度對(duì)降噪性能的影響

在PCD算法的迭代中，需要預(yù)先給出軟收縮函數(shù)S的閾值λ。在確定最優(yōu)的閾值λ時(shí)，可以使用全局閾值規(guī)則和最小最大閾值規(guī)則[7]。本文采用全局閾值規(guī)則，其λ的取值范圍為[7]：

(22)

式中：L是信號(hào)的樣本長(zhǎng)度；σ是Gauss噪聲的標(biāo)準(zhǔn)差。從式(22)可知，要確定閾值λ，就必須先要給出噪聲的強(qiáng)度σ。

音頻信源在被加入Gauss噪聲后生成觀測(cè)信號(hào), 而加入的噪聲強(qiáng)度采用(input)信噪比來(lái)度量。為確定噪聲強(qiáng)度對(duì)降噪性能的影響，本小節(jié)采用小波軟閾值降噪法。小波的類型為db45，當(dāng)觀測(cè)信號(hào)的輸入(input)SNR=10, 15, 20dB時(shí)，利用小波降噪獲得相應(yīng)的輸出SNR值和RMSE值。在這個(gè)仿真中，對(duì)于五個(gè)信源的加噪、降噪過(guò)程重復(fù)進(jìn)行了100次實(shí)驗(yàn)。在每次實(shí)驗(yàn)后，分別計(jì)算并記錄全部五個(gè)觀測(cè)信號(hào)在一次實(shí)驗(yàn)中平均的輸出SNR和RMSE值。仿真實(shí)驗(yàn)的結(jié)果見圖2。

圖2 五個(gè)信號(hào)的平均輸出SNR和RMSE值Fig.2 Average output SNR and RMSE of five signals

從圖2的平均性能指標(biāo)可以看出，隨著觀測(cè)信號(hào)中SNR值的降低(即噪聲強(qiáng)度增大)，小波軟閾值方法的降噪性能變差。特別地，當(dāng)(input)SNR=10時(shí)，降噪后輸出的平均SNR<9.8，即經(jīng)過(guò)降噪處理的信噪比小于降噪前的信噪比。這說(shuō)明由于噪聲強(qiáng)度過(guò)大，信號(hào)被嚴(yán)重污染而不能恢復(fù)。對(duì)音頻信號(hào)來(lái)說(shuō)，(input)SNR=10是一個(gè)非常惡劣的環(huán)境。

為便于分析和比較，對(duì)每一個(gè)輸入的SNR值，分別計(jì)算施加到每個(gè)信源上的噪聲強(qiáng)度σ。其對(duì)應(yīng)的結(jié)果見表1。

表1 輸入SNR與噪聲強(qiáng)度σ的對(duì)應(yīng)關(guān)系Tab.1 Correspondence between input SNR andnoise intensity

從表1的數(shù)據(jù)可以看出，對(duì)于相同的輸入SNR值，不同信源被施加的噪聲強(qiáng)度是不相同的，這是因?yàn)槊總€(gè)信源的幅度(或功率)都不相同所致。對(duì)于(input)SNR=10的最壞情況，施加到信號(hào)x2,x3,x4,x5上的噪聲強(qiáng)度接近0.1，而施加到信號(hào)x1的噪聲強(qiáng)度為0.2。為統(tǒng)一起見，在下面的仿真實(shí)驗(yàn)中，將相同的噪聲強(qiáng)度(σ=0.1)施加于五個(gè)信源形成含噪聲的觀測(cè)信號(hào)，也就是PCD算法的輸入信號(hào)。

4.3 軟閾值對(duì)降噪性能的影響

在給定噪聲強(qiáng)度σ之后，利用式(22)就可以求出軟收縮函數(shù)S的閾值λ。

因?yàn)樾盘?hào)的采樣點(diǎn)數(shù)為：L=216=65 536，因此有(2log2L)1/2=4.7096，從文獻(xiàn)[7]中的Table 2可以知道，在給定σ=1時(shí)，其參考的閾值為λ*=3.310。由于本文中選取σ=0.1，將它代入到式(22)可得結(jié)果：λ≤σ×3.310=0.331 0，顯然這只是λ的一個(gè)取值范圍。為了找出PCD算法的最佳閾值λ，通過(guò)以下的仿真實(shí)驗(yàn)來(lái)進(jìn)行確定。其實(shí)驗(yàn)條件為：每個(gè)觀測(cè)信號(hào)都被分割成長(zhǎng)度為64的音頻幀，設(shè)置PCD算法的最大迭代次數(shù)為M=5。軟閾值λ的值從0.05到0.331 0不斷地變化，其更新步長(zhǎng)為0.01。超完備字典分別采用DCT和Gabor，對(duì)每個(gè)λ值和每個(gè)字典僅進(jìn)行一次實(shí)驗(yàn)。圖3給出了在這種條件下PCD算法輸出的SNR和RMSE值。

圖3 五個(gè)信號(hào)對(duì)不同λ值的輸出SNR和RMSEFig.3 Output SNR and RMSE of five signals for different λ’s

從圖3的結(jié)果可以看出，無(wú)論采用DCT還是Gabor字典，PCD算法對(duì)每個(gè)觀測(cè)信號(hào)的降噪性能，即輸出的SNR和RMSE值隨著閾值λ的增加都有明確的變化趨勢(shì)。顯然地，當(dāng)λ=0.1時(shí)每個(gè)信號(hào)的輸出SNR達(dá)到相應(yīng)的最大值，而每個(gè)信號(hào)的輸出RMSE也達(dá)到相對(duì)應(yīng)的最小值。這說(shuō)明當(dāng)選取閾值λ=0.1時(shí)，PCD算法能獲得最好的降噪效果。

4.4 超完備字典對(duì)降噪性能的影響

從圖3還可看出，對(duì)于不同的超完備字典，PCD算法的去噪性能是不相同的。因此，在給定噪聲強(qiáng)度σ和軟收縮函數(shù)的閾值λ之后，還需要分析超完備字典對(duì)PCD算法降噪性能的影響，以確定最適合于PCD算法的字典。

在這個(gè)仿真中，PCD算法的最大迭代次數(shù)設(shè)置為M=15。分別采用DCT和Gabor字典，對(duì)每個(gè)觀測(cè)信號(hào)利用PCD算法進(jìn)行降噪處理，并計(jì)算和記錄每個(gè)信號(hào)在每次迭代中的性能指標(biāo)。圖4給出了實(shí)驗(yàn)的結(jié)果(輸出的SNR和RMSE值)。

圖4 五個(gè)信號(hào)對(duì)不同字典的輸出SNR和RMSEFig.4 Output SNR and RMSE of five signals for different dictionaries

從圖4的輸出指標(biāo)波形可以看出，隨著迭代次數(shù)的增加，PCD算法的降噪性能逐漸得到提高，即PCD算法隨著迭代次數(shù)的增加而逐漸收斂。觀察第一個(gè)信號(hào)(即signal1)的波形可知，使用Gabor字典所得到的SNR值要大于使用DCT字典對(duì)應(yīng)的SNR值。同樣地，使用Gabor字典得到的RMSE值要小于使用DCT字典的RMSE值。此外，信號(hào)2至信號(hào)5輸出的SNR和RMSE值也都具有相類似的關(guān)系。這意味著使用Gabor字典，PCD算法可以獲得更好的降噪性能。

4.5 算法迭代次數(shù)對(duì)降噪性能的影響

在給定噪聲強(qiáng)度σ，軟閾值λ，并采用Gabor超完備字典的情況下，研究PCD算法的降噪性能與算法迭代次數(shù)的關(guān)系。在這個(gè)仿真中，PCD算法的最大迭代次數(shù)設(shè)為M=50，Gabor字典的尺寸為N×K=64×1 024，它的加權(quán)窗口wd為長(zhǎng)度N=64的對(duì)稱正弦窗口。在PCD算法的每次迭代過(guò)程中，分別記錄每個(gè)信號(hào)的輸出SNR和RMSE值，其結(jié)果見圖5。

圖5 五個(gè)信號(hào)在每次迭代后的輸出SNR和RMSEFig.5 Output SNR and RMSE of five signals after every iteration

從圖5的輸出結(jié)果可以看出，對(duì)于每一個(gè)觀測(cè)信號(hào)，在經(jīng)過(guò)大約25次迭代后，PCD算法降噪的指標(biāo)SNR和RMSE值達(dá)到相對(duì)應(yīng)的穩(wěn)定值。也就是說(shuō)，以降噪為目標(biāo)的PCD算法，其收斂條件為算法需要經(jīng)過(guò)25次的優(yōu)化迭代。

5 結(jié) 語(yǔ)

針對(duì)音頻信號(hào)的降噪問(wèn)題，本文通過(guò)仿真，比較分析了算法的相關(guān)參數(shù)對(duì)PCD降噪性能的影響。仿真結(jié)果表明如下結(jié)論。

1) Gauss白噪聲的標(biāo)準(zhǔn)差為0.1(或輸入的SNR=10)時(shí)會(huì)造成小波降噪性能的嚴(yán)重下降，這是因?yàn)橐纛l信號(hào)本身的幅度(功率)較小所致，此結(jié)論對(duì)音頻信號(hào)處理具有指導(dǎo)性意義。

2) 當(dāng)軟收縮函數(shù)的閾值為0.1、采用Gabor字典時(shí)PCD算法能取得很好的降噪性能。

3) 當(dāng)?shù)螖?shù)為25時(shí)，PCD算法能獲得穩(wěn)定的收斂。