利用心阻抗微分信號的特征點可計算出多個血流動力學參數, 進而判別心功能狀態, 因此特征點的準確提取顯得尤為重要。本文應用實驗室自行設計開發的KF_ICG型無創心功能檢測儀采集了健康人和重慶市大坪醫院22例患者數據, 應用小波閾值法對采集的數據進行降噪處理, 對降噪后的信號采用bior3.7小波進行6層分解后定位特征點。結果表明, 該法無論對健康人還是存在諸多噪聲干擾的臨床患者數據都能有效實現特征點的精確定位, 有助于實現阻抗法無創檢測血流動力學參數的臨床應用。
引用本文： 趙云冬, 季忠, 彭承琳, 霍威. 基于小波變換的心阻抗微分信號去噪及特征點檢測研究. 生物醫學工程學雜志, 2015, 32(2): 284-289. doi: 10.7507/1001-5515.20150052 復制
引言
臨床上，阻抗法是一種無創檢測心輸出量(cardiac output, CO)、每搏輸出量(stroke volume, SV)等多個血流動力學參數的方法[1-2]。目前常用于計算SV的有Kubicek、Bernstein、Sramek等計算公式，其中Kubicek公式如式(1)所示，由該式看出：為得到SV，圖 1心阻抗微分波形中[3]，室縮波C波峰值(dz/dt)max、心室射血始點B點到終點X點間的間隔，亦即左室射血時間(left ventricular ejection time, LVET)必不可少。臨床上采集的信號常常含有諸多噪聲和干擾，幅值較小的特征點在波形中不明顯，造成定位困難，因此有效的去噪方法顯得尤為重要。傳統的頻域去噪法應用頻譜分析技術處理，對于頻率波動大的非平穩信號由于信號本身與噪聲間存在混疊，選取的截止頻率過大達不到去噪目的，過小則將有用信號濾除，導致不能很好地保存有用信號的尖峰和突變部分，這給信號特征點的檢測帶來不利，也就是說傳統的去噪方法存在著保護信號局部性和抑制噪聲之間的矛盾[4]。小波分析具有良好的時頻局部化性質，具有多分辨率帶通濾波特性，在低頻部分有較高的頻率分辨率和較低的時間分辨率，高頻部分有較高的時間分辨率和較低的頻率分辨率，很適合處理正常信號中夾帶的瞬態反常現象[5]。小波分析屬于時頻分析的一種，能同時在時域和頻域對信號進行分析，能有效區分信號的噪聲和突變部分，因而在信噪分離、特征點提取中顯示出獨有的優越性[6]。
| ${\rm{Kubicek}}:{\rm{SV}}={\rho _{\rm{b}}}{\left({\frac{{\rm{L}}}{{{{\rm{Z}}_0}}}} \right)^2}*{\left({\frac{{{\rm{d}}z}}{{{\rm{d}}t}}} \right)_{\max }}*{\rm{LVET}}$ |
圖1
典型的心阻抗血流圖及其相關波形
Figure1.
Typical impedance cardiogram and its related waveforms
本文對采集的心阻抗微分信號應用巴特沃斯帶通濾波器和小波閾值兩種方法進行除噪，對除噪后的信號應用bior3.7小波進行6層分解后成功定位信號特征點，對健康人和患者信號的C波、B點、X點定位的準確率進行了統計，均在97%以上。
1 小波變換法
傳統傅里葉變換將信號展為不同頻率的正弦信號線性疊加，核函數為正弦函數，在時間域上是無限非局部化的，當信號和噪聲在同一頻段疊加時顯得無能為力[7]。而小波變換的核函數是多種多樣的，選擇適當的母小波后，通過對母小波平移、伸縮形成一系列小波，將欲分析的信號投影到小波構成的信號空間中，解決了時域局部化與頻域局部化的矛盾，從而實現在不同的頻率范圍、時間位置對信號的檢測分析[8]。
信號f(x)的小波變換卷積形式為
| ${\rm{W}}{{\rm{T}}_a}{\rm{x}}\left({\rm{t}} \right)={\rm{x}}*{\varphi _{a, \tau }}\left(t \right)=\frac{1}{a}\int {{\rm{x}}\left(\tau \right)\varphi \left({\frac{{t-\tau }}{a}} \right){\rm{d}}t}, $ |
式中a為縮放因子，對應于頻率信息；τ為平移因子，對應于時空信息；φa, τ(t)為母小波。
利用小波變換法檢測信號奇異點的原理是：式(2)中，設θ(t)為一低通平滑函數，且滿足：
| $\int\limits_{-\infty }^\infty {\theta \left(t \right){\rm{d}}t=1, \mathop {\lim }\limits_{\left| {t \to \infty } \right|} } \theta \left(t \right)=0$ |
定義θ(t)的一階、二階導數為：、，則函數φ(1)(t)、φ(2)(t)滿足小波的可容許條件，即φ(1)(t)、φ(2)(t)可作為母小波函數。
函數f(t)在尺度為a、位置為t處的卷積型小波變換為
| $\begin{array}{c} {\rm{W}}_{\rm{a}}^{\left(1 \right)}f\left(t \right)=f*\varphi _a^{\left(1 \right)}\left(t \right)=\frac{1}{a}\int\limits_{-\infty }^{ + \infty } {f\left(\tau \right){\varphi ^{\left(1 \right)}}\left({\frac{{t-\tau }}{a}} \right){\rm{d}}\tau } \\=a\frac{{\rm{d}}}{{{\rm{d}}t}}\left({f*{\theta _a}} \right)\left(t \right) \end{array}$ |
| $\begin{array}{c} {\rm{W}}_{\rm{a}}^{\left(2 \right)}f\left(t \right)=f*\varphi _a^{\left(2 \right)}\left(t \right)=\frac{1}{a}\int\limits_{-\infty }^{ + \infty } {f\left(\tau \right){\varphi ^{\left(2 \right)}}\left({\frac{{t-\tau }}{a}} \right){\rm{d}}\tau } \\={a^2}\frac{{{{\rm{d}}^2}}}{{{\rm{d}}{t^2}}}\left({f*{\theta _a}} \right)\left(t \right) \end{array}$ |
由此看出f(t)關于φ(1)(t)、φ(2)(t)的規范小波變換變成了與光滑函數θa(t)的卷積關于t的一階和二階導數乘以a與a2。這樣Wa(1)f(t)的局部極值點對應Wa(2)f(t)的零交叉點和f(t)*θa(t)的拐點，如圖 2所示，即信號f(t)的特征點可通過小波變換后正負模極大值對間的過零點進行檢測識別[9]。
圖2
信號局部極值與零交叉檢測
Figure2.
Signal location of extrema and zero points
2 實驗處理
2.1 信號來源
本次實驗數據來源于實驗室自行設計開發的KF_ICG型心功能無創檢測分析儀，儀器主界面如圖 3所示。采集的實驗數據構成如下：實驗室健康人6例，年齡范圍22～26歲，平均年齡23.17歲，男性5例，女性1例；重慶市大坪醫院冠心病患者22例，年齡范圍41～80歲，平均年齡62.29歲，男性17例，女性5例。經倫理委員會審查通過，并簽署了患者知情同意書。
圖3
儀器檢測和分析界面
圖中信號自上而下分別為：阻抗變化信號、阻抗微分信號、心電信號、心音信號
Figure3. Instrument interface window for measurement and analysisSignals from top to bottom: impedance change signal, impedance differential signal, electrocardiographic signal, and phonocardiogram signal
2.2 去噪處理
2.2.1 小波閾值法去噪
去噪前，有必要先了解心阻抗信號本身的特點。心阻抗信號對應一個心動周期內主動脈瓣的開閉引起血流量的變化，其能量分布與主動脈瓣開閉的頻率(30～45 Hz)相近[8]，有用信號主要集中在低頻部分，噪聲主要包含：①基線漂移。一般由人體呼吸和電極移動引起，頻率低于10 Hz。②肌電干擾。由人體肌肉顫動所致，頻率范圍廣, 在5 Hz～2 kHz之間。③工頻干擾。頻率固定為50 Hz及其諧波成分。④心電噪聲。位于心臟附近的電極不可避免地采集到心電等其他干擾信號。
應用小波閾值法消除信號噪聲的原理是：原始信號x(n)經小波分解后產生低頻近似分量和高頻細節分量，二者構成一個尺度的分解，分解后的近似分量還可繼續向下分解構成小波分解樹。小波分解后，每一尺度對應各自的小波系數，Donoho提出有用信號產生的小波系數較大且數目較少，而噪聲產生的小波系數幅值較小，這樣選擇合適的閾值對各層小波系數進行收縮處理，再用處理后的小波系數重構信號即可達到去噪目的[10]。具體步驟如下：
(1) 信號的小波分解。選擇合適的母小波并確定分解層次后進行分解運算。
(2) 小波分解高頻系數的閾值量化。
(3) 信號重構。應用處理后的小波系數重構信號，達到消噪目的。
母小波的選取：以上步驟中，母小波的選取和閾值大小的確定是去噪效果的決定因素[11]。母小波具有不唯一性，對同一信號采用不同的母小波分解結果不盡相同，選取的標準有[12]：①支集長度：支集越短，局部化能力越好；②消失矩階數：主要針對數據壓縮；③對稱性：在信號分解與重構中避免相移；④正則性：主要影響信號重構穩定性。本文選取的bior小波是雙正交樣條小波，具有對稱性、緊支撐、線性相位等特性，能有效減少信號誤差，廣泛應用于非線性、非穩定的生理信號處理[12]。結合待處理的心阻抗信號特征，本文最終選取bior3.7小波進行分解重構。
閾值的選取：通常分為硬閾值處理和軟閾值處理[10]，表達式見式(5)，式中w為小波系數、λ為閾值、wλ為閾值處理后的小波系數。軟閾值在閾值處連續，去噪后能產生光滑的結果，連續性好。而硬閾值在閾值處不連續，得到的估計小波系數連續性差, 在某些點上會產生間斷，引起重構信號的吉布斯振蕩[13]，因此本文選擇軟閾值進行處理。
| $\left\{ \begin{array}{l}硬閾值:w\lambda=\left\{ \begin{array}{l} w, \left| w \right| \ge \lambda \\ 0, \left| w \right| < \lambda \end{array} \right.\\軟閾值:w\lambda=\left\{ \begin{array}{l} \left[{{\rm{sign}}\left(w \right)} \right]\left({\left| w \right|-\lambda } \right), \left| w \right| \ge \lambda \\ 0, \left| w \right| < \lambda \end{array} \right. \end{array} \right.$ |
上式中，閾值λ起決定性作用，λ太小噪聲去除不干凈，λ太大則會將有用信號成分去除造成信號失真。常用的閾值計算方法有：固定閾值、基于極大極小原理閾值、自適應閾值。若采用全局固定閾值處理則忽略了各層系數間的差異，導致去噪效果不理想，因此本文采用基于無偏似然估計原理的自適應閾值，即史坦無偏估計(Stein’s unbiased risk estimation, SURE)。對于一個給定閾值，求出其對應的似然估計，在進行非似然最小化，從而得出所選閾值[14]。
設T為一向量，其元素為按從小到大順序排列的小波系數平方，即
| ${\rm{\boldsymbol{T}}}=\left\{ {{w_1}, {w_2}, \cdots, {w_n}} \right\}, $ |
其中w1≤w2…≤wn，再設一風險向量R，其元素如式(7)，求出R元素中最小值rj所對應的下標j，找出wj，則可求的閾值。
| ${r_i}=\left[{{\rm{n}}-2i + \left({n-i} \right){w_i} + \sum\limits_{k=1}^i {{w_k}} } \right]/{\rm{n}}\left({i=1, 2, \cdots, {\rm{n}}} \right)$ |
2.2.2 除噪及效果評價
根據上文步驟，對圖 3所示設備采集到的心阻抗微分信號進行去噪處理。圖 4為利用bior3.7小波對心阻抗微分信號進行6層分解后各尺度下的近似分量和細節分量，由圖看出高頻噪聲分量主要集中在d1、d2、d3上，有用信號集中于d4、d5上。除噪時，根據SURE選取合適的閾值對各層細節分量進行收縮處理，處理時盡可能減小富含噪聲的d1、d2、d3尺度上的小波系數，保留有用信號集中的d4、d5尺度上的系數，這樣利用收縮處理后的各層小波系數加之近似信號a6重構出信號即可達到消噪目的。
圖4
bior3.7小波6層分解后各尺度小波系數
Figure4.
Each scale wavelet coefficients after bior3.7 six layer wavelet decomposition
消噪結果如圖 5所示：(a) 為健康人心阻抗微分信號，與圖 1相比，波形標準、特征點明顯，存在少量噪聲。(b) 為一患者心阻抗微分信號，該患者為男性，57歲，患有冠心病、不穩定性心絞痛。由圖中看出信號受噪聲干擾嚴重，波形在O波處發生變異。(c) 為圖 5(a)除噪后信號，與(a)相比，波形更為平滑干凈，且波峰未被削減。(d) 為圖 5(b)除噪后信號，與(b)相比，除噪后信號變得平滑干凈，特征點凸顯，去噪效果十分明顯。
為了能定量地評價去噪效果，在此采用信噪比(signal to noise ratio, SNR)和均方誤差(mean square error, MSE)進行評價。SNR指有用信號與噪聲的比值，越大表示信號所含噪聲越少；MSE體現了原始信號和去噪后信號間的差異，越小表示去噪效果越好[12]。
| $\left\{ \begin{array}{l} {\rm{SNR}}=10{\log _{10}}\sum\limits_{i=1}^N {{{\left({\frac{{{{\rm{x}}_i}}}{{{{\rm{x}}_i}-{{\rm{y}}_i}}}} \right)}^2}} \\ {\rm{MSE}}=\sqrt {\frac{1}{N}\sum\limits_{i=1}^N {{{\left({{{\rm{y}}_i}-{{\rm{x}}_i}} \right)}^2}} } \end{array} \right.$ |
現實中，完全不含噪聲的心阻抗微分信號是不存在的，因此本文采用圖 5(c)中的信號近似為真實信號。圖 5(c)中信號為實驗室24歲健康男性采用小波閾值法除噪后的信號。式(8)中xi為圖 5(c)中信號，yi為醫院采集的患者去噪后信號，N表示信號長度。
表1
兩種方法去噪后6組信號的信噪比和均方誤差
Table1.
Six groups of signal SNR and MSE after using two de-noising methods
2.3 特征點檢測
準確定位B、C、X各特征點是正確計算SV、CO等血流動力學參數的基礎。信號經小波分解后，得到不同時間尺度上的幅度特性，在某些尺度上，有用信號的特征點十分明顯，同時噪聲被充分抑制，在這些尺度上檢測信號的模極大值對間過零點來定位特征點[7]。心阻抗微分信號的能量分布與主動脈瓣開閉的頻率(30～45 Hz)相近[8]，本文中采樣頻率為1 000 Hz, 根據奈奎斯特采樣定理，采樣率為1 000 Hz的信號最高頻率為500 Hz，6層分解后，尺度4、5集中包含了心阻抗信號的主要能量，高頻系數d5對應的頻帶為15.625～31.25 Hz，與阻抗信號的頻率相近，而肌電干擾等噪聲在S=25尺度下大大減小，結合大量實驗驗證，在S=25尺度下特征點檢測效果最好。
實驗室研制的KF_ICG型無創心功能檢測儀能夠實現信號特征點的自動檢測識別。信號采集完成后(或回放后)，點擊圖 3中的分析按鈕將會自動彈出圖 6所示的特征點檢測界面。圖 6中左下框內顯示了一個周期內完整的心阻抗微分波形，白色虛線即為各特征點檢測結果，由圖看出經小波變換分解后，在感興趣的尺度上成功定位出心阻抗微分信號的特征點。特征點定位完成后，界面右半部分計算出常用的血流動力學參數便于臨床醫生診斷。
圖6
儀器檢測分析界面
Figure6.
Detection and analysis interface of KF_ICG instrument
3 結果與討論
用采集到的實驗室健康人和22例患者數據進行特征點檢測驗證，對同一研究對象分別采用差分閾值法和小波變換法檢測心阻抗微分信號的特征點，最終統計結果如表 2所示。
表2
兩種方法對健康人和患者的C波、B點、X點的檢測結果
Table2.
C wave, B point and X point detection of healthy people and patients by using traditional differential threshold and wavelet threshold
表 2列出了傳統差分閾值法和小波變換法對健康人和患者阻抗微分信號中C波、B點、X點定位的準確率。由表 2得出以下結論：(1)小波法對健康人和患者的定位準確率均高于差分閾值法，但對健康信號二者相差不大。分析原因為健康人年輕(平均23.17歲)，心功能狀態較好，波形幅值大，特征點明顯不易錯檢，因此針對健康信號兩種方法的識別率均較高。(2)針對摻雜了諸多噪聲干擾的臨床信號，小波法的準確率明顯高于閾值法，特別是對幅值小、易被淹沒的B點、X點更為明顯。分析原因：①患者年齡較大(平均62.29歲)，心功能減弱，波形幅值小，B點、X點等特征點在波形圖中不夠明顯；②檢測時要求受試者在安靜的環境下放松靜臥、勻稱呼吸，健康人檢測地點為實驗室，能較好地滿足環境要求；臨床患者檢測時環境嘈雜，特別是某些患者呼吸急促，帶來幅值較大的干擾和噪聲使B點、X點被淹沒其中，而差分閾值法的抗噪能力較差，造成錯檢點增多；③某些患者可能存在影響阻抗檢測的病癥，導致檢測到的波形變異，致使特征點較難定位。但總體說來，小波變換法對健康人和患者的準確率均較高，均在97%以上，結合圖 3儀器檢測界面中列出的臨床常用心功能診斷參數SV、CO值，與重慶市大坪醫院心內科采用的超聲多普勒USCOM1A型儀器的測量值相對比，誤差均在10%以內，證明該法準確有效。
4 結束語
應用阻抗法檢測SV、CO等血流動力學參數是無創、經濟的。為了保證計算機對摻雜了諸多噪聲干擾的臨床數據也能精確地自動識別診斷信號特征點，前期高質量的降噪必不可少。傳統的頻域去噪對頻率波動大的非平穩信號顯得無能為力，而小波變換具有多分辨率分析功能，在生理信號去噪與特征點的檢測中顯示出獨有的優越性[15]。本文應用小波閾值法對采集到的信號進行消噪，并在消噪的基礎上對信號特征點進行自動識別檢測，形成一套完整的心功能信號處理系統。利用該系統采集了28例樣本進行驗證，結果表明該系統無論對健康人還是患者信號都有較好的抗噪識別能力。下一步工作是采集更多的臨床數據，針對某些波形發生變異的數據進行改進，提高設備對不同病癥患者個體差異的魯棒性，從而真正服務于臨床診斷應用。
引言
臨床上，阻抗法是一種無創檢測心輸出量(cardiac output, CO)、每搏輸出量(stroke volume, SV)等多個血流動力學參數的方法[1-2]。目前常用于計算SV的有Kubicek、Bernstein、Sramek等計算公式，其中Kubicek公式如式(1)所示，由該式看出：為得到SV，圖 1心阻抗微分波形中[3]，室縮波C波峰值(dz/dt)max、心室射血始點B點到終點X點間的間隔，亦即左室射血時間(left ventricular ejection time, LVET)必不可少。臨床上采集的信號常常含有諸多噪聲和干擾，幅值較小的特征點在波形中不明顯，造成定位困難，因此有效的去噪方法顯得尤為重要。傳統的頻域去噪法應用頻譜分析技術處理，對于頻率波動大的非平穩信號由于信號本身與噪聲間存在混疊，選取的截止頻率過大達不到去噪目的，過小則將有用信號濾除，導致不能很好地保存有用信號的尖峰和突變部分，這給信號特征點的檢測帶來不利，也就是說傳統的去噪方法存在著保護信號局部性和抑制噪聲之間的矛盾[4]。小波分析具有良好的時頻局部化性質，具有多分辨率帶通濾波特性，在低頻部分有較高的頻率分辨率和較低的時間分辨率，高頻部分有較高的時間分辨率和較低的頻率分辨率，很適合處理正常信號中夾帶的瞬態反常現象[5]。小波分析屬于時頻分析的一種，能同時在時域和頻域對信號進行分析，能有效區分信號的噪聲和突變部分，因而在信噪分離、特征點提取中顯示出獨有的優越性[6]。
| ${\rm{Kubicek}}:{\rm{SV}}={\rho _{\rm{b}}}{\left({\frac{{\rm{L}}}{{{{\rm{Z}}_0}}}} \right)^2}*{\left({\frac{{{\rm{d}}z}}{{{\rm{d}}t}}} \right)_{\max }}*{\rm{LVET}}$ |
圖1
典型的心阻抗血流圖及其相關波形
Figure1.
Typical impedance cardiogram and its related waveforms
本文對采集的心阻抗微分信號應用巴特沃斯帶通濾波器和小波閾值兩種方法進行除噪，對除噪后的信號應用bior3.7小波進行6層分解后成功定位信號特征點，對健康人和患者信號的C波、B點、X點定位的準確率進行了統計，均在97%以上。
1 小波變換法
傳統傅里葉變換將信號展為不同頻率的正弦信號線性疊加，核函數為正弦函數，在時間域上是無限非局部化的，當信號和噪聲在同一頻段疊加時顯得無能為力[7]。而小波變換的核函數是多種多樣的，選擇適當的母小波后，通過對母小波平移、伸縮形成一系列小波，將欲分析的信號投影到小波構成的信號空間中，解決了時域局部化與頻域局部化的矛盾，從而實現在不同的頻率范圍、時間位置對信號的檢測分析[8]。
信號f(x)的小波變換卷積形式為
| ${\rm{W}}{{\rm{T}}_a}{\rm{x}}\left({\rm{t}} \right)={\rm{x}}*{\varphi _{a, \tau }}\left(t \right)=\frac{1}{a}\int {{\rm{x}}\left(\tau \right)\varphi \left({\frac{{t-\tau }}{a}} \right){\rm{d}}t}, $ |
式中a為縮放因子，對應于頻率信息；τ為平移因子，對應于時空信息；φa, τ(t)為母小波。
利用小波變換法檢測信號奇異點的原理是：式(2)中，設θ(t)為一低通平滑函數，且滿足：
| $\int\limits_{-\infty }^\infty {\theta \left(t \right){\rm{d}}t=1, \mathop {\lim }\limits_{\left| {t \to \infty } \right|} } \theta \left(t \right)=0$ |
定義θ(t)的一階、二階導數為：、，則函數φ(1)(t)、φ(2)(t)滿足小波的可容許條件，即φ(1)(t)、φ(2)(t)可作為母小波函數。
函數f(t)在尺度為a、位置為t處的卷積型小波變換為
| $\begin{array}{c} {\rm{W}}_{\rm{a}}^{\left(1 \right)}f\left(t \right)=f*\varphi _a^{\left(1 \right)}\left(t \right)=\frac{1}{a}\int\limits_{-\infty }^{ + \infty } {f\left(\tau \right){\varphi ^{\left(1 \right)}}\left({\frac{{t-\tau }}{a}} \right){\rm{d}}\tau } \\=a\frac{{\rm{d}}}{{{\rm{d}}t}}\left({f*{\theta _a}} \right)\left(t \right) \end{array}$ |
| $\begin{array}{c} {\rm{W}}_{\rm{a}}^{\left(2 \right)}f\left(t \right)=f*\varphi _a^{\left(2 \right)}\left(t \right)=\frac{1}{a}\int\limits_{-\infty }^{ + \infty } {f\left(\tau \right){\varphi ^{\left(2 \right)}}\left({\frac{{t-\tau }}{a}} \right){\rm{d}}\tau } \\={a^2}\frac{{{{\rm{d}}^2}}}{{{\rm{d}}{t^2}}}\left({f*{\theta _a}} \right)\left(t \right) \end{array}$ |
由此看出f(t)關于φ(1)(t)、φ(2)(t)的規范小波變換變成了與光滑函數θa(t)的卷積關于t的一階和二階導數乘以a與a2。這樣Wa(1)f(t)的局部極值點對應Wa(2)f(t)的零交叉點和f(t)*θa(t)的拐點，如圖 2所示，即信號f(t)的特征點可通過小波變換后正負模極大值對間的過零點進行檢測識別[9]。
圖2
信號局部極值與零交叉檢測
Figure2.
Signal location of extrema and zero points
2 實驗處理
2.1 信號來源
本次實驗數據來源于實驗室自行設計開發的KF_ICG型心功能無創檢測分析儀，儀器主界面如圖 3所示。采集的實驗數據構成如下：實驗室健康人6例，年齡范圍22～26歲，平均年齡23.17歲，男性5例，女性1例；重慶市大坪醫院冠心病患者22例，年齡范圍41～80歲，平均年齡62.29歲，男性17例，女性5例。經倫理委員會審查通過，并簽署了患者知情同意書。
圖3
儀器檢測和分析界面
圖中信號自上而下分別為：阻抗變化信號、阻抗微分信號、心電信號、心音信號
Figure3. Instrument interface window for measurement and analysisSignals from top to bottom: impedance change signal, impedance differential signal, electrocardiographic signal, and phonocardiogram signal
2.2 去噪處理
2.2.1 小波閾值法去噪
去噪前，有必要先了解心阻抗信號本身的特點。心阻抗信號對應一個心動周期內主動脈瓣的開閉引起血流量的變化，其能量分布與主動脈瓣開閉的頻率(30～45 Hz)相近[8]，有用信號主要集中在低頻部分，噪聲主要包含：①基線漂移。一般由人體呼吸和電極移動引起，頻率低于10 Hz。②肌電干擾。由人體肌肉顫動所致，頻率范圍廣, 在5 Hz～2 kHz之間。③工頻干擾。頻率固定為50 Hz及其諧波成分。④心電噪聲。位于心臟附近的電極不可避免地采集到心電等其他干擾信號。
應用小波閾值法消除信號噪聲的原理是：原始信號x(n)經小波分解后產生低頻近似分量和高頻細節分量，二者構成一個尺度的分解，分解后的近似分量還可繼續向下分解構成小波分解樹。小波分解后，每一尺度對應各自的小波系數，Donoho提出有用信號產生的小波系數較大且數目較少，而噪聲產生的小波系數幅值較小，這樣選擇合適的閾值對各層小波系數進行收縮處理，再用處理后的小波系數重構信號即可達到去噪目的[10]。具體步驟如下：
(1) 信號的小波分解。選擇合適的母小波并確定分解層次后進行分解運算。
(2) 小波分解高頻系數的閾值量化。
(3) 信號重構。應用處理后的小波系數重構信號，達到消噪目的。
母小波的選取：以上步驟中，母小波的選取和閾值大小的確定是去噪效果的決定因素[11]。母小波具有不唯一性，對同一信號采用不同的母小波分解結果不盡相同，選取的標準有[12]：①支集長度：支集越短，局部化能力越好；②消失矩階數：主要針對數據壓縮；③對稱性：在信號分解與重構中避免相移；④正則性：主要影響信號重構穩定性。本文選取的bior小波是雙正交樣條小波，具有對稱性、緊支撐、線性相位等特性，能有效減少信號誤差，廣泛應用于非線性、非穩定的生理信號處理[12]。結合待處理的心阻抗信號特征，本文最終選取bior3.7小波進行分解重構。
閾值的選取：通常分為硬閾值處理和軟閾值處理[10]，表達式見式(5)，式中w為小波系數、λ為閾值、wλ為閾值處理后的小波系數。軟閾值在閾值處連續，去噪后能產生光滑的結果，連續性好。而硬閾值在閾值處不連續，得到的估計小波系數連續性差, 在某些點上會產生間斷，引起重構信號的吉布斯振蕩[13]，因此本文選擇軟閾值進行處理。
| $\left\{ \begin{array}{l}硬閾值:w\lambda=\left\{ \begin{array}{l} w, \left| w \right| \ge \lambda \\ 0, \left| w \right| < \lambda \end{array} \right.\\軟閾值:w\lambda=\left\{ \begin{array}{l} \left[{{\rm{sign}}\left(w \right)} \right]\left({\left| w \right|-\lambda } \right), \left| w \right| \ge \lambda \\ 0, \left| w \right| < \lambda \end{array} \right. \end{array} \right.$ |
上式中，閾值λ起決定性作用，λ太小噪聲去除不干凈，λ太大則會將有用信號成分去除造成信號失真。常用的閾值計算方法有：固定閾值、基于極大極小原理閾值、自適應閾值。若采用全局固定閾值處理則忽略了各層系數間的差異，導致去噪效果不理想，因此本文采用基于無偏似然估計原理的自適應閾值，即史坦無偏估計(Stein’s unbiased risk estimation, SURE)。對于一個給定閾值，求出其對應的似然估計，在進行非似然最小化，從而得出所選閾值[14]。
設T為一向量，其元素為按從小到大順序排列的小波系數平方，即
| ${\rm{\boldsymbol{T}}}=\left\{ {{w_1}, {w_2}, \cdots, {w_n}} \right\}, $ |
其中w1≤w2…≤wn，再設一風險向量R，其元素如式(7)，求出R元素中最小值rj所對應的下標j，找出wj，則可求的閾值。
| ${r_i}=\left[{{\rm{n}}-2i + \left({n-i} \right){w_i} + \sum\limits_{k=1}^i {{w_k}} } \right]/{\rm{n}}\left({i=1, 2, \cdots, {\rm{n}}} \right)$ |
2.2.2 除噪及效果評價
根據上文步驟，對圖 3所示設備采集到的心阻抗微分信號進行去噪處理。圖 4為利用bior3.7小波對心阻抗微分信號進行6層分解后各尺度下的近似分量和細節分量，由圖看出高頻噪聲分量主要集中在d1、d2、d3上，有用信號集中于d4、d5上。除噪時，根據SURE選取合適的閾值對各層細節分量進行收縮處理，處理時盡可能減小富含噪聲的d1、d2、d3尺度上的小波系數，保留有用信號集中的d4、d5尺度上的系數，這樣利用收縮處理后的各層小波系數加之近似信號a6重構出信號即可達到消噪目的。
圖4
bior3.7小波6層分解后各尺度小波系數
Figure4.
Each scale wavelet coefficients after bior3.7 six layer wavelet decomposition
消噪結果如圖 5所示：(a) 為健康人心阻抗微分信號，與圖 1相比，波形標準、特征點明顯，存在少量噪聲。(b) 為一患者心阻抗微分信號，該患者為男性，57歲，患有冠心病、不穩定性心絞痛。由圖中看出信號受噪聲干擾嚴重，波形在O波處發生變異。(c) 為圖 5(a)除噪后信號，與(a)相比，波形更為平滑干凈，且波峰未被削減。(d) 為圖 5(b)除噪后信號，與(b)相比，除噪后信號變得平滑干凈，特征點凸顯，去噪效果十分明顯。
為了能定量地評價去噪效果，在此采用信噪比(signal to noise ratio, SNR)和均方誤差(mean square error, MSE)進行評價。SNR指有用信號與噪聲的比值，越大表示信號所含噪聲越少；MSE體現了原始信號和去噪后信號間的差異，越小表示去噪效果越好[12]。
| $\left\{ \begin{array}{l} {\rm{SNR}}=10{\log _{10}}\sum\limits_{i=1}^N {{{\left({\frac{{{{\rm{x}}_i}}}{{{{\rm{x}}_i}-{{\rm{y}}_i}}}} \right)}^2}} \\ {\rm{MSE}}=\sqrt {\frac{1}{N}\sum\limits_{i=1}^N {{{\left({{{\rm{y}}_i}-{{\rm{x}}_i}} \right)}^2}} } \end{array} \right.$ |
現實中，完全不含噪聲的心阻抗微分信號是不存在的，因此本文采用圖 5(c)中的信號近似為真實信號。圖 5(c)中信號為實驗室24歲健康男性采用小波閾值法除噪后的信號。式(8)中xi為圖 5(c)中信號，yi為醫院采集的患者去噪后信號，N表示信號長度。
表1
兩種方法去噪后6組信號的信噪比和均方誤差
Table1.
Six groups of signal SNR and MSE after using two de-noising methods
2.3 特征點檢測
準確定位B、C、X各特征點是正確計算SV、CO等血流動力學參數的基礎。信號經小波分解后，得到不同時間尺度上的幅度特性，在某些尺度上，有用信號的特征點十分明顯，同時噪聲被充分抑制，在這些尺度上檢測信號的模極大值對間過零點來定位特征點[7]。心阻抗微分信號的能量分布與主動脈瓣開閉的頻率(30～45 Hz)相近[8]，本文中采樣頻率為1 000 Hz, 根據奈奎斯特采樣定理，采樣率為1 000 Hz的信號最高頻率為500 Hz，6層分解后，尺度4、5集中包含了心阻抗信號的主要能量，高頻系數d5對應的頻帶為15.625～31.25 Hz，與阻抗信號的頻率相近，而肌電干擾等噪聲在S=25尺度下大大減小，結合大量實驗驗證，在S=25尺度下特征點檢測效果最好。
實驗室研制的KF_ICG型無創心功能檢測儀能夠實現信號特征點的自動檢測識別。信號采集完成后(或回放后)，點擊圖 3中的分析按鈕將會自動彈出圖 6所示的特征點檢測界面。圖 6中左下框內顯示了一個周期內完整的心阻抗微分波形，白色虛線即為各特征點檢測結果，由圖看出經小波變換分解后，在感興趣的尺度上成功定位出心阻抗微分信號的特征點。特征點定位完成后，界面右半部分計算出常用的血流動力學參數便于臨床醫生診斷。
圖6
儀器檢測分析界面
Figure6.
Detection and analysis interface of KF_ICG instrument
3 結果與討論
用采集到的實驗室健康人和22例患者數據進行特征點檢測驗證，對同一研究對象分別采用差分閾值法和小波變換法檢測心阻抗微分信號的特征點，最終統計結果如表 2所示。
表2
兩種方法對健康人和患者的C波、B點、X點的檢測結果
Table2.
C wave, B point and X point detection of healthy people and patients by using traditional differential threshold and wavelet threshold
表 2列出了傳統差分閾值法和小波變換法對健康人和患者阻抗微分信號中C波、B點、X點定位的準確率。由表 2得出以下結論：(1)小波法對健康人和患者的定位準確率均高于差分閾值法，但對健康信號二者相差不大。分析原因為健康人年輕(平均23.17歲)，心功能狀態較好，波形幅值大，特征點明顯不易錯檢，因此針對健康信號兩種方法的識別率均較高。(2)針對摻雜了諸多噪聲干擾的臨床信號，小波法的準確率明顯高于閾值法，特別是對幅值小、易被淹沒的B點、X點更為明顯。分析原因：①患者年齡較大(平均62.29歲)，心功能減弱，波形幅值小，B點、X點等特征點在波形圖中不夠明顯；②檢測時要求受試者在安靜的環境下放松靜臥、勻稱呼吸，健康人檢測地點為實驗室，能較好地滿足環境要求；臨床患者檢測時環境嘈雜，特別是某些患者呼吸急促，帶來幅值較大的干擾和噪聲使B點、X點被淹沒其中，而差分閾值法的抗噪能力較差，造成錯檢點增多；③某些患者可能存在影響阻抗檢測的病癥，導致檢測到的波形變異，致使特征點較難定位。但總體說來，小波變換法對健康人和患者的準確率均較高，均在97%以上，結合圖 3儀器檢測界面中列出的臨床常用心功能診斷參數SV、CO值，與重慶市大坪醫院心內科采用的超聲多普勒USCOM1A型儀器的測量值相對比，誤差均在10%以內，證明該法準確有效。
4 結束語
應用阻抗法檢測SV、CO等血流動力學參數是無創、經濟的。為了保證計算機對摻雜了諸多噪聲干擾的臨床數據也能精確地自動識別診斷信號特征點，前期高質量的降噪必不可少。傳統的頻域去噪對頻率波動大的非平穩信號顯得無能為力，而小波變換具有多分辨率分析功能，在生理信號去噪與特征點的檢測中顯示出獨有的優越性[15]。本文應用小波閾值法對采集到的信號進行消噪，并在消噪的基礎上對信號特征點進行自動識別檢測，形成一套完整的心功能信號處理系統。利用該系統采集了28例樣本進行驗證，結果表明該系統無論對健康人還是患者信號都有較好的抗噪識別能力。下一步工作是采集更多的臨床數據，針對某些波形發生變異的數據進行改進，提高設備對不同病癥患者個體差異的魯棒性，從而真正服務于臨床診斷應用。
）