方 震 白忠瑞 陳賢祥 夏 攀 何征嶺 趙榮建
①(中國科學(xué)院空天信息創(chuàng)新研究院 北京100190)
②(中國科學(xué)院大學(xué)電子電氣與通信工程學(xué)院 北京100049)
心沖擊圖(BallistoCardioGram,BCG)是一種無接觸式的心臟活動監(jiān)測手段。它來源于心臟泵血過程中血液流動在人體內(nèi)部產(chǎn)生的一系列機(jī)械沖擊力[1]。與心電圖(ElectroCatdioGram,ECG)等其他生理信息監(jiān)測技術(shù)相比,BCG避免了電極線和電極貼長時間貼在身上帶來的束縛和不適,具有無創(chuàng)、無直接接觸和檢測方便等優(yōu)勢,十分適用于居家養(yǎng)老、病房等需要長時間監(jiān)測生理信號的場景。近年來隨著傳感器和數(shù)字信號處理的發(fā)展,BCG逐漸受到研究人員的重視,研究表明,BCG可以應(yīng)用于心率、心率變異性(Heart Rate Variability,HRV)監(jiān)測[2,3],也可應(yīng)用于心臟收縮性以及心輸出量變化等心血管功能的評估[4]。
但BCG信號在不同測試條件(如測試設(shè)備的安置位置與角度、受試者姿勢等)下表現(xiàn)出較大的不一致性[5],此固有特征在一定程度上限制了分析BCG的形狀特征在評估心血管功能方面的應(yīng)用潛力。相比之下,通過檢測逐拍心動周期估計心率或者計算HRV指標(biāo),是BCG可行性較高的一種應(yīng)用路徑。
HRV定義為逐次心跳周期差異的變化情況,一般通過分析心電信號中R波間期序列獲得。HRV反映了心臟本身竇性心律不齊的程度以及神經(jīng)體液因素與竇房結(jié)之間相互作用的平衡關(guān)系。已有的研究表明,HRV是心源性猝死、冠心病、高血壓病及慢性心力衰竭等心血管疾病及慢性阻塞性肺疾病、糖尿病等疾病預(yù)后的預(yù)測因子[6,7],還能反映出睡眠狀況、精神壓力狀況等多種信息[8,9]。因此使用無感式的BCG信號獲取心率變異性指標(biāo),對于用戶方便舒適地掌握自身健康狀況具有重要意義。
通過BCG獲取HRV的關(guān)鍵在于獲得精準(zhǔn)可靠的逐拍心率。近年來,已有很多學(xué)者在非接觸式獲取逐拍心率方面做了研究:在BCG信號采集裝置方面,有壓電薄膜、壓電陶瓷、光纖、毫米波雷達(dá)、加速度計、壓變電阻以及各類采用MEMS工藝的微傳感器等[10–12];在逐拍心動周期的算法方面,目前有自適應(yīng)閾值法、峰值探測法、倒譜分析法、多示例學(xué)習(xí)法、模板匹配法、混合判決方法[13–16]。值得注意的是,目前此方面研究更多地關(guān)注逐拍心率檢測的覆蓋率和平均心率的準(zhǔn)確度,但事實上,逐拍心動周期計算的精度,即平均絕對誤差(Mean Absolute Error,MAE)對于心率變異性的計算也尤為重要。例如,在5 min短時HRV計算中,若BCG逐拍心動周期序列平均絕對誤差為10 ms,則由其計算得到的頻域特征HF相比ECG真實值的偏差可達(dá)15%~25%,SDNN,p NN50等時域指標(biāo)偏差也可達(dá)10%左右,且此偏差具有隨機(jī)性,這使得BCG所得到的HRV特征指標(biāo)的可靠性降低;而當(dāng)逐拍心動周期的MAE降至4 ms時,逐拍心動周期序列所得到的HF偏差降低到8%以內(nèi),SDNN和p NN50等指標(biāo)偏差降低至2%。
針對目前大多數(shù)方法在逐拍心率計算精度方面的不足,本文完成以下工作:(1)通過采用優(yōu)化的傳感前端設(shè)計,增加傳感器的靈敏度和BCG信號的時間分辨率;(2)通過對心沖擊信號進(jìn)行分析,找到BCG中最適合提取逐拍心動周期的關(guān)鍵成分;(3)提出一種采用聚類的自適應(yīng)模板匹配算法,以準(zhǔn)確提取心動周期信息。最后通過實驗驗證了方法的準(zhǔn)確性。
BCG信號采集裝置的設(shè)計影響逐拍心率提取的準(zhǔn)確度。在各種BCG監(jiān)測方式中,壓電陶瓷傳感器由于其靈敏度高、價格低廉、設(shè)計靈活等優(yōu)勢,具有良好的應(yīng)用前景。本文通過壓電陶瓷非接觸式采集BCG,計算逐拍心率。
本文所設(shè)計的傳感前端外殼為一個ABS材質(zhì)塑料圓盒,壓電陶瓷片被固定在圓盒內(nèi)上壁形變最大位置。外殼上蓋結(jié)構(gòu)設(shè)計為一側(cè)凸起形狀,使其受到震動時響應(yīng)更靈敏。實驗中將傳感前端放置在胸口下方的床墊下,采集到的身體震動信號會經(jīng)放大、濾波以及模數(shù)轉(zhuǎn)換傳輸?shù)綆в械凸乃{(lán)牙(Bluetooth Low Energy,BLE)功能的主控芯片中進(jìn)行下一步的傳輸和處理。BCG信號的采樣頻率也影響逐拍心率計算的精準(zhǔn)性,經(jīng)對比測試,本文設(shè)置BCG采樣頻率為250 Hz,以保證逐拍心率計算的準(zhǔn)確程度,同時最小化資源消耗。
傳感前端采集到的原始的壓電信號可用式(1)表示
其中,r(n)表示呼吸引起的振動成分,b(n)表示BCG成分,g(n)表示噪聲成分。心臟的跳動是一個較為復(fù)雜的機(jī)械運動過程[17],BCG中可以反映心臟機(jī)械運動和心臟泵血引起的身體機(jī)械運動的很多信息。
但對于一般所測得的BCG信號,由心臟原始的泵血引起的信號相對細(xì)微,容易被淹沒在由身體四肢震動造成的信號的主能量中,通過不同的處理方式可以凸顯出BCG中的不同成分。圖1是不同方式處理BCG信號的結(jié)果,虛線是同步采集的ECG的R波所在位置。圖1(a)是經(jīng)過48 Hz低通濾波器預(yù)處理去除高頻噪聲后的壓電信號,為了求得誤差最小的逐拍心動周期誤差,我們嘗試了多種不同的方式對此信號進(jìn)行進(jìn)一步處理,以獲取BCG信號中最適合進(jìn)行精準(zhǔn)逐拍心率提取的“關(guān)鍵成分”,下面優(yōu)選幾種具有代表性的方法具體說明:
方法1使用3~24 Hz的巴特沃斯帶通濾波器進(jìn)行濾波,得到BCG概貌波形,其結(jié)果如圖1(b)。
方法2使用與心沖擊信號最相符的db6小波基對原始壓電信號進(jìn)行9級小波分解,并使用小波分解后的第3~7層細(xì)節(jié)分量重建BCG信號,其結(jié)果如圖1(c)所示。
方法3使用8~24 Hz的巴特沃斯帶通濾波器進(jìn)行濾波,可以得到BCG信號中的由心臟跳動引起的細(xì)微振動成分,其結(jié)果如圖1(d)所示。
方法4對方法1處理的結(jié)果進(jìn)行二次差分操作,求取其加速度信號,結(jié)果如圖1(e)所示。
對于不同的BCG處理方法,每拍心跳在信號上造成的最高點并不是同一位置,為方便起見,本文中統(tǒng)稱這些點為BCG的“J點”。
2.3.1模板學(xué)習(xí)
圖1 BCG不同成分提取結(jié)果
BCG形狀多變,難以提前預(yù)知其波形。我們使用了親和傳播(Affinity Propagation,AP)聚類算法[18],從10 s的BCG信號中進(jìn)行模板的自學(xué)習(xí)。具體步驟如下:
(1)提取2.2節(jié)各方法處理過的10秒BCG信號段中所有的極大值點,以這些極大值點為中心,截取其前后各100采樣點,共約0.8 s長度的信號,得到數(shù)據(jù)樣本集B={b0b1···bn},其中bi代表第i個極值點處截取的信號。
(2)使用AP聚類算法對樣本集B進(jìn)行聚類,并得到模板聚類。AP聚類的具體步驟如下:
步驟1初始化吸引度矩陣R和歸屬度矩陣A為0矩陣。其中R中元素r(i,k)表示數(shù)據(jù)bk適合作為數(shù)據(jù)bi的聚類中心的程度,A的元素a(i,k)表示對象bi適合選擇數(shù)據(jù)對象bk作為其聚類中心的程度。
步驟2分別按照式(2)和式(3)迭代更新吸引度矩陣與歸屬度矩陣。
其中,s(i,k)是相似度矩陣S的元素,為弱化BCG幅值變化的影響,本文沒有使用歐氏距離,而是設(shè)置s(i,k)為Pearson相關(guān)系數(shù)與常數(shù)1的差,如式(4)所示。特殊地,i=k時的相似度s(k,k)稱為參考度,影響聚類的數(shù)目,本文將s(k,k)設(shè)置為相似度的中位數(shù)M(S),可有效完成聚類。最后得到的相似度矩陣S如式(5)所示。
步驟4重復(fù)步驟2、步驟3,直到矩陣穩(wěn)定。按最大相似度規(guī)則將其他點分配到距離它最近的聚類中心點相應(yīng)的聚類中,完成聚類。
(3)對各聚類,計算每個樣本的中點值的平均值,選擇此平均值最大的聚類為理想聚類,計算其算術(shù)平均作為最后的BCG模板,以BCG表示,如圖2所示。
2.3.2心跳位置探測
通過將模板在濾波后BCG數(shù)據(jù)上不斷右滑,形成相關(guān)系數(shù)函數(shù)。對于一段預(yù)處理后的BCG信號,采用以下步驟進(jìn)行模板匹配,確定J點。
(1)通過以下方式構(gòu)造原BCG信號等長的相關(guān)系數(shù)函數(shù)cor(x)。
(2)對前10 s BCG信號做FFT變換,然后進(jìn)行頻域?qū)し宕_定心率近似值HRe。在最初的2 s的cor(x)信號中選取最大值,作為第1個J點的位置,以此為起點,向后搜索[60/(HRe+20)s,60/(HRe–20)s]范圍內(nèi)的最大值點作為下一個J點。當(dāng)連續(xù)5個所選局部最大值小于閾值θ時,說明匹配質(zhì)量下降,重新計算模板,本文中閾值θ設(shè)為0.75。
最后,根據(jù)式(9)和式(10)計算連續(xù)逐拍心動周期和逐拍心率。
其中,IndexJ(n)表示第n個J點的索引位置,F(xiàn)s為采樣頻率。逐拍心動周期和逐拍心率之間可方便地轉(zhuǎn)換,為更清楚地比較計算結(jié)果,本文使用逐拍心動周期進(jìn)行誤差分析。
按照如圖3所示方法進(jìn)行實驗。本文使用壓電陶瓷設(shè)備采集BCG信號,同時為驗證BCG逐拍心率計算的準(zhǔn)確性,在本研究中同步采集單導(dǎo)聯(lián)心電信號進(jìn)行對比。實驗開始時按單導(dǎo)聯(lián)方式在受試者胸腹部佩戴一次性心電電極貼,隨后讓受試者安靜平躺于實驗床墊上,開啟設(shè)備進(jìn)行BCG與ECG的同步采集。通過2.2節(jié)和2.3節(jié)所述方法處理BCG和計算BCG的逐拍心動周期,再使用PT算法[19]定位ECG的R波,計算ECG逐拍心動周期。隨后分析對于2.2節(jié)所述的不同的BCG成分,所提取出的逐拍J-J間期相對于R-R間期的覆蓋率和精準(zhǔn)度。得到每次心跳的R-R間期與J-J間期。所有數(shù)據(jù)分析使用Py Charm與Python3.5進(jìn)行。
圖4(a)和圖4(b)分別展示了受試者1的一段BCG信號在3~24 Hz帶通濾波和8~24 Hz帶通濾波的預(yù)處理方法下,使用本文第2節(jié)所述方法進(jìn)行逐拍心率提取的結(jié)果,圖中使用紅色圓圈標(biāo)出了所定位的J點的位置,虛線位置為ECG的R波位置。圖5則展示了兩種預(yù)處理方式所計算得到的逐拍心動周期與ECG逐拍R-R間期的對比圖??梢娊?jīng)通帶為3~24 Hz的濾波器預(yù)處理之后計算得到的J點定位在原始信號上即為每次心跳所引起的壓電信號最大點位置,但如圖5(a)所示,由此所得到的逐拍心率卻與逐拍R-R間期有較大差別。而經(jīng)通帶為8~24 Hz的濾波器處理后,所提取的J點位置不再是原始壓電信號每拍心跳的最大值,但卻能得到與R-R間期更為吻合的J-J間期,如圖5(b)所示。我們稱這種現(xiàn)象為“關(guān)鍵成分”現(xiàn)象,即8~24 Hz的成分是去除了較低頻的肢體振動等成分、更能反映心室泵血、血液沖擊主動脈引起的心臟原生振動的“關(guān)鍵成分”,能夠獲得與R-R間期相近的J-J間期。
圖2 BCG模板聚類與BCG模板
圖3 非接觸式逐拍心動周期檢測準(zhǔn)確度驗證實驗示意圖
圖4 不同頻帶濾波的BCG模板匹配計算結(jié)果
表1為對受試者1的5 min BCG數(shù)據(jù),共計351次心跳統(tǒng)計的結(jié)果,心跳覆蓋率定義為使用本文算法成功定位到的心跳(即定位J波與R波相隔50 ms內(nèi))的比例,平均絕對誤差定義為逐拍J-J間期與RR間期之差的絕對值的平均值,而平均誤差定義為逐拍J-J間期相對于R-R間期的平均誤差。可見對于受試者1的5 min安靜平躺的BCG數(shù)據(jù),在分別采用4種預(yù)處理方式后,3~24 Hz濾波以及db6小波重構(gòu)都取得了100%的心跳檢出率,但卻有較高的平均絕對誤差。經(jīng)8~24 Hz濾波和二次差分處理后,雖心跳檢出率略有不足,但卻取得了很低的平均絕對誤差。
圖6展示了兩種不同預(yù)處理方式下受試者1的逐拍R-R間期與J-J間期Bland-Altman圖,從中可以看出采用寬帶寬濾波器所計算的結(jié)果(圖6(a))在所有心跳間期分布范圍內(nèi)誤差都明顯高于窄帶濾波器計算結(jié)果(圖6(b))。兩種預(yù)處理方式的95%置信區(qū)間分別為(–32.26 ms,32.26 ms)和(–12.34 ms,12.34 ms),且兩種處理方式所得逐拍心跳間期誤差均值都為0,說明在給定時間段內(nèi)使用BCG和ECG檢測出了相同的心跳次數(shù)。此處需要說明的是,雖表1中8~24 Hz濾波器的處理方法心跳覆蓋率低于100%,但這僅說明BCG有心跳檢測位置與R波位置相差較多,總體依然是檢測到了應(yīng)有的心跳次數(shù)。
本文采集了15名受試者的BCG和ECG數(shù)據(jù),年齡23~34歲,身高156~182 cm,體重45~90 kg。表2展示了所有受試者各5 min安靜平躺狀態(tài)下采集到的數(shù)據(jù)的處理結(jié)果,并且展示了4種BCG成分對逐拍心率提取性能的比較。
由表2可見對于8~24 Hz的帶通濾波和二次差分這兩種方式處理,在檢出心跳覆蓋率方面略有不足,這是因為頻率較高時波形表現(xiàn)得較為復(fù)雜,影響了模板匹配的效果;而同時這兩種方式在平均絕對誤差(MAE)方面有著較大的優(yōu)勢,這是因為高頻成分更能反映心臟跳動的成分而不受整個身體震動的影響。
圖5 ECG-RR間期與不同頻率濾波的逐拍BCG-JJ間期折線圖
表1 受試者1的數(shù)據(jù)處理結(jié)果對比
圖7展示了2.1節(jié)所述各成分所計算得到的J-J間期與R-R間期誤差絕對值的分布,可以看到對于前兩種處理方式,75%的誤差分別在20 ms和16 ms以下,而后兩種處理方式的誤差有75%在4 ms以下,且中位數(shù)、下四分位數(shù)、下限均為0 ms,說明有1/2以上的J-J間期具有與R-R間期相同的采樣點數(shù)。
圖6 受試者1的逐拍RR間期與JJ間期Bland-Altman圖
表2 不同預(yù)處理方式下逐拍心率提取性能
圖7 各成分誤差分布箱型圖
此外,本文還從系統(tǒng)整體的角度對比了同類研究工作在靜息狀態(tài)下測量BCG提取逐拍心率的主要性能指標(biāo),如表3所示??梢姡疚姆椒軌蛟诟采w率和逐拍心動周期的準(zhǔn)確度方面處于領(lǐng)先水平,尤其是本文MAE可達(dá)到4 ms以下,這有助于提高BCG計算HRV的可靠性。
本文提出了一種基于壓電陶瓷傳感器采集BCG并提取精準(zhǔn)逐拍心率的方法,該方法是通過優(yōu)化傳感器性能提高BCG信號質(zhì)量,通過對比BCG信號的不同處理方式,找到最適合提取精準(zhǔn)逐拍心率的成分,最后通過AP聚類和模板匹配的方法確定J點,計算逐拍心率。
表3 BCG逐拍心率提取方法性能對比
本文的逐拍心動周期平均絕對誤差為3.78 ms,心跳檢測覆蓋率為98.5%。本研究結(jié)果表明,通過設(shè)計合適的BCG采集裝置、采用合適的信號處理方法處理BCG信號以及選擇合適的特征點提取算法,可以使BCG提取逐拍心動周期(J-J間期)相對于R-R間期的誤差降至可接受的范圍,這對于提高使用BCG進(jìn)行健康監(jiān)測的效果尤其是提高HRV計算的可靠性具有重要意義。對本文所提出的BCG“關(guān)鍵成分”成因進(jìn)行實驗仿真和驗證,以及對BCG提取逐拍心率精度的影響因素作進(jìn)一步量化分析將是下一步的研究方向。