姜 星 耿讀艷,* 張園園 付志剛
1(河北工業大學省部共建電工裝備可靠性與智能化國家重點實驗室,天津 300130)2(河北工業大學河北省電磁場與電器可靠性重點實驗室,天津 300130)3(中國人民解放軍第254醫院,天津 300142)
心沖擊圖(ballistocardiogram,BCG)可以展示心臟搏動和血液在大動脈流動而引起的人體對外壓力或體表位移的變化,反映心臟的力學特征[1],通過壓力、加速度、電容等傳感器來感知身體的這種微小位移變化,從而獲取心臟搏動引起的體震信息。正常的BCG信號波形及其生理意義如圖1所示,主要包含H、I、J、K、L、M和N等7種波形,其中I、J、K反映心動周期的射血過程,而L、M、N反映舒張過程[2]。相比Holter等心電監測方式,BCG可實現更長期的連續監測,避免電極及其導線束縛造成對人情緒的影響,尤其適合于長時間作息和睡眠監測,比較適用于家庭應用。然而,BCG信號自身十分微弱,且易受到呼吸、體動、工頻噪聲的干擾[3-4],導致直接測量得到的BCG信號含有太多噪聲,不能直接獲得準確的生理信息特征。因此,BCG的降噪處理顯得尤為必要。

圖1 正常BCG信號波形及其生理意義Fig.1 Normal BCG signal and physiological significance
目前對于BCG信號的降噪處理主要是采用小波變換方法,一般是對周期性平穩信號有效,而對包含尖峰或突變的非平穩信號效果不佳[5-6]。小波降噪能同時進行時頻分析,且能有效區分非平穩信號的尖峰與噪聲,但需選取合適的小波基才能達到較好降噪效果,而小波降噪仍基于Fourier變換,尚不能擺脫以Fourier變換為基礎帶來的缺陷[7-9]。經驗模態分解(empirical mode decomposition, EMD)是Huang在1998年提出的一種用于非線性和非平穩時間序列信號的處理方法[10],主要是將信號分解為若干個固有模態函數(intrinsic mode function, IMF)分量,使每個IMF包含不同的頻率成分,從而實現信號特征的有效分解。但是在分解過程中,由于BCG信號中有大量噪聲成分,存在較強的時頻耦合特征,會產生模態混疊現象[11],造成IMF成分中包含不同的時間尺度,從而對各個模態分量的物理意義判別造成困難。獨立分量分析(independent component analysis, ICA)無需信號源先驗知識,在時域和頻域都能將觀測得到的混合信號中的獨立成分分離出來。ICA去噪的關鍵點在于噪聲通道的構建,孫云蓮等利用EMD分解得到的IMF希爾伯特時頻譜來構建噪音通道,在信噪比較高的情況下,基于IMF時頻譜分析的ICA降噪法取得了不錯的效果[12]。本研究在其基礎上,提出一種基于EMD及ICA相結合的BCG信號降噪方法。該方法對含噪BCG信號進行EMD分解,得到頻率從高到低的多個IMF;通過引入模態相關法,找到噪聲與信號的分界IMF;對分界之前的分量構建虛擬通道[13],與原信號一起作為ICA的輸入,得到降噪后BCG信號。
在本研究中,使用的是自主研制的BCG信號采集系統,主要由采集人體信號的壓電薄膜式傳感器、連接傳感器對采集信號做后續處理的電路和上位機軟件組成。將傳感器置于座椅坐墊下,受試者靜止地坐在座椅上,心臟射血時引起身體微小振動,通過傳感器獲取射血引起的重力變化信號[14],經過電路的多級放大和轉換,將其無線傳輸到上位機軟件,最終得到測量的BCG信號。采集電路主要包括電荷放大器、低通濾波器、電壓放大器、電壓抬升電路、信號離散模塊、無線傳輸模塊和串口通信,其總體構成如圖2所示。

圖2 電路整體設計框圖Fig.2 Overall design circuit diagram
經驗模態分解的實質是對一個時間序列信號進行平穩化處理,其結果是將信號中不同尺度的波動或趨勢逐級分解,產生一系列具有不同特征尺度的數據序列,每一個序列稱為IMF。分解得到的IMF必須滿足兩個條件:一是IMF中極值點的個數與過零點的個數相等或不超過1個;二是由極大值與極小值確定的包絡線均值為0。
獨立分量分析是由盲源分離計算發展起來的一種信號處理方法,要求輸入的是多維信號,并且混合信號個數要大于等于源信號個數,才能有效地分離。
由獨立分量分析的原理可知,獨立分量分析要求輸入的是多維信號,并且混合信號個數要大于等于源信號個數,才能有效地分離。基于此,BCG信號經EMD分解后產生若干個基本信號,即
(1)
式中,x(t)表示原始信號,imfi表示第i個IMF分量,rn(t)表示余項。
提取出噪聲所在的IMF分量進行組合,得到噪聲的估計量;與源信號一起作為ICA多維噪聲輸入,構建ICA噪聲通道。本研究采用基于EMD模態相關分選準則[15]判定噪聲與信號的分界層,計算模態函數與原始BCG信號之間的互相關系數,有
(2)



圖3 EMD-ICA降噪原理Fig.3 EMD-ICA de-noising principle
本研究為了驗證所得到的BCG信號的準確性,使用ECG-300 A數字式三道心電圖機和所設計的BCG采集系統,同步采集被試者的ECG和BCG信號,并實時保存。對所采集的BCG信號分析得到心率,并與心電圖機所得心率進行對比,分析BCG的準確性。
為了驗證裝置的有效性,并對人體心沖擊信號的特征進行研究,本研究共采集了10位健康成年人的數據,其中有5位男性、5位女性,年齡均在23~27歲。要求受試者保持平穩的坐姿,坐在放有壓電薄膜傳感器的座椅上,勻速呼吸,將傳感器的輸出端接在信號調理電路上,使用示波器觀察其BCG信號,等待信號穩定后采集2 min的BCG信號,同時采集被試者的ECG信號。將采集到的數字信號經過Matlab處理后得BCG的原始信號,圖4為某一被試者采集的原始BCG信號。

圖4 原始BCG信號Fig.4 Example of the original BCG signal
由圖4可見,BCG信號中可以明顯觀察到J峰,同時可以看到BCG信號的I、J、K波形是有重復性的,所以表示所采集的BCG信號是具有魯棒性的。另外,不同的被試者BCG信號的波形形狀和幅度存在著較大的差別,這些差別來源于每個人不同的心臟活動狀態以及身體結構。
將采集到的BCG信號進行EMD分解,結果如圖5所示。

圖5 EMD分解結果Fig.5 EMD decomposition results
由分解結果可以看出,原始信號經EMD分解得到了幾個固有模態分量和一個趨勢項。由EMD分解過程可知,分解得到的13個IMF分量imf1,imf2,…,imf13中,噪聲對每個IMF分量的支配作用逐漸降低,信號對IMF的支配作用不斷加強,因此需要確定IMF分量中信號與噪聲的分界[16]。根據式(1)計算出各階固有模態分量與原始BCG信號的相關系數,得到結果后繪制相應的曲線,如圖6所示。

圖6 相關系數曲線Fig.6 The correlation coefficient curve table
按照前述信號與噪聲的分界點判定方法,由圖6可以看出,imf4為分界IMF。根據噪聲與信號的分界,選擇第1~3個IMF分量作為虛擬噪聲通道。由于原始BCG信號進行EMD分解后imf4之前的IMF分量大都是噪聲成分,所以借鑒信號時序平移思想構造ICA的多維輸入[17],其基本思想是:把原始BCG信號分解后得到的由第1~3個IMF分量構建的虛擬噪聲通道時序序列向左平移l個位置,向左溢出部分拼接到該時序序列的右端,這樣就可以得到新的含噪信號作為ICA的噪聲通道輸入。這個新的含噪信號與原信號相比信噪比幾乎一致,不僅保留了原信號中的有效成分,而且信號中含噪聲的成分也沒有大的改變,僅是噪聲形態發生了改變。
這樣就構造了多個輸入信號,包括兩個虛擬噪聲通道和原始BCG信號。在獨立分量分析中,常用的算法是Fast ICA和Infomax算法,從計算速度方面考慮,本研究采用一種基于負熵的判據和一種非常高效的FastI CA[18-19]算法進行分析,其框圖如圖7所示。負熵定義如下:

(3)
式中,E為均值運算,g為非線性函數,Ygauss是高斯變量。Fast ICA其步驟如下:
1) 對觀測信號進行中心化處理,使其均值為0。
2) 利用特征值分解法對中心化后的數據進行白化,得到標準化數據Z(t)。
3) 選擇要估計的獨立成分個數m,迭代次數為p,初始值為p=1。
4) 選擇一個隨機的初始分離矩陣W。
5) 迭代更新分離矩陣,令
(4)
(5)
Wp=Wp/‖Wp‖
(6)
式(4)中,g(y)=tanh(a1,y),a1=1。
6) 判斷Wp是否收斂,如果不收斂,返回步驟5)。

圖7 Fast ICA算法框圖Fig.7 Fast ICA algorithm block diagram
7) 令p=p+1,假如p≤m,返回步驟4),重復執行步驟4)~7),直到將所有需要的獨立分量提取出來為止,算法結束。
為了定量說明本方法的降噪效果,同時采用EMD方法和小波去噪方法[20],分別對本研究的含噪BCG信號進行降噪處理,比較不同方法的降噪效果。為了定量分析3種方法的降噪效果,計算其降噪后信號與原始信號的相關系數(R值)、信噪比(Rsn值)和能量比(Esn值)。相關系數值說明了降噪后的信號與原始信號之間整體波形的相似度,信噪比說明了信號的去噪情況,而能量比值反映了降噪后信號占原信號的能量百分比。其中,信噪比Rsn定義[21]為
(7)

若降噪后信號具有較高的信噪比,說明降噪效果好。
信號能量計算公式[16]如下:
(8)
式中,x(t)為信號幅值,T為采樣周期,m為采樣點數。
降噪后的信號占原信號的能量百分比Esn定義為
Esn=E0/E
(9)
式中,E表示原始信號能量,E0表示降噪后信號能量。
Esn越大,說明降噪后的信號越高,越能保持原信號特征,也越接近原信號。
為了驗證信號采集系統的準確性,經過心沖擊圖和心電圖方法計算得出心率[22],如表1所示。對其用配對t檢驗方法進行檢驗,結果顯示兩種方法計算得出的心率不具有顯著性差異(P>0.05),說明本研究的BCG采集系統具有一定的準確性。

表1 被試者的心率計算結果

圖8 BCG信號降噪前后對比。(a)降噪前BCG信號;(b)降噪后BCG信號Fig.8 BCG signal de-noising before and after the comparison. (a) Noise reduction before BCG signal;(b) Noise reduction after BCG signal
根據上述的分析方法,由Fast ICA最終得到降噪后的BCG信號。圖8所示為BCG信號降噪前后的波形,可以觀察到降噪處理后BCG的信號特征得到了明顯的保留,信號波形較之前也趨于平滑。為了進一步證明降噪效果,計算了降噪前后BCG信號的功率譜密度(power spectral density,PSD),如圖9所示。可以看出,在高頻部分BCG信號功率譜明顯趨于平滑,且在低頻部分沒有明顯的衰減,進一步說明本研究采用的方法在去除噪聲的同時,很好地保留了BCG信號的低頻特性。
利用算法評價指標,分別計算10名被試的原始BCG信號降噪前后的相關系數、信噪比、能量百分比,得到的結果如圖10所示。對其進行單因素ANOVA檢驗,結果顯示3種降噪方法計算得出的相關系數、信噪比、能量比具有顯著性差異(P<0.05);并計算其平均值,結果如表2所示。

表2 不同降噪算法效果對比Tab.2 Comparison of different noise reduction algorithms

圖9 降噪前后BCG信號的PSD。(a)降噪前;(b)降噪后Fig.9 PSD before and after de-noising BCG signal. (a) Before de-noising;(b) After de-noising

圖10 不同降噪方法處理前后的對比。(a)相關系數對比;(b)信噪比對比;(c)能量比對比Fig.10 Comparison of quantitative figure before and after different noise reduction methods. (a) Comparison of correlation coefficient;(b) Comparison of signal-to-noise ratio; (c) Comparison of energy percentage
由圖10可以看出,EMD-ICA方法得到的相關系數、信噪比和能量百分比的各項指標均優于小波和EMD方法的指標。如表2所示,對其進行配對t檢驗,EMD-ICA相比EMD方法信噪比顯著提高(5.19±1.29vs14.87±3.04,P<0.001),相比小波方法能量百分比顯著提高 (88.81±2.81vs96.64±2.92,P<0.001)。雖然小波降噪方法的信噪比較高,但是相關系數和能量百分比相對于其他兩種方法最低。EMD降噪方法雖然相關系數和能量百分比相對較高,但是信噪比相比其他方法最低。3種方法比較來看,小波降噪方法的降噪效果雖較好,但也損失了部分原信號的信號特征;EMD降噪方法能較好地保持原信號特征,接近原信號,但降噪效果不明顯;EMD-ICA降噪方法的降噪效果最明顯,最大限度地保留了原信號的信號特征。
消除噪聲是BCG信號預處理的重要內容,由于采集系統和外部環境等因素的影響,采集的BCG信號數據含有不同類別的噪聲,具有非線性、非平穩性,嚴重影響心率等生理信息測量的準確性。綜合考慮降噪效果和計算量,本方法避免了文獻[23]中需要大量實驗確定小波基與分解尺度的問題。從信號分解角度,小波變換中基函數是固定的,不能匹配信號的全部,用自適應EMD分解方法代替復雜的小波變化,大大降低了算法的復雜度。文獻[14]中EMD分解BCG信號確定IMF階數后,直接對各分量進行重構,不能有效消除IMF中由于高頻噪聲的干擾以及時間尺度信號變化等因素引起的模態混疊噪聲。本研究提出信號時序平移方法構造多個觀測信號,通過ICA可以實現混疊頻率的有效分離,從而估計出IMF中的原始信號,能有效消除EMD分解過程中模態混疊現象對降噪結果的影響,提高降噪效果。
EMD方法因在處理非線性、非平穩信號上的優越性而得到了廣泛的應用,能夠自適應地分解信號,但其分解過程中面臨的端點效應問題阻礙了其發展。對于該問題,目前提出的很多方法都是從算法上入手來抑制端點效應的影響,忽略了噪聲信號在EMD分解過程中對端點效應的影響。本研究選擇增加BCG信號長度,以減少端點效應對信號降噪的影響。此外,分離信號的幅度和順序上存在不確定性,但是BCG信號的特征存在于相對幅頻特征中,因此對BCG信號的研究影響不大。
由于壓電薄膜傳感器的測量原理,BCG信號極易受到身體晃動的干擾,咳嗽、扭動身體等引起的輕微運動都會影響信號的質量,使得特征難以辨認,不能有效地提取信號。本系統目前適用于安靜狀態下獲取BCG信號,比較適合睡眠時檢測,特別適用于家庭應用。因此,未來工作首先要解決BCG信號在各種條件下的采集與降噪問題,以保證BCG信號的質量。其次,隨著對BCG信號生理意義研究的不斷深入,希望探求一種既能代替心電圖實現多參數無感檢測又可以將BCG信號包含的心血管系統生理信息用于輔助醫生的方法,從而對各種心血管疾病進行診斷與分析。這是今后研究的方向,其實現將有力于促進人類健康水平的提高。
本研究綜合考慮EMD與ICA的技術優點,提出一種基于EMD-ICA的BCG信號降噪方法。該方法利用EMD將BCG信號分解,根據模態相關準則將分解信號進行信號層與噪聲層的判定,得到重構信號并作為ICA的輸入,最后去除了原始信號中的噪聲成分。該方法克服了EMD分解過程中產生的模態混疊問題,同時解決了ICA方法中多維輸入的要求,實現了多周期BCG信號的降噪。通過與小波和EMD降噪方法進行對比實驗研究,表明該方法能夠更有效地識別心臟的動力學信號,這對后期研究日常作息條件下與人體健康狀態相關的生理病理信息有重要意義。