研究基于希爾伯特黃變換(HHT)的心音包絡提取系統在LabVIEW上的完整實現。首先使用聲卡采集心音信號, 并在LabVIEW上實現了完整的基于HHT的心音采集、預處理和包絡提取功能的系統程序, 最后使用案例證明了該系統能夠簡便地實現心音信號采集、預處理和包絡提取。該系統較好地保留和顯示了心音包絡特征, 并且其程序和方法對振動、語音等的研究有重要的參考價值。
引用本文： 譚志向, 張懿, 曾德平, 王華. 基于希爾伯特-黃變換的心音包絡提取在LabVIEW上的實現. 生物醫學工程學雜志, 2015, 32(2): 263-268. doi: 10.7507/1001-5515.20150048 復制
引言
心音信號是一種非線性非平穩性生理信號，能反映心臟的多種生理病理信息，因此心音的檢測和分析是了解心臟、血管等健康狀態的一種重要手段。心音包絡不僅更明顯地突出了心音主要成分的特征，而且是心音特征識別的基礎，因此，準確地提取心音包絡具有重要意義[1]。
心音包絡提取的方法很多。Ruffo等[2]使用Teager能量算子和自相關算法提取胎心音包絡，該方法比較簡單，但是它對信號的信噪比要求較高。全海燕等[3]提出了實現較為復雜的基于小波多分辨率分析的心音分段算法。周靜[4]提出了歸一化平均香農能量分布法，該方法雖然很好地得出了心音的能量分布圖，但是在對心音分段計算時有重疊，影響了信號的時間關系。姚曉帥等[5]使用數學形態學算法做心音分段，該方法實現較簡單，但是它很難找到合適的結構元素來適合復雜的心音信號。趙治棟等[6]使用小波閾值法和Hilbert變換提取心音包絡，該方法中小波閾值選取的過程較復雜。
本文在研究和分析了目前信號包絡提取的發展現狀的基礎上，利用希爾伯特-黃變換(Hilbert-Huang transform，HHT)[7]的研究成果和虛擬儀器的優勢，開發了一個基于LabVIEW和聲卡的心音包絡提取系統，該系統為探索心音信號新的處理方法、提取心音特征等提供了一個廉價、快速和開放的平臺。
1 希爾伯特-黃變換
HHT使用信號的尺度特征對信號進行分解，是一種新的處理非線性非平穩信號的時頻分析方法。它不僅具有良好的局部適應特性，在引入瞬時頻率后，實現了信號的時頻分析，提高了信號處理的靈活性和有效性，而且HHT具有小波變換多分辨的優勢，同時克服了小波變換中小波基函數選擇的缺陷，所以HHT對非線性非平穩信號的處理比小波變換更有優勢[8]。
HHT理論中，經驗模態分解(empirical mode decomposition, EMD)是核心，圖 1是EMD的實現流程圖[9]。實現HHT時，先用EMD分解信號，獲得若干代表瞬時頻率分量的本征模態函數(intrinsic mode function，IMF)和一個殘余分量，再用Hilbert變換求各個IMF的解析函數，最后對各個解析函數進行處理得到Hilbert譜和邊際譜，即可得到信號的時頻域性質。
圖1
EMD實現流程圖
Figure1.
EMD realization flow chart
在圖 1中，IMF有兩個條件：①在整個時間序列中，所有極大值和極小值與過零點的數量相等或者最多只能差一個；②在任意時刻上，由信號的上包絡線和下包絡線所求的均值為零。第一個條件類似傳統平穩高斯過程的分布，而第二個條件要求信號的局部均值為零。為了解決非線性非平穩信號難以定義局部時間的缺陷，HHT用信號的上下包絡來近似局部對稱，并且限定IMF的條件，使得IMF只包含一個基本函數模式的振蕩波形，所以，任意信號在完成EMD分解后，得到了多個只包含一個基本函數模式的IMF，即實現了信號的局部表征。這個過程，我們不僅可以得到信號的時頻域性質，而且可以依據需要，去掉不需要的分量，對剩下的分量進行重構，即可實現信號除噪[10]。
使用HHT對心音處理時，先使用EMD將心音信號分解為多個IMF分量，然后按照心音的頻譜性質選擇包含心音主要部分的分量，去掉包含噪聲的分量，最后對選擇的分量進行Hilbert變換，即可實現對心音信號的預處理和包絡提取[11]。
2 基于LabVIEW的心音包絡提取系統
2.1 系統組成
如圖 2所示，系統組成的心音信號采集端是將一個3M聽診器(Littmann 2450)的耳件去掉，然后把一個聲傳感器(Panasonic WM-64C微型麥克風，響應頻率的范圍20~16 000 Hz)塞進聽診器的膠管內并密封，聲傳感器通過屏蔽導線連接一個3.5 mm四極音頻插頭，最后將該音頻插頭插入電腦的Line In接口[12]，在電腦上打開預先編寫好的LabVIEW軟件并調試可用，即可實現心音的采集和其它功能。
圖2
系統組成
Figure2.
Components of the system
系統基于HHT的心音包絡提取在LabVIEW上實現的功能包括：
(1)基于聲卡的心音采集、WAV格式文件存儲和回放。
(2)基于LabVIEW實現HHT中的EMD過程，實現對心音信號的預處理。
(3)選擇包含心音主成分的IMF，進行Hilbert變換，提取信號包絡。
心音頻段主要在500 Hz以下[13]，本論文實驗使用的電腦聲卡響應頻率范圍25~16 000 Hz，而通過聲卡采集的心音信號含有大量噪聲，在提取信號包絡之前先用EMD將采集信號分解，選擇出合適的IMF，去掉噪聲分量，實現采集信號的預處理。接著，在使用Hilbert變換提取信號包絡時，由于EMD中用樣條插值構造上下包絡有端點效應問題，在選定IMF后，先對該IMF進行Hilbert變換，然后尋找該Hilbert變換后的信號的極大值點，用三次樣條差值曲線擬合并順次連接這些極大值點，得到一條新的上包絡線，最后將這個上包絡線確定為最終所需要的信號包絡[14]。
參考李國成[15]和李博[16]的研究，將系統分為幾個主要模塊，第一個是基于聲卡的心音采集、存儲和回放模塊；第二個是EMD模塊，它通過調用IMF模塊和判斷剩余分量模塊實現，而IMF函數模塊又依次調用局部均值模塊和包絡線構造模塊。最后使用主函數將各個模塊組合在一起。
2.2 系統實現
(1)基于聲卡的心音采集、存儲和回放
如圖 3所示，程序(a)通過使用配置聲音輸入控件和讀取聲音輸入控件實現聲卡對聲音的采集，然后使用寫入聲音文件控件、打開聲音文件控件和文件對話框實現聲音文件的存儲，并且聲音的采集和存儲是同步的。程序(b)使用簡易讀取文件控件和文件選擇路徑實現存儲音頻文件的讀取回放。最后使用一個條件結構將基于聲卡的心音采集、文件存儲和回放結合在一起，默認使用(a)，當點擊布爾按鈕后，將關閉(a)，啟用(b)。
圖3
基于聲卡的心音采集、存儲和回放程度框圖
(a)基于聲卡的心音采集、WAV格式文件存儲程序框圖；(b)WAV格式文件讀取程序框圖
Figure3. Block diagram of heart sound collection, storage and playback based on the sound card(a) the block diagram of heart sounds collection and WAV format file storage based on sound card; (b) the block diagram of WAV format file reading
(2)三次樣條插值構造包絡線
如圖 4所示，原始信號的大小確定了程序中循環運算的次數，輸入極值位置和極值幅值后利用三次樣條插值，即可實現整個信號包絡線的構造。
圖4
使用三次樣條插值構造包絡線程序框圖
Figure4.
Block diagram of constructing the envelope using the cubic spline interpolation
(3)局部均值
如圖 5所示，對原信號分別使用Peak DetectorⅥ尋找信號的波峰和波谷，然后調用(1)的三次樣條子Ⅵ確定局部信號極值的包絡，最后對信號的上下包絡求均值，輸出局部均值。
圖5
局部均值程序框圖
Figure5.
Block diagram of local mean
(4)IMF
在圖 1的EMD實現流程圖中可知，將原始數據減去局部均值得到殘差H(t)，當H(t)滿足IMF的兩個條件時即可被認為是一個IMF。現實中我們常使用仿柯西收斂準則限定方差的方法來判定，公式為
| ${S_d}=\sum\limits_{t=0}^n {M_{k-2}^2\left(t \right)/\sum\limits_{t=0}^n {h_{k-1}^2\left(t \right), } } $ |
式中Mk-1(t)和hk-1(t)表示第k-1次產生的均值包絡和殘差，并且方差一般在0.2~0.3之間取值。如圖 6所示，方差選擇為0.25，當方差小于或等于0.25時，循環結束，輸出一個IMF。
圖6
IMF程序框圖
Figure6.
Block diagram of IMF
(5)EMD
實現EMD時，對剩余分量的判斷十分關鍵。每得到一個IMF，就會對應產生一個剩余分量，將此剩余分量當做輸入信號重復EMD的過程，直到最終的剩余分量滿足人為設定的結束條件即可。一般設定的結束條件主要包括：①剩余分量小于預定的誤差；②剩余分量為一單調函數，此時不可再從其中提取出IMF。所以，如圖 7所示，原始信號X(t)經過EMD分解后就變成了n個IMF和一個剩余分量Rn(t)的和。在圖 8的判斷剩余分量方法中，為了提高準確率，在利用Peak Detector控件對剩余分量采取極值點的限制時，也限制了剩余分量的均方誤差。
圖7
EMD程序框圖
Figure7.
Block diagram of EMD
圖8
判斷剩余分量程序框圖
Figure8.
Block diagram of determining the remaining components
(6)主函數
最后將各個模塊組合在主函數里，如圖 9所示，LabVIEW的程序是并行運行，所以程序中將條件選擇結構和選擇信號Ⅵ的控制按鈕連接在一起，保證了輸入的待處理信號與條件結構中選擇的模式的一致性，同時程序中使用順序結構將信號處理部分放在最后執行，使得各個模塊都能按照順序執行，避免了沖突。
圖9
主函數
Figure9.
Main function
3 使用案例
以采集一個25周歲的成人心音為例，將信號采集端置于被測試者胸壁的合適位置，如圖 10所示，啟動軟件配置好參數后開始采集。
圖10
心音采集
Figure10.
Heart sound acquisition
在采集信號時，單擊軟件前面板的“IMF與信號包絡”選板，在選板上通過選擇“選擇IMF分量”按鈕，可以查看各個IMF的波形圖及其包絡圖，這樣的好處主要是不僅可以清晰觀察各個IMF的波形圖和對應的包絡圖，還可以依據心音的性質選定效果最好的IMF和對應的包絡圖，甚至可以依據IMF和對應包絡圖的質量來判斷信號采集端的情況，可通過重新選定采集位置或改進信號采集端等方法獲得質量更好的心音信號。
如圖 11所示，在采集心音時，通過選擇獲得最好的IMF分量波形圖及其對應的信號包絡圖。
圖11
心音IMF分量圖及其信號包絡圖
Figure11.
IMF of heart sound and its signal envelope
4 結論
論文實現了基于LabVIEW的HHT方法，給出了詳細的實現程序，并經過實例驗證了該系統的功能。該系統的優點首先是使用HHT實現了心音預處理和包絡提取，很好地顯示了心音信號的能量分布，其次是將信號采集、存儲、預處理和包絡提取等通過電腦聲卡和LabVIEW組合成了一個快速的原型化系統，代替了昂貴的傳統儀器，加快了研究進度。而且，相比于歸一化平均香農能量法[17]，本論文的方法沒有破壞信號的時間特性，更好地保留了信號特征；相比于數學形態學方法，本論文的方法不需要大量的實驗確定結構元素的形狀和長度等，節省了研究精力和時間。
同時，從圖 11的結果中發現，信號包絡仍有缺陷，可能是心音信號中仍有噪聲干擾，EMD分解時并沒有完全濾掉所有噪聲，所以在信號預處理時可考慮加入其它預處理方法。
本系統快速地實現了心音的采集、預處理和包絡提取等，為心音的其它處理提供了一個廉價的原型化平臺，而且，該系統開放地允許使用者依據自己的需求嘗試更多方法來處理心音信號。此外，本系統實現的只是HHT一方面的應用，系統的程序和方法可以被拓展到更多的應用，例如應用于振動、語音和其它生物醫學信號的預處理、性質分析和包絡提取等，還可以對該系統進行更深入的探索，例如算法優化等，所以本論文的程序和方法對其它研究有重要的參考價值。
引言
心音信號是一種非線性非平穩性生理信號，能反映心臟的多種生理病理信息，因此心音的檢測和分析是了解心臟、血管等健康狀態的一種重要手段。心音包絡不僅更明顯地突出了心音主要成分的特征，而且是心音特征識別的基礎，因此，準確地提取心音包絡具有重要意義[1]。
心音包絡提取的方法很多。Ruffo等[2]使用Teager能量算子和自相關算法提取胎心音包絡，該方法比較簡單，但是它對信號的信噪比要求較高。全海燕等[3]提出了實現較為復雜的基于小波多分辨率分析的心音分段算法。周靜[4]提出了歸一化平均香農能量分布法，該方法雖然很好地得出了心音的能量分布圖，但是在對心音分段計算時有重疊，影響了信號的時間關系。姚曉帥等[5]使用數學形態學算法做心音分段，該方法實現較簡單，但是它很難找到合適的結構元素來適合復雜的心音信號。趙治棟等[6]使用小波閾值法和Hilbert變換提取心音包絡，該方法中小波閾值選取的過程較復雜。
本文在研究和分析了目前信號包絡提取的發展現狀的基礎上，利用希爾伯特-黃變換(Hilbert-Huang transform，HHT)[7]的研究成果和虛擬儀器的優勢，開發了一個基于LabVIEW和聲卡的心音包絡提取系統，該系統為探索心音信號新的處理方法、提取心音特征等提供了一個廉價、快速和開放的平臺。
1 希爾伯特-黃變換
HHT使用信號的尺度特征對信號進行分解，是一種新的處理非線性非平穩信號的時頻分析方法。它不僅具有良好的局部適應特性，在引入瞬時頻率后，實現了信號的時頻分析，提高了信號處理的靈活性和有效性，而且HHT具有小波變換多分辨的優勢，同時克服了小波變換中小波基函數選擇的缺陷，所以HHT對非線性非平穩信號的處理比小波變換更有優勢[8]。
HHT理論中，經驗模態分解(empirical mode decomposition, EMD)是核心，圖 1是EMD的實現流程圖[9]。實現HHT時，先用EMD分解信號，獲得若干代表瞬時頻率分量的本征模態函數(intrinsic mode function，IMF)和一個殘余分量，再用Hilbert變換求各個IMF的解析函數，最后對各個解析函數進行處理得到Hilbert譜和邊際譜，即可得到信號的時頻域性質。
圖1
EMD實現流程圖
Figure1.
EMD realization flow chart
在圖 1中，IMF有兩個條件：①在整個時間序列中，所有極大值和極小值與過零點的數量相等或者最多只能差一個；②在任意時刻上，由信號的上包絡線和下包絡線所求的均值為零。第一個條件類似傳統平穩高斯過程的分布，而第二個條件要求信號的局部均值為零。為了解決非線性非平穩信號難以定義局部時間的缺陷，HHT用信號的上下包絡來近似局部對稱，并且限定IMF的條件，使得IMF只包含一個基本函數模式的振蕩波形，所以，任意信號在完成EMD分解后，得到了多個只包含一個基本函數模式的IMF，即實現了信號的局部表征。這個過程，我們不僅可以得到信號的時頻域性質，而且可以依據需要，去掉不需要的分量，對剩下的分量進行重構，即可實現信號除噪[10]。
使用HHT對心音處理時，先使用EMD將心音信號分解為多個IMF分量，然后按照心音的頻譜性質選擇包含心音主要部分的分量，去掉包含噪聲的分量，最后對選擇的分量進行Hilbert變換，即可實現對心音信號的預處理和包絡提取[11]。
2 基于LabVIEW的心音包絡提取系統
2.1 系統組成
如圖 2所示，系統組成的心音信號采集端是將一個3M聽診器(Littmann 2450)的耳件去掉，然后把一個聲傳感器(Panasonic WM-64C微型麥克風，響應頻率的范圍20~16 000 Hz)塞進聽診器的膠管內并密封，聲傳感器通過屏蔽導線連接一個3.5 mm四極音頻插頭，最后將該音頻插頭插入電腦的Line In接口[12]，在電腦上打開預先編寫好的LabVIEW軟件并調試可用，即可實現心音的采集和其它功能。
圖2
系統組成
Figure2.
Components of the system
系統基于HHT的心音包絡提取在LabVIEW上實現的功能包括：
(1)基于聲卡的心音采集、WAV格式文件存儲和回放。
(2)基于LabVIEW實現HHT中的EMD過程，實現對心音信號的預處理。
(3)選擇包含心音主成分的IMF，進行Hilbert變換，提取信號包絡。
心音頻段主要在500 Hz以下[13]，本論文實驗使用的電腦聲卡響應頻率范圍25~16 000 Hz，而通過聲卡采集的心音信號含有大量噪聲，在提取信號包絡之前先用EMD將采集信號分解，選擇出合適的IMF，去掉噪聲分量，實現采集信號的預處理。接著，在使用Hilbert變換提取信號包絡時，由于EMD中用樣條插值構造上下包絡有端點效應問題，在選定IMF后，先對該IMF進行Hilbert變換，然后尋找該Hilbert變換后的信號的極大值點，用三次樣條差值曲線擬合并順次連接這些極大值點，得到一條新的上包絡線，最后將這個上包絡線確定為最終所需要的信號包絡[14]。
參考李國成[15]和李博[16]的研究，將系統分為幾個主要模塊，第一個是基于聲卡的心音采集、存儲和回放模塊；第二個是EMD模塊，它通過調用IMF模塊和判斷剩余分量模塊實現，而IMF函數模塊又依次調用局部均值模塊和包絡線構造模塊。最后使用主函數將各個模塊組合在一起。
2.2 系統實現
(1)基于聲卡的心音采集、存儲和回放
如圖 3所示，程序(a)通過使用配置聲音輸入控件和讀取聲音輸入控件實現聲卡對聲音的采集，然后使用寫入聲音文件控件、打開聲音文件控件和文件對話框實現聲音文件的存儲，并且聲音的采集和存儲是同步的。程序(b)使用簡易讀取文件控件和文件選擇路徑實現存儲音頻文件的讀取回放。最后使用一個條件結構將基于聲卡的心音采集、文件存儲和回放結合在一起，默認使用(a)，當點擊布爾按鈕后，將關閉(a)，啟用(b)。
圖3
基于聲卡的心音采集、存儲和回放程度框圖
(a)基于聲卡的心音采集、WAV格式文件存儲程序框圖；(b)WAV格式文件讀取程序框圖
Figure3. Block diagram of heart sound collection, storage and playback based on the sound card(a) the block diagram of heart sounds collection and WAV format file storage based on sound card; (b) the block diagram of WAV format file reading
(2)三次樣條插值構造包絡線
如圖 4所示，原始信號的大小確定了程序中循環運算的次數，輸入極值位置和極值幅值后利用三次樣條插值，即可實現整個信號包絡線的構造。
圖4
使用三次樣條插值構造包絡線程序框圖
Figure4.
Block diagram of constructing the envelope using the cubic spline interpolation
(3)局部均值
如圖 5所示，對原信號分別使用Peak DetectorⅥ尋找信號的波峰和波谷，然后調用(1)的三次樣條子Ⅵ確定局部信號極值的包絡，最后對信號的上下包絡求均值，輸出局部均值。
圖5
局部均值程序框圖
Figure5.
Block diagram of local mean
(4)IMF
在圖 1的EMD實現流程圖中可知，將原始數據減去局部均值得到殘差H(t)，當H(t)滿足IMF的兩個條件時即可被認為是一個IMF。現實中我們常使用仿柯西收斂準則限定方差的方法來判定，公式為
| ${S_d}=\sum\limits_{t=0}^n {M_{k-2}^2\left(t \right)/\sum\limits_{t=0}^n {h_{k-1}^2\left(t \right), } } $ |
式中Mk-1(t)和hk-1(t)表示第k-1次產生的均值包絡和殘差，并且方差一般在0.2~0.3之間取值。如圖 6所示，方差選擇為0.25，當方差小于或等于0.25時，循環結束，輸出一個IMF。
圖6
IMF程序框圖
Figure6.
Block diagram of IMF
(5)EMD
實現EMD時，對剩余分量的判斷十分關鍵。每得到一個IMF，就會對應產生一個剩余分量，將此剩余分量當做輸入信號重復EMD的過程，直到最終的剩余分量滿足人為設定的結束條件即可。一般設定的結束條件主要包括：①剩余分量小于預定的誤差；②剩余分量為一單調函數，此時不可再從其中提取出IMF。所以，如圖 7所示，原始信號X(t)經過EMD分解后就變成了n個IMF和一個剩余分量Rn(t)的和。在圖 8的判斷剩余分量方法中，為了提高準確率，在利用Peak Detector控件對剩余分量采取極值點的限制時，也限制了剩余分量的均方誤差。
圖7
EMD程序框圖
Figure7.
Block diagram of EMD
圖8
判斷剩余分量程序框圖
Figure8.
Block diagram of determining the remaining components
(6)主函數
最后將各個模塊組合在主函數里，如圖 9所示，LabVIEW的程序是并行運行，所以程序中將條件選擇結構和選擇信號Ⅵ的控制按鈕連接在一起，保證了輸入的待處理信號與條件結構中選擇的模式的一致性，同時程序中使用順序結構將信號處理部分放在最后執行，使得各個模塊都能按照順序執行，避免了沖突。
圖9
主函數
Figure9.
Main function
3 使用案例
以采集一個25周歲的成人心音為例，將信號采集端置于被測試者胸壁的合適位置，如圖 10所示，啟動軟件配置好參數后開始采集。
圖10
心音采集
Figure10.
Heart sound acquisition
在采集信號時，單擊軟件前面板的“IMF與信號包絡”選板，在選板上通過選擇“選擇IMF分量”按鈕，可以查看各個IMF的波形圖及其包絡圖，這樣的好處主要是不僅可以清晰觀察各個IMF的波形圖和對應的包絡圖，還可以依據心音的性質選定效果最好的IMF和對應的包絡圖，甚至可以依據IMF和對應包絡圖的質量來判斷信號采集端的情況，可通過重新選定采集位置或改進信號采集端等方法獲得質量更好的心音信號。
如圖 11所示，在采集心音時，通過選擇獲得最好的IMF分量波形圖及其對應的信號包絡圖。
圖11
心音IMF分量圖及其信號包絡圖
Figure11.
IMF of heart sound and its signal envelope
4 結論
論文實現了基于LabVIEW的HHT方法，給出了詳細的實現程序，并經過實例驗證了該系統的功能。該系統的優點首先是使用HHT實現了心音預處理和包絡提取，很好地顯示了心音信號的能量分布，其次是將信號采集、存儲、預處理和包絡提取等通過電腦聲卡和LabVIEW組合成了一個快速的原型化系統，代替了昂貴的傳統儀器，加快了研究進度。而且，相比于歸一化平均香農能量法[17]，本論文的方法沒有破壞信號的時間特性，更好地保留了信號特征；相比于數學形態學方法，本論文的方法不需要大量的實驗確定結構元素的形狀和長度等，節省了研究精力和時間。
同時，從圖 11的結果中發現，信號包絡仍有缺陷，可能是心音信號中仍有噪聲干擾，EMD分解時并沒有完全濾掉所有噪聲，所以在信號預處理時可考慮加入其它預處理方法。
本系統快速地實現了心音的采集、預處理和包絡提取等，為心音的其它處理提供了一個廉價的原型化平臺，而且，該系統開放地允許使用者依據自己的需求嘗試更多方法來處理心音信號。此外，本系統實現的只是HHT一方面的應用，系統的程序和方法可以被拓展到更多的應用，例如應用于振動、語音和其它生物醫學信號的預處理、性質分析和包絡提取等，還可以對該系統進行更深入的探索，例如算法優化等，所以本論文的程序和方法對其它研究有重要的參考價值。
