胎心音是一種非線性非平穩信號,采集時含有大量噪聲,信號處理時需要一種良好的除噪方法。本文先用截止頻率為200 Hz的巴特沃斯低通濾波器和重采樣做預處理,再使用經驗模態分解法(EMD)做信號分解,選擇含有目標信號的分量進行類小波軟閾值自適應算法除噪,最后組合重構得到除噪后的胎心音信號。模態分解時,使用添加掩模信號等方法消除模態混疊,用鏡像延拓法消除端點效應,并引用Rilling的研究設定停止規則。該除噪方法一次性地消除了基線漂移和噪聲,同小波變換(WT)、數學形態學(MM)和傅里葉變換(FT)比較,信噪比(SNR)明顯改善,均方根誤差(RMSE)最小,能夠滿足實際應用的需要。
引用本文： 劉俏俏, 譚志向, 張懿, 王華. 基于經驗模態分解的胎心音信號除噪方法. 生物醫學工程學雜志, 2015, 32(4): 740-745,772. doi: 10.7507/1001-5515.20150135 復制
引言
在胎兒心臟中，心肌本身的舒張和瓣膜的關閉以及血流沖擊所產生的聲音叫做胎心音(fetal heart sound,FHS)[1]。胎心音包含較響的第一心音(S1)和第二心音(S2)，以及很微弱的第三心音(S3)和第四心音(S4)。Nagel[2]研究發現，胎心音采集信號主要含有母體心音、內部噪聲和胎心音S1與S2，其中內部噪聲和母體心音的頻譜主要在20 Hz以下，S1在23、34和75 Hz有三個譜峰，S2在30 Hz處有一個譜峰，并且噪聲頻段與胎心音頻段重合。同時，他還通過統計研究發現，隨著懷孕周數增加，S2的頻率逐漸降低，到一定階段趨于穩定，而S1的頻率先上升，在34周時達到最大，然后又逐漸下降，最后趨于穩定，即S1與S2的頻率與懷孕周數有較明顯的復雜函數關系。這都說明胎心音采集信號是一種非線性非平穩信號，要從中提取出有用的胎心音信號有很大困難[3]。
目前，因為胎心音檢測技術比超聲、超聲多普勒和胎心電等技術更安全，也更容易實現，故而受到越來越多的關注，但是在母體腹部采集得到的信號中含有大量噪聲，胎心音主要成分都淹沒在噪聲中，所以胎心音信號除噪方法是胎心音檢測技術的一個研究重點。周禮杲等[4]分析了胎心音的頻譜特點，提出了使用頻率濾波方法分離S1和S2，通過帶通濾波和匹配濾波等方法實現胎心音除噪，但是胎心音信號非線性非平穩的特點，即不同的胎兒在不同的時間或者環境，胎心音頻率均會發生變化，使該方法難以普遍適用。Varady等[5]使用帶通濾波器、小波變換(wavelet transform，WT)和互相關等方法做信號處理。 Chourasia等[6]使用小波工具箱和自適應濾波等方法對采集信號除噪。小波變換的本質是基于傅里葉的，對非平穩信號處理能力有限，并且它存在小波母函數選擇困難等問題。Ruffo等[7]使用匹配濾波、帶通濾波、Teager能量算子(Teager energy operator,TEO)和自相關等方法做信號處理和計算胎心率，該方法雖然較容易實現，但是在使用TEO能量算子時對處理信號的信噪比(signal noise ratio，SNR)要求較高。希爾伯特黃變換(Hilbert-Huang transform,HHT)是一種新的非線性非平穩信號處理方法，它在分解信號時能得到復雜信號的瞬時頻率，這與小波變換不同，而且它在地震、雷達、振動和成人心音等信號處理中得到了廣泛的應用，并取得了良好的效果，但是目前還沒有人將該方法應用于胎心音處理研究，所以本論文在比較和分析了該方法的理論和優點后，將HHT的核心——經驗模態分解(empirical mode decomposition method,EMD)應用于胎心音信號除噪。
1 經驗模態分解
由文獻[8]，HHT包含經驗模態分解EMD和希爾伯特變換(Hilbert transform,HT)，其中EMD能夠對復雜信號進行多層次分解，獲得若干代表瞬時頻率分量的本征模態函數(intrinsic mode function,IMF)和一個殘余分量。EMD分解信號時[9]，可以依據IMF的兩個條件[8]設置分解條件和停止規則，使得每個IMF只包含一個基本函數模式的振蕩波形。所以，任一復雜信號被EMD分解后，均能得到多個局部表征信號特點的IMF，這時可以依據需要，去掉不需要的IMF，組合重構剩下的IMF，即可實現信號除噪的目的[10]。
EMD在分解信號時，存在三個問題——端點效應[11]、停止規則和模態混疊[12]，需要采用合適的方法處理，否則信號分解可能不充分或者不準確。針對本論文，使用EMD分解胎心音采集信號時，端點效應問題是難以確定信號端點是否是極值點，論文中使用鏡像延拓法，在端點附近尋找合適的位置鏡像延拓一段信號，這樣在信號極值點包絡插值擬合時，就可以避免判斷端點是否是極值點的問題；EMD分解信號是為了得到單一分量，所以需要按照IMF的定義設定停止規則，同時，在實際使用中，為了提高信號處理效率，沒必要將信號徹底分解，這時需要對停止規則做些改變；模態混疊問題是因為噪聲頻率與胎心音有用信號頻率有重疊，可能造成信號分解不充分，使得某些分量中還含有其它頻率成分的信號，造成后期信號除噪困難，論文將采用加入掩模信號的方法消除模態混疊。雖然EMD分解信號需要注意一些問題，但是它具有小波變換多分辨率的優勢，同時克服了小波變換中小波基選擇的缺陷，對胎心音等非線性非平穩信號的處理更有優勢[13]。
2 胎心音除噪方法
本論文的胎心音采集信號均是采用參考文獻[14]的裝置和方法采集，使用一個高靈敏度的3M改裝聽診器、NI myDAQ數據采集卡和LabVIEW(Laboratory Virtual Instrumentation Engineering Workbench)控制軟件組成采集系統，實現了在母體腹部采集信號并以WAV格式保存原始音頻信號。實驗對象是在超聲醫療國家工程研究中心的員工及員工家屬內征集的10個孕婦志愿者(孕期從25周到35周)，受試者本人與家屬均簽訂知情同意書。按照文獻[14]的裝置和方法共采集200例數據，經過篩選，最后得到150例有用的胎心音數據。
該方法采集到的胎心音信號含有大量噪聲，包括采集電路引入的基線漂移噪聲，為了除去這些噪聲，如圖 1所示，主除噪部分采用EMD對信號進行模態分解，在分解時考慮消除模態混疊和端點效應的影響，并設置合適的停止規則，接著將EMD分解出的含有有用信號信息的分量采用類小波軟閾值自適應除噪，最后將除噪后的選擇單分量重構組合，即可得到除噪信號。
圖1
胎心音信號處理流程圖
Figure1.
Processing flow chart of FHS
2.1 預處理
為了節省資源，提高后期軟件程序的執行效率，針對胎心音信號頻段主要在20~200 Hz的特點，預處理中先采用截止頻率為200 Hz的巴特沃斯低通濾波器，然后將信號進行重采樣，以降低數據量。
胎心音采集信號采樣率為8 000 Hz，時長為5 s，圖 2是原胎心音采集信號經低通濾波和分頻重采樣后，得到信號采樣率2 000 Hz，并且由圖中可以看到相同時間長度上，信號的點數有所減少。因為胎心音的主頻段在200 Hz以內，遠遠小于重采樣后信號頻率2 000 Hz，所以重采樣后的信號依然包含胎心音信號的全部信息，頻譜沒有變化，而且數據量減少后，程序運行效率得到了提高。
圖2
信號預處理
Figure2.
Signal preprocessing
2.2 經驗模態分解
使用EMD對信號進行分解得到多個IMF和一個殘余分量的過程是一個篩選過程，篩選時要保證信號盡可能地分解為單分量信號，即分解的單分量信號不僅需要滿足IMF的條件，還需要消除EMD存在的三個問題。本研究將含噪信號按照以下步驟進行分解：
(1)使用EMD分解之前考慮端點效應影響，本論文使用鏡像延拓法處理信號端點效應。
(2)為了抑制模態混疊，在采集信號x分解時引入一個正弦掩模信號y，如公式(1)，它的頻率一般大于x的最高頻率，因為EMD分解時，模態混疊主要是由于頻率重合引起的，并且主要在高頻段，當加入一個高頻信號后，即可在EMD的低階分解時將高頻和低頻分開，從而達到消除模態混疊的作用，因為胎心音信號有用頻段在200 Hz以下，所以這里選擇正弦掩模信號的頻率為200 Hz，如圖 3為所選擇的正弦掩模信號波形圖。
| $y = mean\left( x \right) * \sin \left( {2\pi ft} \right)$ |
圖3
掩模信號
Figure3.
Masking signal
(3)加入掩模信號后，分解方法為：
| $\begin{array}{l} imf1 = emd\left( {x + y} \right)\\ imf2 = emd\left( {x - y} \right)\\ imf = \left( {imf1 + imf2} \right)/2 \end{array}$ |
上式中，imf1、imf2和imf均為單分量信號組成的矩陣。并且，使用x+y和x-y做EMD分解后求均值，這樣處理的原因是如果信號x不存在模態混疊，加入的y相當于一個噪聲，經過這樣的處理，可以消除y的影響，雖然這樣處理的弊端是增加了計算量，但是這種方法能夠最大限度地降低模態混疊的影響，在資源允許的情況下，是可以接受的。
(4)EMD分解時還需要設定合適的停止規則，按照IMF的條件設定判定規則，按照EMD分解結果設定算法結束規則，參考Rilling等[15]的研究，本論文設定的停止規則為：
①待處理數據的極值點個數不能少于3個，這是判定EMD結束的規則之一。
②假設三個變量T1、T2和α，且T1<T2，三個變量取值均在[0, 1]內，令m(t)、emax(t)和emin(t)分別表示待處理信號的包絡均值和上、下包絡振幅。然后要求滿足關系：
A、x(t)分解中的分量中任意一點值都小于T2，這是為了保證不出現局部過于不對稱的點。
B、x(t)分解中的分量的點超過T1的比例不大于α，這是為了保證對稱性，防止EMD分解時篩選次數太多，降低算法效率。
C、極值點和過零點數目相差不大于1，這是簡單中值規則，也是EMD分解結束判據之一。
經過前人大量的實驗和分析，并反復驗證掩模信號和停止規則，一般選擇T1=0.05，T2=0.5，α=0.05，可得到信號所有IMF分量和殘余分量。
使用該分解過程，不僅能很好地抑制模態混疊的影響，而且保證了信號的恰當分解，對后面的信號處理有重要意義。圖 4是預處理后的信號的EMD分解圖，信號被分解為9層，得到8個IMF分量和一個剩余分量，其中目標信號主要在IMF2、IMF3和IMF4分量中，而噪聲主要在IMF1以及IMF5之后的分量中。
圖4
EMD結果圖
Figure4.
Result of EMD
2.3 類小波軟閾值自適應算法
當含噪信號經過EMD分解得到多個單分量后,閾值函數的選擇以及閾值選取是各個單分量除噪的兩個關鍵問題。軟閾值函數是將信號絕對值小于閾值的數置為零，對大于閾值的數進行閾值收縮處理，適合用于心音的處理。軟閾值函數表達式為：
| ${\eta _T}\left( {\rm{x}} \right) = \left\{ \begin{array}{l} {\rm{sign}}\left( x \right) * \left( {\left| x \right| - T} \right),\left| x \right| \ge T\\ 0,\left| x \right| \le T \end{array} \right.$ |
其中T為閾值，ηT(x)是除噪函數，sign(x)是符號函數，x為輸入信號。
在閾值選擇時，若閾值太小，則可能除噪后信號中依然有大量噪聲；若閾值太大，除噪后信號的方差減小造成大偏差，有可能使信號過于平滑。考慮到EMD分解信號的特點，閾值采用基于高斯白噪聲的自適應閾值，閾值函數公式為[16]：
| ${T_k} = c * \sqrt {{E_k} * 2 * \ln N} ,k = 1,2, \cdots $ |
其中c表示閾值系數，一般取值范圍為(0.5，0.7)，經過大量實驗分析，這里選擇c=0.6，N表示信號長度，Ek(k=1,2,…)表示第k個IMF的噪聲強度，其計算公式為[17]：
| ${E_k} = \left\{ \begin{array}{l} \sigma _1^2,k = 1\\ \frac{{\sigma _1^2}}{{0.719}} * {2.01^{ - k}},k = 2,3, \cdots \end{array} \right.$ |
其中σ12是第一個IMF分量的方差。
如圖 4所示，信號經過EMD分解后，IMF1分量主要是噪聲，而后面幾個分量和剩余分量是低頻部分或趨勢信號，可以直接刪除，目標信號主要在IMF2、IMF3和IMF4分量中，因此在使用類小波軟閾值自適應函數除噪時，沒必要將所有分量都進行除噪，僅需要將含有目標信號主要成分的分量選擇出來，并對這些分量分別進行自適應除噪，然后將除噪后的分量組合重構，便可得到除噪信號。
如圖 5所示，對IMF2、IMF3和IMF4分量進行除噪，進行重構組合后，得到的干凈信號相比于原始信號，不僅消除了基線漂移的干擾，也很好地濾除了噪聲。
圖5
除噪結果圖
Figure5.
Result of denoising
3 總結
為了分析該方法的除噪性能，這里將本論文的除噪方法與其它除噪方法如小波變換、數學形態學(mathematical morphology,MM)和傅里葉變換(Fourier transform,FT)相比較，評價除噪效果采用SNR、均方根誤差(root-mean-square error,RMSE)來表示。其中，SNR是除噪后的信號與噪聲的功率之比，單位為dB，SNR越大說明除噪效果越好；公式(7)中RMSEnoise表示噪聲的均方根誤差，Original(n)為原始采集信號，Denoise(n)為去噪后的信號，RMSEnoise越小，除噪方法越好。
| ${\rm{SNR}} = 10 * \lg \left( {{P_{{\rm{signal}}}}/{P_{{\rm{noise}}}}} \right),P = \frac{1}{N}\sum\limits_{n = 0}^{N - 1} {s{{\left( n \right)}^2}} $ |
| ${\rm{RMS}}{{\rm{E}}_{{\rm{noise}}}} = \sqrt {\frac{1}{N}{{\sum\limits_n {\left[ {{\rm{Original}}\left( n \right) - {\rm{Denoise}}\left( n \right)} \right]} }^2}} $ |
如圖 6所示，在SNR改善方面，EMD+軟閾值法提升SNR較為明顯，最大可提升10 dB，傅里葉變換法最差，數學形態學法和小波變換效果相仿；在RMSE比較中，EMD+軟閾值法除噪后的RMSE最小。
圖6
除噪效果比較
Figure6.
Comparison of denoising result
本論文的除噪方法，使用一次流程，不僅矯正了基線漂移，還除去了噪聲，對信號SNR改善明顯，除噪信號的RMSE也很小，所得除噪信號極大地保留了目標信號的有用信息，并且EMD和類小波軟閾值自適應算法容易實現，對同類信號均可實現類似的處理，在實際使用中具有重要價值和意義。
4 討論
論文中針對EMD分解信號時需要注意的三個問題采用了合理簡便的處理，從圖 6的結果可以看出，這些處理很好地解決了除噪問題，取得了較好的效果。
胎心音檢測技術可以得到胎心率等多種胎兒生理參數，是一種價廉簡便的胎兒健康監護技術，有光明的應用前景。目前，胎心音檢測技術的研究重點之一是采集的胎心音除噪。相比于其它方法，本論文的除噪方法和流程一次性地去除了采集信號中的噪聲，保留了信號中的有用成分。本研究為胎心音的除噪提供了新的方法，也為胎心音檢測技術的應用起到了推動作用。同時，該方法不僅適用于胎心音處理，對其它信號處理研究也具有借鑒意義，例如語音、雷達等信號的除噪。
本論文的實驗數據較少，未來的主要工作一是選擇醫院臨床實際采集的大量胎心音數據，用于算法的驗證和分析；二是在采集信號除噪的基礎上，計算胎心率等生理信息，并建立實際可用的系統，實現胎兒健康監護。
引言
在胎兒心臟中，心肌本身的舒張和瓣膜的關閉以及血流沖擊所產生的聲音叫做胎心音(fetal heart sound,FHS)[1]。胎心音包含較響的第一心音(S1)和第二心音(S2)，以及很微弱的第三心音(S3)和第四心音(S4)。Nagel[2]研究發現，胎心音采集信號主要含有母體心音、內部噪聲和胎心音S1與S2，其中內部噪聲和母體心音的頻譜主要在20 Hz以下，S1在23、34和75 Hz有三個譜峰，S2在30 Hz處有一個譜峰，并且噪聲頻段與胎心音頻段重合。同時，他還通過統計研究發現，隨著懷孕周數增加，S2的頻率逐漸降低，到一定階段趨于穩定，而S1的頻率先上升，在34周時達到最大，然后又逐漸下降，最后趨于穩定，即S1與S2的頻率與懷孕周數有較明顯的復雜函數關系。這都說明胎心音采集信號是一種非線性非平穩信號，要從中提取出有用的胎心音信號有很大困難[3]。
目前，因為胎心音檢測技術比超聲、超聲多普勒和胎心電等技術更安全，也更容易實現，故而受到越來越多的關注，但是在母體腹部采集得到的信號中含有大量噪聲，胎心音主要成分都淹沒在噪聲中，所以胎心音信號除噪方法是胎心音檢測技術的一個研究重點。周禮杲等[4]分析了胎心音的頻譜特點，提出了使用頻率濾波方法分離S1和S2，通過帶通濾波和匹配濾波等方法實現胎心音除噪，但是胎心音信號非線性非平穩的特點，即不同的胎兒在不同的時間或者環境，胎心音頻率均會發生變化，使該方法難以普遍適用。Varady等[5]使用帶通濾波器、小波變換(wavelet transform，WT)和互相關等方法做信號處理。 Chourasia等[6]使用小波工具箱和自適應濾波等方法對采集信號除噪。小波變換的本質是基于傅里葉的，對非平穩信號處理能力有限，并且它存在小波母函數選擇困難等問題。Ruffo等[7]使用匹配濾波、帶通濾波、Teager能量算子(Teager energy operator,TEO)和自相關等方法做信號處理和計算胎心率，該方法雖然較容易實現，但是在使用TEO能量算子時對處理信號的信噪比(signal noise ratio，SNR)要求較高。希爾伯特黃變換(Hilbert-Huang transform,HHT)是一種新的非線性非平穩信號處理方法，它在分解信號時能得到復雜信號的瞬時頻率，這與小波變換不同，而且它在地震、雷達、振動和成人心音等信號處理中得到了廣泛的應用，并取得了良好的效果，但是目前還沒有人將該方法應用于胎心音處理研究，所以本論文在比較和分析了該方法的理論和優點后，將HHT的核心——經驗模態分解(empirical mode decomposition method,EMD)應用于胎心音信號除噪。
1 經驗模態分解
由文獻[8]，HHT包含經驗模態分解EMD和希爾伯特變換(Hilbert transform,HT)，其中EMD能夠對復雜信號進行多層次分解，獲得若干代表瞬時頻率分量的本征模態函數(intrinsic mode function,IMF)和一個殘余分量。EMD分解信號時[9]，可以依據IMF的兩個條件[8]設置分解條件和停止規則，使得每個IMF只包含一個基本函數模式的振蕩波形。所以，任一復雜信號被EMD分解后，均能得到多個局部表征信號特點的IMF，這時可以依據需要，去掉不需要的IMF，組合重構剩下的IMF，即可實現信號除噪的目的[10]。
EMD在分解信號時，存在三個問題——端點效應[11]、停止規則和模態混疊[12]，需要采用合適的方法處理，否則信號分解可能不充分或者不準確。針對本論文，使用EMD分解胎心音采集信號時，端點效應問題是難以確定信號端點是否是極值點，論文中使用鏡像延拓法，在端點附近尋找合適的位置鏡像延拓一段信號，這樣在信號極值點包絡插值擬合時，就可以避免判斷端點是否是極值點的問題；EMD分解信號是為了得到單一分量，所以需要按照IMF的定義設定停止規則，同時，在實際使用中，為了提高信號處理效率，沒必要將信號徹底分解，這時需要對停止規則做些改變；模態混疊問題是因為噪聲頻率與胎心音有用信號頻率有重疊，可能造成信號分解不充分，使得某些分量中還含有其它頻率成分的信號，造成后期信號除噪困難，論文將采用加入掩模信號的方法消除模態混疊。雖然EMD分解信號需要注意一些問題，但是它具有小波變換多分辨率的優勢，同時克服了小波變換中小波基選擇的缺陷，對胎心音等非線性非平穩信號的處理更有優勢[13]。
2 胎心音除噪方法
本論文的胎心音采集信號均是采用參考文獻[14]的裝置和方法采集，使用一個高靈敏度的3M改裝聽診器、NI myDAQ數據采集卡和LabVIEW(Laboratory Virtual Instrumentation Engineering Workbench)控制軟件組成采集系統，實現了在母體腹部采集信號并以WAV格式保存原始音頻信號。實驗對象是在超聲醫療國家工程研究中心的員工及員工家屬內征集的10個孕婦志愿者(孕期從25周到35周)，受試者本人與家屬均簽訂知情同意書。按照文獻[14]的裝置和方法共采集200例數據，經過篩選，最后得到150例有用的胎心音數據。
該方法采集到的胎心音信號含有大量噪聲，包括采集電路引入的基線漂移噪聲，為了除去這些噪聲，如圖 1所示，主除噪部分采用EMD對信號進行模態分解，在分解時考慮消除模態混疊和端點效應的影響，并設置合適的停止規則，接著將EMD分解出的含有有用信號信息的分量采用類小波軟閾值自適應除噪，最后將除噪后的選擇單分量重構組合，即可得到除噪信號。
圖1
胎心音信號處理流程圖
Figure1.
Processing flow chart of FHS
2.1 預處理
為了節省資源，提高后期軟件程序的執行效率，針對胎心音信號頻段主要在20~200 Hz的特點，預處理中先采用截止頻率為200 Hz的巴特沃斯低通濾波器，然后將信號進行重采樣，以降低數據量。
胎心音采集信號采樣率為8 000 Hz，時長為5 s，圖 2是原胎心音采集信號經低通濾波和分頻重采樣后，得到信號采樣率2 000 Hz，并且由圖中可以看到相同時間長度上，信號的點數有所減少。因為胎心音的主頻段在200 Hz以內，遠遠小于重采樣后信號頻率2 000 Hz，所以重采樣后的信號依然包含胎心音信號的全部信息，頻譜沒有變化，而且數據量減少后，程序運行效率得到了提高。
圖2
信號預處理
Figure2.
Signal preprocessing
2.2 經驗模態分解
使用EMD對信號進行分解得到多個IMF和一個殘余分量的過程是一個篩選過程，篩選時要保證信號盡可能地分解為單分量信號，即分解的單分量信號不僅需要滿足IMF的條件，還需要消除EMD存在的三個問題。本研究將含噪信號按照以下步驟進行分解：
(1)使用EMD分解之前考慮端點效應影響，本論文使用鏡像延拓法處理信號端點效應。
(2)為了抑制模態混疊，在采集信號x分解時引入一個正弦掩模信號y，如公式(1)，它的頻率一般大于x的最高頻率，因為EMD分解時，模態混疊主要是由于頻率重合引起的，并且主要在高頻段，當加入一個高頻信號后，即可在EMD的低階分解時將高頻和低頻分開，從而達到消除模態混疊的作用，因為胎心音信號有用頻段在200 Hz以下，所以這里選擇正弦掩模信號的頻率為200 Hz，如圖 3為所選擇的正弦掩模信號波形圖。
| $y = mean\left( x \right) * \sin \left( {2\pi ft} \right)$ |
圖3
掩模信號
Figure3.
Masking signal
(3)加入掩模信號后，分解方法為：
| $\begin{array}{l} imf1 = emd\left( {x + y} \right)\\ imf2 = emd\left( {x - y} \right)\\ imf = \left( {imf1 + imf2} \right)/2 \end{array}$ |
上式中，imf1、imf2和imf均為單分量信號組成的矩陣。并且，使用x+y和x-y做EMD分解后求均值，這樣處理的原因是如果信號x不存在模態混疊，加入的y相當于一個噪聲，經過這樣的處理，可以消除y的影響，雖然這樣處理的弊端是增加了計算量，但是這種方法能夠最大限度地降低模態混疊的影響，在資源允許的情況下，是可以接受的。
(4)EMD分解時還需要設定合適的停止規則，按照IMF的條件設定判定規則，按照EMD分解結果設定算法結束規則，參考Rilling等[15]的研究，本論文設定的停止規則為：
①待處理數據的極值點個數不能少于3個，這是判定EMD結束的規則之一。
②假設三個變量T1、T2和α，且T1<T2，三個變量取值均在[0, 1]內，令m(t)、emax(t)和emin(t)分別表示待處理信號的包絡均值和上、下包絡振幅。然后要求滿足關系：
A、x(t)分解中的分量中任意一點值都小于T2，這是為了保證不出現局部過于不對稱的點。
B、x(t)分解中的分量的點超過T1的比例不大于α，這是為了保證對稱性，防止EMD分解時篩選次數太多，降低算法效率。
C、極值點和過零點數目相差不大于1，這是簡單中值規則，也是EMD分解結束判據之一。
經過前人大量的實驗和分析，并反復驗證掩模信號和停止規則，一般選擇T1=0.05，T2=0.5，α=0.05，可得到信號所有IMF分量和殘余分量。
使用該分解過程，不僅能很好地抑制模態混疊的影響，而且保證了信號的恰當分解，對后面的信號處理有重要意義。圖 4是預處理后的信號的EMD分解圖，信號被分解為9層，得到8個IMF分量和一個剩余分量，其中目標信號主要在IMF2、IMF3和IMF4分量中，而噪聲主要在IMF1以及IMF5之后的分量中。
圖4
EMD結果圖
Figure4.
Result of EMD
2.3 類小波軟閾值自適應算法
當含噪信號經過EMD分解得到多個單分量后,閾值函數的選擇以及閾值選取是各個單分量除噪的兩個關鍵問題。軟閾值函數是將信號絕對值小于閾值的數置為零，對大于閾值的數進行閾值收縮處理，適合用于心音的處理。軟閾值函數表達式為：
| ${\eta _T}\left( {\rm{x}} \right) = \left\{ \begin{array}{l} {\rm{sign}}\left( x \right) * \left( {\left| x \right| - T} \right),\left| x \right| \ge T\\ 0,\left| x \right| \le T \end{array} \right.$ |
其中T為閾值，ηT(x)是除噪函數，sign(x)是符號函數，x為輸入信號。
在閾值選擇時，若閾值太小，則可能除噪后信號中依然有大量噪聲；若閾值太大，除噪后信號的方差減小造成大偏差，有可能使信號過于平滑。考慮到EMD分解信號的特點，閾值采用基于高斯白噪聲的自適應閾值，閾值函數公式為[16]：
| ${T_k} = c * \sqrt {{E_k} * 2 * \ln N} ,k = 1,2, \cdots $ |
其中c表示閾值系數，一般取值范圍為(0.5，0.7)，經過大量實驗分析，這里選擇c=0.6，N表示信號長度，Ek(k=1,2,…)表示第k個IMF的噪聲強度，其計算公式為[17]：
| ${E_k} = \left\{ \begin{array}{l} \sigma _1^2,k = 1\\ \frac{{\sigma _1^2}}{{0.719}} * {2.01^{ - k}},k = 2,3, \cdots \end{array} \right.$ |
其中σ12是第一個IMF分量的方差。
如圖 4所示，信號經過EMD分解后，IMF1分量主要是噪聲，而后面幾個分量和剩余分量是低頻部分或趨勢信號，可以直接刪除，目標信號主要在IMF2、IMF3和IMF4分量中，因此在使用類小波軟閾值自適應函數除噪時，沒必要將所有分量都進行除噪，僅需要將含有目標信號主要成分的分量選擇出來，并對這些分量分別進行自適應除噪，然后將除噪后的分量組合重構，便可得到除噪信號。
如圖 5所示，對IMF2、IMF3和IMF4分量進行除噪，進行重構組合后，得到的干凈信號相比于原始信號，不僅消除了基線漂移的干擾，也很好地濾除了噪聲。
圖5
除噪結果圖
Figure5.
Result of denoising
3 總結
為了分析該方法的除噪性能，這里將本論文的除噪方法與其它除噪方法如小波變換、數學形態學(mathematical morphology,MM)和傅里葉變換(Fourier transform,FT)相比較，評價除噪效果采用SNR、均方根誤差(root-mean-square error,RMSE)來表示。其中，SNR是除噪后的信號與噪聲的功率之比，單位為dB，SNR越大說明除噪效果越好；公式(7)中RMSEnoise表示噪聲的均方根誤差，Original(n)為原始采集信號，Denoise(n)為去噪后的信號，RMSEnoise越小，除噪方法越好。
| ${\rm{SNR}} = 10 * \lg \left( {{P_{{\rm{signal}}}}/{P_{{\rm{noise}}}}} \right),P = \frac{1}{N}\sum\limits_{n = 0}^{N - 1} {s{{\left( n \right)}^2}} $ |
| ${\rm{RMS}}{{\rm{E}}_{{\rm{noise}}}} = \sqrt {\frac{1}{N}{{\sum\limits_n {\left[ {{\rm{Original}}\left( n \right) - {\rm{Denoise}}\left( n \right)} \right]} }^2}} $ |
如圖 6所示，在SNR改善方面，EMD+軟閾值法提升SNR較為明顯，最大可提升10 dB，傅里葉變換法最差，數學形態學法和小波變換效果相仿；在RMSE比較中，EMD+軟閾值法除噪后的RMSE最小。
圖6
除噪效果比較
Figure6.
Comparison of denoising result
本論文的除噪方法，使用一次流程，不僅矯正了基線漂移，還除去了噪聲，對信號SNR改善明顯，除噪信號的RMSE也很小，所得除噪信號極大地保留了目標信號的有用信息，并且EMD和類小波軟閾值自適應算法容易實現，對同類信號均可實現類似的處理，在實際使用中具有重要價值和意義。
4 討論
論文中針對EMD分解信號時需要注意的三個問題采用了合理簡便的處理，從圖 6的結果可以看出，這些處理很好地解決了除噪問題，取得了較好的效果。
胎心音檢測技術可以得到胎心率等多種胎兒生理參數，是一種價廉簡便的胎兒健康監護技術，有光明的應用前景。目前，胎心音檢測技術的研究重點之一是采集的胎心音除噪。相比于其它方法，本論文的除噪方法和流程一次性地去除了采集信號中的噪聲，保留了信號中的有用成分。本研究為胎心音的除噪提供了新的方法，也為胎心音檢測技術的應用起到了推動作用。同時，該方法不僅適用于胎心音處理，對其它信號處理研究也具有借鑒意義，例如語音、雷達等信號的除噪。
本論文的實驗數據較少，未來的主要工作一是選擇醫院臨床實際采集的大量胎心音數據，用于算法的驗證和分析；二是在采集信號除噪的基礎上，計算胎心率等生理信息，并建立實際可用的系統，實現胎兒健康監護。
