王 娜,呂東澔,張 勇

(內(nèi)蒙古科技大學(xué)信息工程學(xué)院,內(nèi)蒙古 包頭 014010)

1 引言

隨著科技的進(jìn)步和社會(huì)的發(fā)展,數(shù)字化的信號(hào)已經(jīng)成為當(dāng)前時(shí)代人們獲取信息的最主要來源。由于設(shè)備內(nèi)部電子元器件老化、環(huán)境中雜散的電磁波輻射干擾等因素,信號(hào)在采集、傳輸和處理的過程中必然會(huì)夾雜一定的噪聲,噪聲的出現(xiàn)會(huì)使原始信號(hào)和接收到的信號(hào)之間產(chǎn)生誤差,導(dǎo)致在后續(xù)分析時(shí)原始信號(hào)出現(xiàn)失真的現(xiàn)象。因此想要從信號(hào)中提取有用的信息,首先需要對(duì)其進(jìn)行降噪處理,信號(hào)降噪是信號(hào)處理的一個(gè)核心任務(wù),其目的是抑制含噪信號(hào)中噪聲成分并增強(qiáng)有用的信號(hào)成分,有效的降噪方法往往對(duì)更高層任務(wù)的成敗具有重要影響。相關(guān)學(xué)者在信號(hào)去噪處理方面進(jìn)行了大量的研究,例如,Gabor[1]提出了短時(shí)傅里葉變換算法,信號(hào)經(jīng)過傅里葉變換得到頻率信號(hào),可在頻域分析信號(hào)的頻譜,但是不能提取時(shí)變、非平穩(wěn)信號(hào)的局部特征;模極大值去噪[2]、小波閾值去噪[3]是基于小波思想提出的典型去噪算法,主要是通過設(shè)置小波系數(shù)來達(dá)到抑制噪聲的目的,然而在小波域中執(zhí)行軟閾值處理時(shí)會(huì)產(chǎn)生偽像現(xiàn)象。現(xiàn)有的降噪算法雖具有良好的降噪性能。但是利用這些算法在對(duì)稀疏信號(hào)進(jìn)行降噪時(shí)會(huì)產(chǎn)生一些階梯偽像噪聲,如寄生噪聲尖峰和Gibbs振蕩[4,5],導(dǎo)致原始信號(hào)中重要信號(hào)丟失。而Rudin Osher和Fatemi[6]基于變分法的思想提出基于全變分(Total Variation, TV)正則化的降噪方法,由于該方法主要是利用了信號(hào)在時(shí)域內(nèi)的光滑性,在對(duì)稀疏信號(hào)進(jìn)行降噪時(shí)能保留信號(hào)的邊緣信息并準(zhǔn)確還原原始信號(hào)。

TV去噪是一種基于基礎(chǔ)信號(hào)為分段常數(shù)信號(hào)(包含信號(hào)的導(dǎo)數(shù)是常數(shù)這一情況)的非線性濾波方法[7],這種信號(hào)出現(xiàn)在地球科學(xué)、生物物理學(xué)和其它領(lǐng)域[8],TV去噪技術(shù)還與其它方法結(jié)合使用,以處理更一般類型的信號(hào)[9-11],但是TV去噪方法在處理稀疏信號(hào)時(shí)更為突出。TV降噪算法是通過對(duì)一個(gè)特定代價(jià)函數(shù)進(jìn)行優(yōu)化,目前算法研究大多是基于凸懲罰函數(shù)提出的代價(jià)函數(shù)。例如迭代軟閾值算法(Iterating Soft Threshold Algorithm, ISTA)[12]用于解決1范數(shù)正則化反問題,ISTA具有線性收斂性能夠?qū)崿F(xiàn)對(duì)稀疏信號(hào)的降噪目的,但是需要經(jīng)過多次迭代后才能收斂且降噪精度較差。Selesnick I W等人提出的優(yōu)化最小化(Majorization Minimization, MM)[13]算法、迭代裁剪(Clipping, CLIP)[14]算法等能夠有效的找到函數(shù)的最小值,達(dá)到降噪的目的,但是不能精確估計(jì)分段常數(shù)信號(hào)峰值處的幅度。這些降噪算法的代價(jià)函數(shù)都是以凸懲罰函數(shù)為前提進(jìn)行的,因?yàn)橥箲土P中的1范數(shù)能夠提高信號(hào)的稀疏性。

研究表明非凸罰函數(shù)可以更有效的使信號(hào)稀疏化,相關(guān)學(xué)者[15-18]利用非凸懲罰函數(shù)的特點(diǎn)來估計(jì)稀疏信號(hào),但是此類算法中很少能保持代價(jià)函數(shù)的凸性。本文基于非凸懲罰項(xiàng)構(gòu)建了算法的代建函數(shù),它推廣了標(biāo)準(zhǔn)的懲罰函數(shù),通過選擇合適的優(yōu)化參數(shù)使其滿足凸性條件,使去噪問題的代價(jià)函數(shù)保持凸性。新的懲罰可以更準(zhǔn)確的估計(jì)分段常數(shù)信號(hào)中峰值處的幅值,降噪精度高于基于凸懲罰的MM優(yōu)化算法以及CLIP算法。

2 全變分降噪

TV去噪算法在降噪領(lǐng)域中應(yīng)用廣泛,該方法的具體思路是把需要降噪信號(hào)的代價(jià)函數(shù)構(gòu)建成一個(gè)凸函數(shù),對(duì)這個(gè)凸函數(shù)進(jìn)行求導(dǎo),找到凸函數(shù)的最小值,這個(gè)最小值就是利用全變分降噪得到的去噪信號(hào)。

假設(shè)噪聲數(shù)據(jù)y(n)為

y(n)=x(n)+w(n),n=0,…,N-1

(1)

其中x(n)是近似分段常數(shù)信號(hào),w(n)是高斯白噪聲。TV去噪通過解決以下優(yōu)化問題來估計(jì)信號(hào)x(n)。

(2)

其中x*為降噪后信號(hào),λ為正則化參數(shù)且λ>0用于平衡函數(shù)的權(quán)重,φn(xn)為正則項(xiàng)且φn(xn)>0,可以通過設(shè)計(jì)不同的正則化項(xiàng)設(shè)計(jì)降噪算法可以實(shí)現(xiàn)含噪信號(hào)的精準(zhǔn)降噪。

N點(diǎn)信號(hào)x由以下矢量表示

x=[x(0),…,x(N-1)]T

(3)

矩陣D被定義為

(4)

N點(diǎn)信號(hào)x的一階差分矩陣由Dx給出,其中矩陣D為(N-1)×N的一階差分矩陣。

使用以上表示方法,TV去噪問題(2)可簡(jiǎn)潔的寫成

(5)

代價(jià)函數(shù)(5)為全變分去噪模型的標(biāo)準(zhǔn)函數(shù),本文算法是在此標(biāo)準(zhǔn)函數(shù)的基礎(chǔ)上進(jìn)行的推廣設(shè)計(jì),結(jié)合非凸懲罰項(xiàng)提出的全變分(Non-Convex Penalty Total Variation, NCP-TV)降噪算法。

3 本文算法

3.1 可微分凸函數(shù)

本文所提出的非凸懲罰函數(shù)是從標(biāo)準(zhǔn)的懲罰函數(shù)中減去一個(gè)可微的凸函數(shù),其中可微的凸函數(shù)是可分離的或二元函數(shù)的和,首先定義可微的凸函數(shù)如下

令0<α<1,定義Kα:RN→R

(6)

D為式(4)中提到的一階差分矩陣。

其中Kα為2范數(shù),2范數(shù)可以尋求最優(yōu)問題且為二范數(shù)可微易求解。通過定義可知Kα存在二階導(dǎo)且二階導(dǎo)數(shù)為α,又因?yàn)?<α<1,由此可知函數(shù)Kα是可微的凸函數(shù)。

命題1:.函數(shù)Kα可定義如下

K0(x)=0

(7)

(8)

證明:

若α=0:設(shè)α=0,b=0代入式(6)可得式(7)。

若0<α<1,通過對(duì)函數(shù)(6)進(jìn)行最小化可得函數(shù)b的定義即

(9)

=[(Db)′Db]′+α(b-x)

=D′Db+b′D′D+α(b-x)

=2D′Db+αb-αx

(10)

令:s′(b)=0則

2D′Db+αb=αx

(2D′D+α)b=αx

(11)

b=(2D′D+αI)-1αx

3.2 非凸懲罰函數(shù)

研究證明非凸懲罰函數(shù)可以更有效的提高信號(hào)的稀疏性,所以本文改進(jìn)的算法是基于非凸懲罰的全變分去噪算法,此函數(shù)是標(biāo)準(zhǔn)全變分函數(shù)的一種推廣形式。

命題2:φα(x)設(shè)計(jì)為非凸函數(shù),其定義如下

(12)

D為一階差分矩陣。

命題3:當(dāng)0<α<1,非凸懲罰函數(shù)φα(x)的非凸區(qū)間定義如下

(13)

證明:

當(dāng)一個(gè)凸函數(shù)減去一個(gè)凸函數(shù)時(shí),得到的函數(shù)很可能在某一定義域內(nèi)為負(fù)值,所以定義如(13)所示的非凸區(qū)間:

由式(6)定義可知函數(shù)Kα(x)為二次函數(shù)所以

Kα≥0

(14)

又因?yàn)楹瘮?shù)φα(x)為正則化函數(shù)所以

(15)

可得

(16)

由此可推導(dǎo)出式(13)。

3.3 代價(jià)函數(shù)

命題4:.由全變分降噪定義可得基于非凸懲罰的代價(jià)函數(shù)為:

(17)

對(duì)此函數(shù)進(jìn)行降噪本質(zhì)上是求其最優(yōu)解的問題即

(18)

證明:

令n=1式(18)可寫成以下形式

(19)

構(gòu)成關(guān)于x的二項(xiàng)式函數(shù)其中c為關(guān)于b的常數(shù)項(xiàng),二次函數(shù)中二次項(xiàng)系數(shù)大于0,函數(shù)開口向上即為凸函數(shù)由此可得

(20)

得

(21)

3.4 降噪算法

命題5:設(shè)y∈RN,λ>0,算法中產(chǎn)生迭代,其迭代式x(k)定義如下所示

B(k)=y+λα[x-(2D′D)+αI]-1αx

(22)

x(k+1)=tvd(zk,λ)

(23)

證明:

NCP-TV降噪算法使用前后向分裂方法(Forward Backward Splitting, FBS)[19]對(duì)代價(jià)函數(shù)進(jìn)行優(yōu)化,FBS方法是一種迭代優(yōu)化算法,該算法最小化的函數(shù)形式如下所示

F(x)=f1(x)+f2(x)

(24)

其中f1(x)和f2(x)都是凸函數(shù),FBS方法如下

B(k)=x(k)-?f1(x(k))

(25)

(26)

通過x(k+1)進(jìn)行迭代至收斂到F(x)的最小值。

對(duì)式(17)使用FBS算法可得

Aα(x)=Fα(x)

=f1(x)+f2(x)

(27)

其中

(28)

(29)

f1(x)的梯度為

?f1(x)=x-y-?Kα(x)=x-y-α(x-b)

(30)

將式(11)中b代入(27)可得

?f1(x)=x-y-αλ[x-(2D′D+αI)-1αx]

(31)

將式(26)(28)代入上式(22)(23)可得(19)(20)形式的FBS算法。

基于全變分模型的降噪算法步驟如下:

1)輸入含噪信號(hào)y,初始化x,λ和α輸入?yún)?shù);

2)設(shè)置k=1:Nit,其中Nit為需要迭代的次數(shù);

3)輸入中間變量B(k);

4)對(duì)x(k)進(jìn)行全變分降噪;

5)k=k+1;

返回步驟2),反復(fù)迭代降噪,直到收斂到代價(jià)函數(shù)的最優(yōu)值,迭代結(jié)束,最優(yōu)降噪次數(shù)為Nit;

4 數(shù)據(jù)分析及實(shí)驗(yàn)結(jié)果

4.1 仿真信號(hào)降噪

本文實(shí)驗(yàn)通過對(duì)兩種不同的信號(hào)進(jìn)行降噪處理,第一種為仿真信號(hào)由Matlab生成的分段常數(shù)信號(hào),在此信號(hào)的基礎(chǔ)上人為的加入一組高斯噪聲,通過使用MM、CLIP降噪算法和NCP-TV降噪算法對(duì)此含噪信號(hào)進(jìn)行降噪處理,對(duì)比三種算法的降噪效果;第二種是NCP-TV算法對(duì)實(shí)際采集的肌電信號(hào)進(jìn)行降噪處理,分析降噪效果并驗(yàn)證NCP-TV降噪算法在實(shí)際應(yīng)用中的有效性。

對(duì)如圖1所示的含噪分段常數(shù)信號(hào)進(jìn)行降噪處理,這是帶有加性高斯白噪聲(σ=0.5)的信號(hào)(長(zhǎng)度N=256)。

圖1 含噪信號(hào)

如圖2顯示了CLIP、MM算法和NCP-TV算法的去噪結(jié)果,其中CLIP算法是基于TV降噪的一種降噪方法,由仿真結(jié)果可看出信號(hào)平滑區(qū)域降噪結(jié)果較差,MM算法是標(biāo)準(zhǔn)懲罰函數(shù)的一種全變分去噪方法,由圖可以看出MM算法降噪后的信號(hào)基本可以跟隨原始信號(hào),但是峰值的降噪效果有一定偏差,從NCP-TV算法的去噪結(jié)果可以看出本文算法可以更精確的估計(jì)信號(hào)峰值及平滑區(qū)域信號(hào)的幅值,降噪精度優(yōu)于上述兩種算法。

圖2 三種降噪算法對(duì)比

為了更加直觀的分析三種降噪算法的降噪效果,計(jì)算了三種算法降噪后的誤差曲線如圖3所示,讓噪聲標(biāo)準(zhǔn)差滿足0.2≤σ≤1.0,分別計(jì)算每一個(gè)σ值對(duì)應(yīng)的噪聲方差,方差用來度量隨機(jī)變量和數(shù)學(xué)期望(即均值)之間的偏離程度。由圖3可知NCP-TV算法降噪后信號(hào)的方差明顯低于CLIP算法和MM算法,說明NCP-TV算法降噪后信號(hào)與原始信號(hào)之間的誤差較小降噪精度較高。

圖3 兩種算法誤差對(duì)比

分析信號(hào)降噪常用的分析方法有信噪比(SNR)、均方根誤差(RMSE)、平均絕對(duì)偏差(MAE)等。

其中SNR為有用信號(hào)功率與噪聲功率的比值,以dB為單位,SNR越高,代表信號(hào)中噪聲值越小,故信號(hào)精度越高。其公式為

(32)

其中Ps表示有用信號(hào)功率,Pn表示噪聲的有效功率。

RMSE、MAE可以反映預(yù)測(cè)模型與實(shí)際數(shù)據(jù)誤差值大小,數(shù)值越低代表兩者越接近精度越高其計(jì)算公式如下

(32)

(33)

其中yi表示真實(shí)值即原始信號(hào)fi表示預(yù)測(cè)值即降噪后的信號(hào)。

三種算法的分析結(jié)果如表1所示:由表可知三種降噪算法的SNR值都高于含噪信號(hào),說明三種算法都達(dá)到了降噪的目的,其中NCP-TV算法的SNR值較其它兩種算法最高說明其信號(hào)中噪聲值最小降噪效果較好;由RMSE、MAE值對(duì)比可知,NCP-TV算法的RMSE和MAE值最低說明其降噪后的數(shù)據(jù)與原始信號(hào)誤差值最小,降噪精度較高。

表1 兩種不同算法特征指標(biāo)對(duì)比

4.2 肌電信號(hào)降噪

在實(shí)際采樣中,生物組織受到高壓電脈沖后期的電流信號(hào)均為方波信號(hào),求導(dǎo)信號(hào)值為零,通過判斷零信號(hào)和求導(dǎo)后的零信號(hào)占總體信號(hào)一半以上,這樣確定的信號(hào)是稀疏的。如圖4所示為實(shí)際采集到肌電信號(hào)的部分放大數(shù)據(jù)多數(shù)為方波信號(hào)可知肌電信號(hào)也為稀疏信號(hào)。

圖4 肌電稀疏信號(hào)

圖5為實(shí)際采集到肌電信號(hào)的完整數(shù)據(jù),采集過程中由于環(huán)境或設(shè)備等干擾因素肌電信號(hào)中夾雜大量干擾信號(hào),所以在對(duì)信號(hào)進(jìn)行分析之前需要對(duì)其進(jìn)行降噪處理,可采用NCP-TV算法驗(yàn)證在實(shí)際應(yīng)用中的有效性。

圖5 肌電含噪信號(hào)

通過上一章節(jié)的介紹確定NCP-TV降噪算法的可行性,為了提高肌電信號(hào)的質(zhì)量采用NCP-TV算法對(duì)實(shí)際采取的肌電信號(hào)進(jìn)行降噪,降噪實(shí)驗(yàn)如圖6所示。

圖6 肌電信號(hào)降噪

為了更加清楚地觀察降噪效果,提取信號(hào)中的前3500個(gè)數(shù)據(jù)如圖7所示為原始肌電信號(hào),通過NCP-TV算法進(jìn)行降噪,降噪后的信號(hào)如圖8所示,從降噪效果來看,所采用的降噪算法能夠消弱噪聲對(duì)肌電信號(hào)的影響。

圖7 肌電含噪信號(hào)

圖8 NCP-TV算法降噪

為了測(cè)試算法的適用性,對(duì)信號(hào)進(jìn)行特征提取并分析,肌肉疲勞標(biāo)準(zhǔn)是隨著疲勞程度的加深,即平均功率頻率(MPF)值的降低,MPF反映信號(hào)頻譜含量的變化,一些影響信號(hào)頻譜的因素會(huì)影響該參數(shù)對(duì)疲勞的表征效果,有些疲勞過程中MPF可能下降不明顯,其計(jì)算公式如下所示

(34)

式中p(f)為功率譜密度,對(duì)原始肌電信號(hào)段作頻譜變換得到,f1和f2是硬件濾波的截止頻率。

表2是肌電信號(hào)靜止?fàn)顟B(tài)與疲勞狀態(tài)降噪前后的MPF值,可以發(fā)現(xiàn)利用NCP-TV算法去噪后肌肉由靜止到疲勞狀態(tài)肌電信號(hào)MPF值降低的速率增加,從而說明降噪后的信號(hào)能更好的反映肌肉疲勞的標(biāo)準(zhǔn)。

表2 疲勞程度MPF值

表3為肌電信號(hào)靜止?fàn)顟B(tài)與疲勞狀態(tài)降噪前后的MAE值前后對(duì)比,隨著疲勞程度的加深信號(hào)的MAE值呈增長(zhǎng)趨勢(shì),由表可以看出信號(hào)降噪后其MAE值較降噪前增長(zhǎng)的速率明顯提升更好的反映了肌肉的疲勞標(biāo)準(zhǔn)。

表3 疲勞程度MAE值

由于肌電信號(hào)可以反映肌肉的活動(dòng)特征和肌肉疲勞等信息,但因其包含多種噪聲干擾無法直接進(jìn)行應(yīng)用,通過NCP-TV降噪算法對(duì)其進(jìn)行降噪處理,降噪后信號(hào)才可以進(jìn)一步進(jìn)行分析處理應(yīng)用于疲勞分析等領(lǐng)域。實(shí)驗(yàn)結(jié)果表明,NCP-TV算法在實(shí)際應(yīng)用中具有有效性。

5 結(jié)束語

本文研究了稀疏信號(hào)的精確降噪問題,提出了NCP-TV降噪算法,該算法通過構(gòu)建新的代價(jià)函數(shù)利用前向-后向方法對(duì)稀疏信號(hào)進(jìn)行降噪,并將全變分模型建為最優(yōu)模型,通過迭代尋優(yōu)獲得該模型的最優(yōu)解,該算法能夠有效的提高信號(hào)幅值的準(zhǔn)確性,仿真結(jié)果表明,NCP-TV降噪算法是有效和可行的,降噪效果優(yōu)于其它對(duì)比算法并且在實(shí)際應(yīng)用中具有有效性。