郭興明 袁志會
(重慶大學生物工程學院,重慶 400044)
心音信號是一種典型的非平穩(wěn)、非線性信號。對于心音信號的研究,傳統(tǒng)的處理方法一般是在假設(shè)心音信號為平穩(wěn)或分段平穩(wěn)的前提下進行的,這就需要用虛假的高頻成分去補償信號的非平穩(wěn)性,容易引起能量的擴散,以致得到錯誤的信息。經(jīng)驗?zāi)J椒纸夥椒?empirical mode decomposition,EMD)是1998年由 HUANG等提出的[1],它是基于信號的局部特征時間尺度,把復(fù)雜的信號分解為固有模態(tài)函數(shù)(intrinsic mode functions,IMF)之和,每一個IMF所包含的頻率成分不僅與頻率有關(guān)還隨信號本身變化而變化。因此,該方法具有自適應(yīng)性,非常適合處理非線性與非平穩(wěn)信號?;谶@樣的特性,EMD也被廣泛用于心音信號特征信息的提取[2-3]。
在獲取IMF分量時,由于采用三次樣條插值的方法獲取信號上下包絡(luò)線的瞬時平均,這使得EMD分解過程中存在特殊的邊緣效應(yīng),即端點效應(yīng)問題,嚴重的會影響整個信號的分析,這是 EMD分解中的一個難點。目前常見的鏡像延拓法[4]、基于神經(jīng)網(wǎng)絡(luò)的延拓法[5]、基于多項式擬合的延拓方法等雖然在一定程度上抑制了端點效應(yīng)[6],但這些方法的本質(zhì)都是在邊界兩端預(yù)測數(shù)據(jù)來抑制端點效應(yīng)對信號分解的影響,而沒有考慮信號中噪聲對三次樣條插值的影響。由于信號在采集過程中容易疊加大量的噪聲信號,上述方法并不能從根本上解決端點效應(yīng)問題,如果噪聲信號在EMD分解過程中一直存在,就會增加大量的無用頻率分量,隨著分解層數(shù)的增加,邊界處的波動逐漸向內(nèi)“污染”整個數(shù)據(jù)序列,從而加劇端點效應(yīng)。文獻[2]在用EMD方法對心音的處理中只注意到了噪聲對產(chǎn)生端點效應(yīng)的影響,而忽略了樣條插值的影響,文獻[3]則和文獻[2]相反,只考慮了樣條插值的影響。因此只有同時考慮噪聲和樣條插值的影響才能更好的解決端點效應(yīng),提高EMD的分解質(zhì)量。
本研究在經(jīng)驗?zāi)J椒纸獾幕A(chǔ)上,結(jié)合近年來本科研團隊的相關(guān)研究工作提出一種新的心音信號分析方法[7-8],即首先采用基于小波的方法去除心音信號中的噪聲干擾,再用改進的EMD端點效應(yīng)的新方法即自適應(yīng)波形匹配延拓法處理心音信號,通過對得到的IMF分量進行Hilbert變換,為后續(xù)心音的定位和分析提供參考。實驗表明該方法減少了EMD的分解層數(shù)以及端點效應(yīng)對信號分析的影響,大大提高信號特征提取的精度和時效性。
在心音信號的采集過程中,總有噪音干擾伴隨信號出現(xiàn),心音信號中的噪聲主要來自環(huán)境噪聲、工頻噪聲、儀器本身的噪聲以及采集設(shè)備與對象皮膚的摩擦等[9]。在進行 EMD分解時,如果信號被噪聲污染,將會極大的影響EMD過程中對極大值和極小值的有效提取,使要取的包絡(luò)線出現(xiàn)不同程度的失真,同時在對原信號運用EMD進行波形匹配延拓時,尋找到的匹配波形必然也夾雜著噪聲干擾,端點處的延拓同樣受到噪聲的干擾,影響IMF分量的精度。因此,先利用小波去處噪聲對心音信號的干擾,減少端點效應(yīng)產(chǎn)生的可能性,這樣可以提高IMF的精度,減少EMD分解過程的分解層數(shù)。
在充分考慮小波函數(shù)的性質(zhì)和心音信號特點的基礎(chǔ)上,選用 coif3小波函數(shù)rigrsure閾值處理方法對心音信號的消噪,去噪結(jié)果如圖1所示,(a)為原始心音信號,(b)為小波去除噪聲后的信號,從圖中可以看出消噪后的信號較好的保留心音中的主要成分,減少毛刺等高頻信號對波形的干擾,效果較好,從而減少噪聲對IMF分量的影響。
圖1 心音小波去噪后的圖形。(a)原始信號;(b)去噪后的信號Fig.1 Wavelet denoising diagram.(a)Original heart sound signal;(b)Denoised heart sound signal
EMD方法基于這樣的假設(shè)[1]:任何信號都是由一系列幅度和相位都隨時間變化的基本模式分量構(gòu)成,它必須滿足兩個條件,一是它的零點數(shù)與極點數(shù)相等或至多相差1個;二是由它的極大值和極小值確定的上下包絡(luò)線均值為0。EMD將這種基本模式分量定義為固有模態(tài)函數(shù),即IMF。EMD方法的一個重要步驟是構(gòu)造信號的上下包絡(luò)線,得到信號的瞬時平均,Huang等提出以信號的極大值點和極小值點擬合三次樣條曲線的方式來構(gòu)造上下包絡(luò)線,由于信號在端點處往往并非極值點,因此三次樣條曲線容易在數(shù)據(jù)序列的兩端出現(xiàn)發(fā)散現(xiàn)象,并且這種發(fā)散的結(jié)果會伴隨著各個 IMF的篩選過程,逐漸向內(nèi)“污染”整個數(shù)據(jù)序列而使所得結(jié)果嚴重失真,這就是EMD分解的端點效應(yīng)。
文中第一部分論述了噪聲對EMD分解過程中端點效應(yīng)的影響,消除噪聲的干擾從一定程度上能抑制端點效應(yīng),但端點效應(yīng)產(chǎn)生的本質(zhì)原因是對信號進行三次樣條插值時,由于信號在端點處并非極值點而造成包絡(luò)線在信號兩端發(fā)散產(chǎn)生的,因此,要想更好的抑制端點效應(yīng),僅僅考慮消除噪聲干擾還是不夠的,還要充分考慮信號的內(nèi)部特性。
端點延拓是解決EMD端點效應(yīng)問題的有效手段。目前,普遍采用的抑制端點效應(yīng)的方法有鏡像延拓法[4]、基于神經(jīng)網(wǎng)絡(luò)的延拓法[5]、基于多項式擬合的延拓方法等[6]。這些方法對端點效應(yīng)的抑制都有一定的效果,卻又都存在著各自的問題:鏡像延拓方法因為它可能要截去部分數(shù)據(jù),在處理短數(shù)據(jù)時效果比較差;對于神經(jīng)網(wǎng)絡(luò)方法而言,處理速度比較慢;多項式擬合方法適應(yīng)性較差。本研究充分考慮信號自身的內(nèi)在相似特性,采用一種基于波形匹配的自適應(yīng)端點延拓方法。
延拓出的波形盡可能地符合原信號的變化趨勢,這樣的延拓才有意義。如果在原始信號內(nèi)部尋找一段與邊緣處變化趨勢非常相似的子波,在每次提取IMF前,用這段子波來延拓信號邊緣處的波形,使原端點的延展最大限度的維持原來的變化趨勢?;诖艘胱赃m應(yīng)波形配備延拓法[11]:從原始信號內(nèi)部找出最符合信號邊緣處變化趨勢的波形對信號進行延拓,最大限度地維護信號的內(nèi)在變化趨勢,而對內(nèi)在規(guī)律較弱、邊界數(shù)據(jù)變化異常的信號可以只考慮信號邊緣的極值點信息,具體通過設(shè)定一閾值來實現(xiàn),這樣使該方法具有自適應(yīng)性。
信號的延拓包括左右兩端,下面以左端的延拓為例來說明該方法,設(shè)原始信號為x(t)、pi、qi(i=1、2、3、…)分別為信號 x(t)的極大值和極小值,其對應(yīng)的時間分別為 tpi、tqi。設(shè)信號的左端數(shù)據(jù)為 x(1),p1,q1為左端起信號x(t)的第一個最大值和最小值,定義x(1)、p1、q1三點構(gòu)成的三角波形為特征波形,沿著x(t)尋找與特征波形最為相似的三角波形x(i)、pi、qi為匹配波形,將匹配波形前的數(shù)據(jù)作為x(t)的延拓波形,這樣即得到符合信號變化得自然趨勢,具體步驟如下:
步驟1:找出除特征波形外,所有匹配三角波形的起始點值x(i),其對應(yīng)的時間點定義為
步驟2:計算所有匹配三角波形與特征波形的匹配誤差,誤差公式為
步驟3:找出最小的匹配誤差值 mine(i),并設(shè)定一閾值α,如果mine(i)<α,則將 mine(i)對應(yīng)的三角波形作為匹配波形,將匹配波形前的數(shù)據(jù)延拓到原始信號的前面,如果 mine(i)≥α,跳轉(zhuǎn)至4。
步驟4:直接設(shè)置信號端點處的極值,即分別求出最靠近信號端的相鄰的N個極大值點的平均值和M個極小值點的平均值,將其分別作為信號x(t)的極大值和極小值,具體N和M的大小可根據(jù)信號在邊界處的跳變情況而定。
右端點也按上述方法延拓,將延拓后的信號做EMD分解,對分解得到的 IMF分量按原信號x(t)對應(yīng)的時刻進行截取,這樣便能獲得改善了端點效應(yīng)的分解結(jié)果。上述方法不僅充分考慮了信號的內(nèi)在規(guī)律與變化趨勢,減小了分解結(jié)果在端點處的振蕩,并且對內(nèi)部規(guī)律較弱或邊界數(shù)據(jù)異常的信號,也能進行較好的處理,與通常的端點處理算法比較,有更強的自適應(yīng)性。
為了驗證所提出算法的可靠性和有效性,對40例心音(正常心音10例,異常心音30例,每例的記錄時間為20s,采樣頻率為11025Hz)進行 Hilbert變換提取包絡(luò)后用雙閾值分割法實現(xiàn)心音的定位[8]。由于文獻[8]中通過舍棄信號前后的數(shù)據(jù)點來減少端點效應(yīng),但本算法不需要舍棄信號前后的數(shù)據(jù)就能很好的抑制端點效應(yīng),而此對雙閾值分割定位策略做如下修正。
1)S1持續(xù)時間約為100~110 ms,S2持續(xù)時間約為70~80 ms,用一個包含110個數(shù)據(jù)點(時限為150 ms左右)的窗,以平移的方式對數(shù)據(jù)段進行處理,這樣保證完整的S1或S2出現(xiàn)在窗內(nèi)。找出每個窗中的最大值,依次存入到數(shù)組M當中,求得平均值Ma,舍棄M中小于Ma的值,再對剩余的數(shù)據(jù)求平均值Mb。利用第一次平均值來剔除基線過大擾動的影響能取得較好的效果。令H=0.3Mb,L=0.05Mb,確定雙閾值的高閾值H和低閾值L。
2)平移窗,以 H為臨界值,找出每次平移后窗中第一個大于等于該值的點,如圖2中的 A點,把其數(shù)據(jù)點位置存入到數(shù)組A中。在H臨界值下,不僅保證夠檢測出幅值較大的有用信號,同時能夠去除幅值較小的噪音及雜音的干擾。
3)依次讀取A數(shù)組中數(shù)據(jù)點的位置,找出向前搜索和向后搜索時各自的第一個小于等于L的數(shù)據(jù)點,如圖2中的B和C點,存入數(shù)組B和數(shù)組C。
4)依次從數(shù)據(jù)組B和數(shù)據(jù)組C中取出一個數(shù),組成一對數(shù)據(jù),用C中的數(shù)據(jù)減去B中的數(shù)據(jù),得到的數(shù)據(jù)存入數(shù)組D中,則D中存儲的是檢測到的每個包絡(luò)的寬度(可以轉(zhuǎn)化為時限)。
5)S1和S2交替出現(xiàn)的,所以取出D中的奇數(shù)項和偶數(shù)項,分別存入數(shù)組E和數(shù)組F中,對E和F中的數(shù)值求平均值并比較,因為S1到S2的間隔都要小于 S2到下一個S1的間隔,判定該間隔的起點為 S2,終點為S1,至此,心音分段成功。
圖2 雙閾值分割Fig.2 Segmentation with two thresholds
圖3是原始心音信號的基于波形匹配延拓的EMD結(jié)果,其中第一層為未去噪的心音信號。圖4是去噪后心音信號的基于波形匹配延拓的EMD結(jié)果,第一層為去噪后的心音信號。采用自適應(yīng)波形配備延拓法來從算法本身改善EMD的端點效應(yīng)問題,從圖中可以看出,端點效應(yīng)得到很好的抑制,但是圖3中由于噪聲的影響,不必要的干擾信號參與分解,將每層分解的包絡(luò)線的誤差進行積累,傳遞到后面的分解層,增加分解的層數(shù),從而影響 IMF分量的精度,圖4中由于消除了噪聲的干擾,減少了無用的頻率分量,第一層和第二層IMF已很好地將心音信號分離出來。
對比圖3和圖4可以看出經(jīng)過小波去噪再用改進端點效應(yīng)的 EMD方法分解心音信號,不僅可以減少噪聲對邊緣效應(yīng)的影響,還可以有效地減少累計誤差,減少計算量。從仿真結(jié)果可以看到效果是比較理想的,去噪后心音信號的 EMD分解得到的IMF分量明顯好于直接對信號進行EMD分解得到的IMF分量。隨著IMF分量階數(shù)的增加,其分解出信號的頻率隨之降低,低頻段的IMF分量所含的心音信號的有效信息越來越少,通過圖4觀察各階IMF分量,第一和第二階 IMF分量已經(jīng)能很好地表達原信號,因此在后續(xù)的分析中可用這兩階信號之和代替原信號。圖5中(a)是進行小波去噪后的心音信號,(b)是用IMF1和IMF2之和代替原信號,從圖中對比可以發(fā)現(xiàn),通過小波去噪后信號仍存在不少毛刺和不規(guī)整之處,而EMD分解過后,信號變得更加規(guī)則,研究中用這兩階IMF之和來替代原始信號,利于心音信號的后續(xù)處理。
圖3 原始心音信號延拓后的EMD結(jié)果Fig.3 The EMD result of original heart sound signal after extension
圖4 去噪心音信號延拓后的EMD結(jié)果Fig.4 The EMD result of denoising heart sound signal after extension
經(jīng)過本算法對心音信號進行EMD分解重構(gòu)后,用Hilbert變換來提取心音信號的包絡(luò)并用雙閾值分割法對心音進行定位分段。心音信號的包絡(luò)定位圖如圖6所示,從圖中可以看出已經(jīng)很好的把S1,S2定位出來。定位之后對心音信號進行分段,根據(jù)定位信息確定出S1,S2的邊界,如圖7所示。
圖5 心音信號EMD分解重構(gòu)結(jié)果。(a)小波去噪后的信號;(b)重構(gòu)后的信號Fig.5 TheEMD decomposition reconstruction result of heart sound signal.(a)Wavelet denoising diagram;(b)Signal of the reconstruction
圖6 心音信號的包絡(luò)與定位Fig.6 The envelope and localization ofheart sound signal
圖7 分段后心音信號(“o”表示第一心音;“*”表示第二心音)Fig.7 Segmentation of heart sound signal(“o”marks the first heart sound;“*”marks the second heart sound)
對40例心音定位結(jié)果如表1所示,S1的檢出率為97.05%,S2的檢出率為97.12%,與文獻[8]中的結(jié)果相比檢出率有所提高。在準確定位的基礎(chǔ)上,對以上40例心音信號共計1080個心動周期進行分段,其分段結(jié)果如表2所示,從表2中可以看出,該算法對正常心音有較好的分段效果,準確率達到97.03%,異常心音的分段準確率為86.54%與文獻[11-12]中所用方法結(jié)果相符。
表1 定位正確率Tab.1 Correct rate of location
表2 分段正確率Tab.2 Correct rate of segmentation
文獻[8]中對第一心音和第二心音的識別率不太高的主要原因是有些噪聲的幅值高于了鄰近的S1或S2,從而產(chǎn)生了誤判。本研究采用小波變換結(jié)合改進后的EMD算法處理的心音信號減少了信號的不規(guī)整之處和噪聲的影響,選取第一階和第二階固有模態(tài)函數(shù)之和替代原始心音信號,減少了EMD的分解層數(shù),提高了 IMF分量的精度,S1和 S2的檢出率得到提高,為心音信號的準確分段提供了參考。在準確定位的基礎(chǔ)上對心音信號進行分段,結(jié)果表明該算法對正常心音有良好的分段效果,對異常心音的分段也能達到較好的效果,進一步分析表明,由于異常心音中出現(xiàn)寬分裂,以及連續(xù)性的高幅度心雜音導(dǎo)致分段效果不是十分理想,如何對上述心音進行準確的分段是下一步要研究的內(nèi)容。
由于不同類型的信號所需要的小波基不同,且由于EMD本身也具有去噪性,自適應(yīng)波形匹配延拓法中,閾值α也要根據(jù)具體信號而定,設(shè)定不當也會影響EMD分解的準確性,所以本研究提出的算法存在一定的局限性。在實際應(yīng)用中小波基的選擇與自適應(yīng)波形匹配中閾值α的選擇達到一個良好的統(tǒng)一,以此來提高IMF分量的精度。怎樣使小波基與自適應(yīng)波形匹配中閾值α達到良好的統(tǒng)一是下一步研究的重點,也是提高心音分段正確率的關(guān)鍵。
EMD方法自產(chǎn)生以來因其在處理非線性非平穩(wěn)信號上的優(yōu)越性,得到廣泛的應(yīng)用,但是 EMD分解過程中面臨的端點效應(yīng)問題阻礙了其發(fā)展。對于該問題目前提出的很多方法都是從算法上入手來抑制端點效應(yīng),忽略了噪聲信號在EMD分解過程中對端點效應(yīng)的影響。本研究從心音信號的噪聲特征和內(nèi)在特性出發(fā),提出先用小波去噪,去掉噪聲干擾,減少不必要的分解層帶來的邊緣效應(yīng)誤差,之后再進行 EMD分解,提高分解的精度。針對EMD的端點延拓,提出了一種自適應(yīng)波形匹配延拓法來解決,它從原理上減小了延拓信號端點包絡(luò)線的發(fā)散程度,有效地抑制了EMD分解的端點效應(yīng),提高IMF分量的準確性。在此基礎(chǔ)上,根據(jù)心音信號生理上的時限關(guān)系,提出了用于雙閾值分割的具體閾值,準確定位分段心音信號中S1和S2。
實驗結(jié)果表明本研究為心音的定位與分段提供準確的信息,為進一步的心音的分類識別以及輔助診斷奠定了良好基礎(chǔ)。
[1]Huang Wentao,Shen Zhen,Steven R.Long,et al.The empirica mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis[J].Proc R Soc Lond A,1998,454:903 - 995.
[2]許曉飛,林勇,嚴彬彬.基于希爾伯特-黃變換的心音包絡(luò)提取[J].航天醫(yī)學與醫(yī)學工程,2008,21(2):134-136.
[3]林勇,許曉飛.基于經(jīng)驗?zāi)J椒纸獾男囊糇詣臃侄嗡惴ǎ跩].中國生物醫(yī)學工程學報,2008,27(4):485-489.
[4]Zhao Jinping,Huang daji.Mirror extending and circular spline function for empirical mode decomposition method[J].Journal of Zhejiang University,2001,2(3):247-252.
[5]鄧擁軍,王偉,錢成春,等.EMD方法及 HILBERT變換中邊界問題的處理[J].科學通報,2001,46(3):257 - 263.
[6]劉慧婷,張旻,程家興,等.基于多項式擬合算法的 EMD端點問題的處理[J].計算機工程與應(yīng)用,2004,40(16):84-86,100.
[7]郭興明,湯麗平,陳麗珊,等.經(jīng)驗?zāi)J椒纸庠赒RS波群和T波檢測中的應(yīng)用[J].電子科技大學學報,2011,40(1):142-145.
[8]郭興明,林輝杰,肖守中.心音中醫(yī)學指標的提?。跩].計算機工程與應(yīng)用,2011,47(3):214 -216.
[9]元秀華,謝定,吳承德.心音信號測量中的噪音干擾分析與濾除方法[J].中國現(xiàn)代醫(yī)學雜志,1999,9(6):65,67.
[10]張亢,程軍圣,楊宇.基于自適應(yīng)波形匹配延拓的局部均值分解端點效應(yīng)處理方法[J].中國機械工程,2010,21(4):459-460.
[11]Liang H,Lukkarinen S,Hartimo I,Heart sound segmentation algorithm based on heart sound envelogram[J].Computers in Cardiology,1997,24:105-108.
[12]Liang Huiying,Lukkarinen Sakari,Hartimo Iiro.A heart sound segmentation algorithm using wavelet decomposition and reconstruction[C]// Proceedings of19th International Conference,IEEE-EMBS.Chicago,USA:IEEE,1997:1630 -1633.