本文為解決傳統心沖擊圖（BCG）提取心率不足的問題，提出了一種改進的心沖擊圖檢測心率的方法，首先通過嵌入式傳感器實時采集微弱的心臟活動信號，對信號濾波、轉換獲得局部 BCG 節拍，然后不經由心跳模板而直接由 BCG 節拍估計心率。與其他方法相比，本文所提方法在心率數據準確性和抗干擾性方面具有較強優勢，也實現了非接觸式的在線檢測，最終通過分析 13 名受試者約 20 000 次的心率數據，達到了平均節拍誤差為 0.86%，覆蓋率 96.71% 的滿意結果，或可為醫院臨床和居家養老提供一種估計心率的新方法。
引用本文： 任彧, 劉穩, 高志剛. 一種改進的用于心率估計的峰值提取方法. 生物醫學工程學雜志, 2019, 36(5): 834-840. doi: 10.7507/1001-5515.201810041 復制
引言
目前心臟類疾病是導致人類死亡的常見原因之一，如何預防和發現這類疾病是人們最關心的問題，而檢測和診斷的關鍵是如何得到長時間心臟活動的數據。傳統 24 h 心電圖（electrocardiogram，ECG）可以獲取一段時間內的心臟活動數據，但很難找到便捷且費用適宜的儀器進行較長時間的心臟活動數據監測和提取。本文使用嵌入式傳感器可以長時間、實時監視受試者的心率，給醫院臨床和居家養老提供了一種長時間獲取心率數據并且有助于心臟類疾病在線檢測的方法。
心沖擊圖（ballistocardiogram，BCG）是一種記錄心臟活動引起的身體震動的方法[1]，將其用于測試心率時，相比較目前常用的 ECG 信號而言，BCG 信號檢測不需要將電極附著到受試者身體上，不易受到受試者的心情和情緒影響，且可在受試者不知情的狀況下進行測試并記錄數據，通過對大量心率數據的特征提取和匹配，有利于早期發現心臟類疾病，對臨床診斷有一定的參考價值。
近年來，已有學者研究了利用 BCG 數據進行心率統計和預測[2-4]。例如：Brüser 等[2]利用聚類方法為給定的 BCG 數據記錄確定合適的心跳模板，并利用該模板估計隨后的 BCG 節拍。當模板匹配一致時，該方法能獲得較好的結果，但缺陷是模板的建立嚴重依賴于訓練階段，如果算法在無監督訓練期間學習了錯誤的模版，則所有后續心跳檢測都將失敗。此外，當受試者的位置移動、肢體狀態發生變化時，會導致心跳模板發生變化，需要重新訓練模板。
對具有標準波形的 BCG 數據來說，可以通過傳統的峰值檢測方法來識別 J 峰[5]，如圖 1 所示。但 BCG 信號提取時易受到其他信號的干擾，干擾源包括：運動偽影、BCG 信號振蕩和工頻干擾等。因此使用 BCG 數據進行心率提取時需要解決：① 外界因素造成的噪聲：由于受試者與傳感器非直接接觸的方式，BCG 信號非常微弱，當受到外界噪聲干擾時更難以確定 J 峰的確切位置，如圖 1 所示，藍色波形表示了受試者的一個心跳周期在受到外界噪聲（如工頻噪聲）干擾時的 BCG 波形。② 人為因素干擾：在測試期間，如果受試者肢體狀態變化造成其與傳感器之間的接觸不穩定，會產生嚴重的噪聲，即偽影。這些偽影能通過數據分析來獲取，也可以利用數據過濾進行消除，但會丟失一段時間的 BCG 數據。因此，傳統峰值檢測算法在這些情況下無法準確地提取心率參數。如圖 1 所示，黑色線條波形為受試者肢體狀態改變時 BCG 信號振幅的變化。
圖1
標準 BCG 波形、非標準 BCG 波形和噪聲
Figure1.
Standard BCG waveform, non-standard BCG waveform and noise
為了提高心率測量的準確性，研究人員們提出了許多方法來消除運動引起的噪聲，其中包括：標準數字濾波器[6]、線性相位濾波器[7]和小波變換[8]等。然而，現有大部分方法是針對呼吸或細微的身體運動等產生的微弱噪聲設計，當采集到的傳感器信號含有大量的肢體活動信號時，BCG 信號成分變得微弱，這些算法都會表現不佳。盡管上述方法在一定程度上能夠消除部分噪聲，但需增加復雜的濾波、變換算法對信號進行處理，這樣會提高對嵌入式儀器設備的硬件資源要求，導致成本增加，并不適合進一步推廣應用。
針對以上問題，本文提出了一種新型的峰值檢測方法，該算法不需要提前獲得心跳模板，不需要定義檢測匹配模板波形，直接從信號數據中估計心跳間隔并獲得心率數據，以期實現在有局部干擾的 BCG 信號中識別有效的心率數據。另外，本文也設計了在疊加偽影的 BCG 信號中提取心率數據的方法，該方法通過信號濾波處理和構造偽影模型來估計 BCG 信號在偽影段的 BCG 節拍，以提高在存有偽影的 BCG 信號中提取心率數據的準確率。
1 方法
1.1 數據來源
本次實驗分析使用的數據來自本課題組招募的兩個年齡段的 13 名志愿者：① 18～25 歲，8 男、2 女，編號為：1～10；② 55～65 歲，3 男，編號為：11～13。其中，編號為 13 的受試者患有心臟病和高血壓，其余受試者身體健康，無心臟相關疾病史，均自愿參加本研究測試，并簽署了知情同意書和承諾書。本文數據采集過程均未違反任何倫理道德準則。
每位受試者均使用 3 導聯 ECG 信號監測儀（AD8232，Analog Devices，美國）采集其 ECG 信號，并經頻率為 1 000 Hz 的二階高通濾波器濾波后，提取 R-R 間期數據，將該數據作為金標準以備參考[9]。
數據采集試驗于杭州電子科技大學杭電 TI 聯合實驗室進行，將高精度單軸傾角傳感器（SCA61T-FAHH1G，VTI，日本）［靈敏度：70 mV/（°），后文簡稱：傾角傳感器］安裝于普通床墊底部靠近人體胸部的位置，以便記錄心臟的振動和呼吸運動。
如圖 2 所示，受試者呈仰臥位，使用嵌入在床墊底部的傾角傳感器獲取 BCG 數據，數據通過采樣率為 1 kHz 的 12 位模數轉換器獲取。傾角傳感器所獲取的信號表示受試者體內血液置換對傾角傳感器的逐搏反作用力，并且它與人體生理活動產生的沿水平方向的作用力成正比。
圖2
BCG 系統
Figure2.
System of BCG
1.2 預處理
受試者對傾角傳感器的逐搏反作用力經采集系統轉換為一個帶有直流偏移的交流電壓輸出。首先對傾角傳感器采集的原始信號進行高通濾波，截止頻率為 3 Hz，去除由床和身體產生的直流分量，經差分放大后進行 50 Hz 陷波和低通濾波，截止頻率為 11 Hz，去除工頻干擾和高頻噪聲，最終得到 BCG 信號，如圖 2 所示。
1.3 基本算法
本文采用的算法為局部相鄰波峰和波谷的最大差值檢測，對應于單個心臟跳動的 BCG 信號峰值模式的心跳檢測算法（以下簡稱：改進的峰值檢測方法）。設預處理后的信號表示為 x［n］，在連續的信號上迭代移動分析窗口，提取特征信息。如式（1）所示：
 |
其中，HR 為心率，Thi 為心跳間隔。因人體正常心率范圍為 40～140 次/min，由式（1）可得心跳間隔的最大值和最小值，定義心跳間隔區間為［Tmin，Tmax］。
將分析窗口區間定為［Ni，Ni + 2·Tmin·Fs］，以確保在此區間內至少檢測到一次心跳，其中 Ni 表示樣本起始位置，Fs 表示采樣頻率。算法實現分為 6 個步驟，如下所示，在步驟（1）和步驟（2）中，對窗口 x［Ni，Ni + 2·Tmin·Fs］進行 2·Tmin·Fs 次迭代。
步驟（1）：找出區間內的所有波峰 P，即極大值并記錄其索引，如式（2）所示：
 |
其中，val 表示波峰的幅度大小，idx 表示波峰在 x［n］出現的索引。
找出區間內的所有波谷 T，即極小值并記錄其索引，如式（3）所示：
 '/> |
其中，val′表示波谷的幅度大小，idx′表示波谷在 x［n］出現的索引。
步驟（2）：計算分析窗口的最大振幅 
，并檢查 
是否在特定的閾值范圍內。偽影的最高幅值遠超出正常 BCG 信號幅值，因此，如果 
位于閾值內，則此分析窗口中的信號被認為是可用的 BCG 信號。否則，視其為無效 BCG 信號，算法跳到步驟（5）。
步驟（3）：設結構體數組 P 和 T 的長度為 M，對 P 和 T 進行 M 次迭代。
找出相鄰波峰和波谷最大差值 Max 并記錄其索引，如式（4）所示：
 '/> |
其中，j 表示數組中的索引位置，val″ 表示幅值，idx″ 表示為在 x［n］中的索引。
步驟（4）：計算 Max［j］.idx″ 與 Max［j ? 1］.idx″ 的差值，即 BCG 節拍，判斷其是否在區間［Tmin，Tmax］范圍內。如果是，算法跳到步驟（6）。如果不是，視其為無效 BCG 節拍，算法繼續執行步驟（5）。
步驟（5）：當受試者翻動或側身時，通過傳感器采集到的肢體活動信號已經覆蓋大部分 BCG 信號，信噪比顯著降低，使得幾秒鐘內不可能進行可靠的心率估計。而在靜息狀態條件下，在短暫時間內心率變化不會太大，在大多數情況下連續的心臟跳動存在一定規律。本算法利用這些特征較準確地估計了在偽影下的 BCG 節拍，其目的是在數據不可用時將測量誤差盡量降低到最小。
如圖 3 所示，圖中紅點代表在 3 min 靜息條件下，檢測到的隨機選取的編號為 1 的受試者其連續標準 ECG 信號的心跳間隔。將其繪制成曲線并將其平滑處理，如圖中黑色波形線所示，可以觀察到，在短時間內（大約 10～20 次心跳間隔），曲線近似于一條微微傾斜的直線和一條含絕對值的正弦曲線的疊加，如式（5）所示：
圖3
在靜息下的心跳間隔
Figure3.
Heartbeat interval at rest
 |
對于無效窗口，本文記錄最近的 20 次 BCG 節拍數據，并近似地擬合式（5）的曲線，跟據此曲線實現對無效 BCG 信號的 BCG 節拍估計。
步驟（6）：由步驟（4）或步驟（5）估計得到 BCG 節拍，由式（1）計算心率參數 HR。
2 實驗驗證
本文先使用傾角傳感器采集 BCG 原始信號，然后對信號濾波、整合，最后利用商業數學軟件 MATLAB 2016a（MathWorks inc，美國）分析處理過的 BCG 信號。
首先根據所提方法計算得出 J-J 間期（即 BCG 信號中相鄰兩個 J 峰之間的時間間隔），再由 Hamilton 等[9]的方法獲取標準的 R-R 間期（即 ECG 信號中相鄰的兩個 R 波之間的時間間隔），其次分別統計 J-J 間期和 R-R 間期的數量（分別以符號 J?J 和 R?R 表示），并計算覆蓋率（即 J-J 間期與 R-R 間期的數量百分比，以符號 CR 表示），統計假陽性（即 J 峰相對 R 波的延時時間超過 100 ms）和假陰性（即未檢測到的 J 峰數量），再計算假陽率（以符號 FP 表示）和假陰率（以符號 FN 表示），然后計算 J-J 間期相對 R-R 間期的平均逐拍間隔誤差（表示為：平均誤差 E）和絕對平均間隔誤差（表示為：絕對誤差 Eabs）。最后，通過這些計算量和誤差統計來評估算法的性能。
2.1 測試示例
如圖 4 所示，為隨機選取的編號 1 的受試者（以下簡稱為受試者 1）部分測試示例，圖中展示了藍框和紅框內的放大數據，綠色波形上明顯突出的部分波形代表受試者在有肢體活動時產生的偽影，紅點表明了算法匹配的心跳位置，綠色波形上重疊的藍色波形表明為所提方法對有偽影疊加的 BCG 數據模擬重建。由于 BCG 信號相對 ECG 信號有一定的時間延遲，為了實現數據同步和消除計算誤差，本文也設定了 0～200 ms 的 BCG 信號和 ECG 信號的延遲補償時間。
圖4
受試者 1 的部分測試數據
Figure4.
A partial data of the subject 1
如圖 5 所示，展示了受試者 1 的 J-J 間期和 R-R 間期，可以看出大部分情況下，藍色波形和紅色波形幾乎重合在一起，而藍色波形中部分突出的毛剌可認為是所估計偽影段的 J-J 間期；圖 5 還展示了受試者 1 逐搏間隔誤差的 Bland-Altman 圖[10]，從中可以看出 99% 的間隔誤差位于 ? 71～71 ms 之間，95% 的間隔誤差位于 ? 15～15 ms 之間；子圖平均心率圖展示的是受試者 1 平均心率的變化（即按分鐘統計心跳數量），所估計的平均心率與參考 ECG 信號的誤差范圍在 ? 2～2 beat/min 內。總的來說，從圖 5 可以觀察到 BCG 信號和 ECG 信號逐搏間隔之間具有良好的一致性。
圖5
受試者 1 的 J-J 間隔和 R-R 間隔、Bland-Altman 圖和平均心率
Figure5.
J-J interval and R-R interval, Bland-Altman plot and mean heart rate of Subject 1
2.2 數據分析
如表 1 所示，表中記錄了改進的峰值檢測方法的結果。假陰率和覆蓋率是相對的，隨著不同受試者疊加的偽影長度變化，錯過的假陰性隨之發生相應變化，覆蓋率也會發生改變。表 1 的性能數據代表不同受試者之間的平均誤差、絕對誤差、覆蓋率、假陰率和假陽率。算法對 BCG 信號疊加偽影的敏感度會影響覆蓋率。當受試者在床上有大量肢體活動時，其測試結果不是很理想，原因是：由于考慮到要保持與參考數據的同步，該方法只適用于處理含短暫偽影的 BCG 數據，對于較長偽影的 BCG 數據，不能正確地進行 BCG 節拍估計。通常情況下，該方法可以準確估計大部分 BCG 信號的心率，但如果由于 BCG 信號失真等原因，導致算法估計的當前心跳信號位置和實際位置相差較大時，估計的 J-J 間期發生偏移，從而產生累加誤差，影響下次 J-J 間期的估計，進而假陽率增加，其結果不可靠。但重要的是，該方法的誤差水平仍然維持在較低的范圍內，而且其絕對誤差也保持在可接受的限度內。
表1
改進的峰值檢測方法的性能
Table1.
The performance of improved peak detection method
此外，本文還使用相同的數據集評估了 Choi 等[5]所提出傳統峰值檢測方法的性能，如表 2 所示，和 Brüser 等[2]所提出層次聚類算法的性能，如表 3 所示。Choi 等[5]提出的方法可以在一段時間內找到 4 個子區間內的局部峰值，并從這些局部峰值中選擇最大峰值作為 J 峰，同時還附加了一些額外的規則；Brüser 等[2]使用層次聚類，先從記錄中推斷出心跳模板，然后根據該模板選擇出信號中最適的位置作為心跳位置，進而確定心跳間隔。從平均假陰率來看，改進的峰值檢測方法為 3.29%，而忽略掉偽影的 Choi 等[5]為 5.92%，表明偽影占比大概在 5.92% 左右，這說明了改進的峰值檢測方法通過重建部分偽影數據用于心跳間隔的估計。從平均假陽率來看，Brüser 等[2]幾乎為 0，但平均覆蓋率只有 89.39%，可以看出這種對波形要求較高的算法，提高了對心跳位置的識別率而降低了覆蓋率。改進的峰值檢測方法的覆蓋率、平均假陽率、平均誤差和絕對誤差分別為 96.71%、1.01%、0.86% 和 6.79 ms，這表明與 Choi 等[5]和 Brüser 等[2]相比，該方法獲得的 J-J 間期偏差更小，具有更高的覆蓋率和更好的整體準確性。總體來說，改進的峰值檢測方法要優于 Choi 等[5]和 Brüser 等[2]的方法。
表2
傳統 J 峰檢測算法的性能
Table2.
Performance of traditional J-peak detection algorithms
表3
層次聚類算法的性能
Table3.
Performance of hierarchical clustering algorithm
3 結束語
本文提出了一種改進 BCG 信號估計心率的方法，該方法是利用 BCG 信號的 J 峰提取心跳節拍，通過尋找局部 BCG 信號中相鄰波峰、波谷的最大差值的位置，將此位置設定為 J 峰，最后獲得 BCG 節拍。對于有偽影疊加的 BCG 數據段，本文設計了一種偽影模型可以近似估計該數據段的 BCG 節拍。
從 13 名受試者獲得 260 min 測試數據，得到約 20 000 次心跳的 BCG 信號樣本，本文對方法進行了驗證，可以看到實驗的數據輸出與黃金標準參考[9]之間具有良好一致性。實驗結果還表明，即使存在復雜的偽影和其他弱噪聲的情況下，該方法也能保證較低的誤差水平和較好的抗干擾能力。同時，該方法也與已有的其他方法進行了比較，實驗結果同樣具有較高的覆蓋率和較好的準確率。
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
引言
目前心臟類疾病是導致人類死亡的常見原因之一，如何預防和發現這類疾病是人們最關心的問題，而檢測和診斷的關鍵是如何得到長時間心臟活動的數據。傳統 24 h 心電圖（electrocardiogram，ECG）可以獲取一段時間內的心臟活動數據，但很難找到便捷且費用適宜的儀器進行較長時間的心臟活動數據監測和提取。本文使用嵌入式傳感器可以長時間、實時監視受試者的心率，給醫院臨床和居家養老提供了一種長時間獲取心率數據并且有助于心臟類疾病在線檢測的方法。
心沖擊圖（ballistocardiogram，BCG）是一種記錄心臟活動引起的身體震動的方法[1]，將其用于測試心率時，相比較目前常用的 ECG 信號而言，BCG 信號檢測不需要將電極附著到受試者身體上，不易受到受試者的心情和情緒影響，且可在受試者不知情的狀況下進行測試并記錄數據，通過對大量心率數據的特征提取和匹配，有利于早期發現心臟類疾病，對臨床診斷有一定的參考價值。
近年來，已有學者研究了利用 BCG 數據進行心率統計和預測[2-4]。例如：Brüser 等[2]利用聚類方法為給定的 BCG 數據記錄確定合適的心跳模板，并利用該模板估計隨后的 BCG 節拍。當模板匹配一致時，該方法能獲得較好的結果，但缺陷是模板的建立嚴重依賴于訓練階段，如果算法在無監督訓練期間學習了錯誤的模版，則所有后續心跳檢測都將失敗。此外，當受試者的位置移動、肢體狀態發生變化時，會導致心跳模板發生變化，需要重新訓練模板。
對具有標準波形的 BCG 數據來說，可以通過傳統的峰值檢測方法來識別 J 峰[5]，如圖 1 所示。但 BCG 信號提取時易受到其他信號的干擾，干擾源包括：運動偽影、BCG 信號振蕩和工頻干擾等。因此使用 BCG 數據進行心率提取時需要解決：① 外界因素造成的噪聲：由于受試者與傳感器非直接接觸的方式，BCG 信號非常微弱，當受到外界噪聲干擾時更難以確定 J 峰的確切位置，如圖 1 所示，藍色波形表示了受試者的一個心跳周期在受到外界噪聲（如工頻噪聲）干擾時的 BCG 波形。② 人為因素干擾：在測試期間，如果受試者肢體狀態變化造成其與傳感器之間的接觸不穩定，會產生嚴重的噪聲，即偽影。這些偽影能通過數據分析來獲取，也可以利用數據過濾進行消除，但會丟失一段時間的 BCG 數據。因此，傳統峰值檢測算法在這些情況下無法準確地提取心率參數。如圖 1 所示，黑色線條波形為受試者肢體狀態改變時 BCG 信號振幅的變化。
圖1
標準 BCG 波形、非標準 BCG 波形和噪聲
Figure1.
Standard BCG waveform, non-standard BCG waveform and noise
為了提高心率測量的準確性，研究人員們提出了許多方法來消除運動引起的噪聲，其中包括：標準數字濾波器[6]、線性相位濾波器[7]和小波變換[8]等。然而，現有大部分方法是針對呼吸或細微的身體運動等產生的微弱噪聲設計，當采集到的傳感器信號含有大量的肢體活動信號時，BCG 信號成分變得微弱，這些算法都會表現不佳。盡管上述方法在一定程度上能夠消除部分噪聲，但需增加復雜的濾波、變換算法對信號進行處理，這樣會提高對嵌入式儀器設備的硬件資源要求，導致成本增加，并不適合進一步推廣應用。
針對以上問題，本文提出了一種新型的峰值檢測方法，該算法不需要提前獲得心跳模板，不需要定義檢測匹配模板波形，直接從信號數據中估計心跳間隔并獲得心率數據，以期實現在有局部干擾的 BCG 信號中識別有效的心率數據。另外，本文也設計了在疊加偽影的 BCG 信號中提取心率數據的方法，該方法通過信號濾波處理和構造偽影模型來估計 BCG 信號在偽影段的 BCG 節拍，以提高在存有偽影的 BCG 信號中提取心率數據的準確率。
1 方法
1.1 數據來源
本次實驗分析使用的數據來自本課題組招募的兩個年齡段的 13 名志愿者：① 18～25 歲，8 男、2 女，編號為：1～10；② 55～65 歲，3 男，編號為：11～13。其中，編號為 13 的受試者患有心臟病和高血壓，其余受試者身體健康，無心臟相關疾病史，均自愿參加本研究測試，并簽署了知情同意書和承諾書。本文數據采集過程均未違反任何倫理道德準則。
每位受試者均使用 3 導聯 ECG 信號監測儀（AD8232，Analog Devices，美國）采集其 ECG 信號，并經頻率為 1 000 Hz 的二階高通濾波器濾波后，提取 R-R 間期數據，將該數據作為金標準以備參考[9]。
數據采集試驗于杭州電子科技大學杭電 TI 聯合實驗室進行，將高精度單軸傾角傳感器（SCA61T-FAHH1G，VTI，日本）［靈敏度：70 mV/（°），后文簡稱：傾角傳感器］安裝于普通床墊底部靠近人體胸部的位置，以便記錄心臟的振動和呼吸運動。
如圖 2 所示，受試者呈仰臥位，使用嵌入在床墊底部的傾角傳感器獲取 BCG 數據，數據通過采樣率為 1 kHz 的 12 位模數轉換器獲取。傾角傳感器所獲取的信號表示受試者體內血液置換對傾角傳感器的逐搏反作用力，并且它與人體生理活動產生的沿水平方向的作用力成正比。
圖2
BCG 系統
Figure2.
System of BCG
1.2 預處理
受試者對傾角傳感器的逐搏反作用力經采集系統轉換為一個帶有直流偏移的交流電壓輸出。首先對傾角傳感器采集的原始信號進行高通濾波，截止頻率為 3 Hz，去除由床和身體產生的直流分量，經差分放大后進行 50 Hz 陷波和低通濾波，截止頻率為 11 Hz，去除工頻干擾和高頻噪聲，最終得到 BCG 信號，如圖 2 所示。
1.3 基本算法
本文采用的算法為局部相鄰波峰和波谷的最大差值檢測，對應于單個心臟跳動的 BCG 信號峰值模式的心跳檢測算法（以下簡稱：改進的峰值檢測方法）。設預處理后的信號表示為 x［n］，在連續的信號上迭代移動分析窗口，提取特征信息。如式（1）所示：
|
其中，HR 為心率，Thi 為心跳間隔。因人體正常心率范圍為 40～140 次/min，由式（1）可得心跳間隔的最大值和最小值，定義心跳間隔區間為［Tmin，Tmax］。
將分析窗口區間定為［Ni，Ni + 2·Tmin·Fs］，以確保在此區間內至少檢測到一次心跳，其中 Ni 表示樣本起始位置，Fs 表示采樣頻率。算法實現分為 6 個步驟，如下所示，在步驟（1）和步驟（2）中，對窗口 x［Ni，Ni + 2·Tmin·Fs］進行 2·Tmin·Fs 次迭代。
步驟（1）：找出區間內的所有波峰 P，即極大值并記錄其索引，如式（2）所示：
|
其中，val 表示波峰的幅度大小，idx 表示波峰在 x［n］出現的索引。
找出區間內的所有波谷 T，即極小值并記錄其索引，如式（3）所示：
| '/> |
其中，val′表示波谷的幅度大小，idx′表示波谷在 x［n］出現的索引。
步驟（2）：計算分析窗口的最大振幅 ，并檢查 是否在特定的閾值范圍內。偽影的最高幅值遠超出正常 BCG 信號幅值，因此，如果 位于閾值內，則此分析窗口中的信號被認為是可用的 BCG 信號。否則，視其為無效 BCG 信號，算法跳到步驟（5）。
步驟（3）：設結構體數組 P 和 T 的長度為 M，對 P 和 T 進行 M 次迭代。
找出相鄰波峰和波谷最大差值 Max 并記錄其索引，如式（4）所示：
| '/> |
其中，j 表示數組中的索引位置，val″ 表示幅值，idx″ 表示為在 x［n］中的索引。
步驟（4）：計算 Max［j］.idx″ 與 Max［j ? 1］.idx″ 的差值，即 BCG 節拍，判斷其是否在區間［Tmin，Tmax］范圍內。如果是，算法跳到步驟（6）。如果不是，視其為無效 BCG 節拍，算法繼續執行步驟（5）。
步驟（5）：當受試者翻動或側身時，通過傳感器采集到的肢體活動信號已經覆蓋大部分 BCG 信號，信噪比顯著降低，使得幾秒鐘內不可能進行可靠的心率估計。而在靜息狀態條件下，在短暫時間內心率變化不會太大，在大多數情況下連續的心臟跳動存在一定規律。本算法利用這些特征較準確地估計了在偽影下的 BCG 節拍，其目的是在數據不可用時將測量誤差盡量降低到最小。
如圖 3 所示，圖中紅點代表在 3 min 靜息條件下，檢測到的隨機選取的編號為 1 的受試者其連續標準 ECG 信號的心跳間隔。將其繪制成曲線并將其平滑處理，如圖中黑色波形線所示，可以觀察到，在短時間內（大約 10～20 次心跳間隔），曲線近似于一條微微傾斜的直線和一條含絕對值的正弦曲線的疊加，如式（5）所示：
圖3
在靜息下的心跳間隔
Figure3.
Heartbeat interval at rest
|
對于無效窗口，本文記錄最近的 20 次 BCG 節拍數據，并近似地擬合式（5）的曲線，跟據此曲線實現對無效 BCG 信號的 BCG 節拍估計。
步驟（6）：由步驟（4）或步驟（5）估計得到 BCG 節拍，由式（1）計算心率參數 HR。
2 實驗驗證
本文先使用傾角傳感器采集 BCG 原始信號，然后對信號濾波、整合，最后利用商業數學軟件 MATLAB 2016a（MathWorks inc，美國）分析處理過的 BCG 信號。
首先根據所提方法計算得出 J-J 間期（即 BCG 信號中相鄰兩個 J 峰之間的時間間隔），再由 Hamilton 等[9]的方法獲取標準的 R-R 間期（即 ECG 信號中相鄰的兩個 R 波之間的時間間隔），其次分別統計 J-J 間期和 R-R 間期的數量（分別以符號 J?J 和 R?R 表示），并計算覆蓋率（即 J-J 間期與 R-R 間期的數量百分比，以符號 CR 表示），統計假陽性（即 J 峰相對 R 波的延時時間超過 100 ms）和假陰性（即未檢測到的 J 峰數量），再計算假陽率（以符號 FP 表示）和假陰率（以符號 FN 表示），然后計算 J-J 間期相對 R-R 間期的平均逐拍間隔誤差（表示為：平均誤差 E）和絕對平均間隔誤差（表示為：絕對誤差 Eabs）。最后，通過這些計算量和誤差統計來評估算法的性能。
2.1 測試示例
如圖 4 所示，為隨機選取的編號 1 的受試者（以下簡稱為受試者 1）部分測試示例，圖中展示了藍框和紅框內的放大數據，綠色波形上明顯突出的部分波形代表受試者在有肢體活動時產生的偽影，紅點表明了算法匹配的心跳位置，綠色波形上重疊的藍色波形表明為所提方法對有偽影疊加的 BCG 數據模擬重建。由于 BCG 信號相對 ECG 信號有一定的時間延遲，為了實現數據同步和消除計算誤差，本文也設定了 0～200 ms 的 BCG 信號和 ECG 信號的延遲補償時間。
圖4
受試者 1 的部分測試數據
Figure4.
A partial data of the subject 1
如圖 5 所示，展示了受試者 1 的 J-J 間期和 R-R 間期，可以看出大部分情況下，藍色波形和紅色波形幾乎重合在一起，而藍色波形中部分突出的毛剌可認為是所估計偽影段的 J-J 間期；圖 5 還展示了受試者 1 逐搏間隔誤差的 Bland-Altman 圖[10]，從中可以看出 99% 的間隔誤差位于 ? 71～71 ms 之間，95% 的間隔誤差位于 ? 15～15 ms 之間；子圖平均心率圖展示的是受試者 1 平均心率的變化（即按分鐘統計心跳數量），所估計的平均心率與參考 ECG 信號的誤差范圍在 ? 2～2 beat/min 內。總的來說，從圖 5 可以觀察到 BCG 信號和 ECG 信號逐搏間隔之間具有良好的一致性。
圖5
受試者 1 的 J-J 間隔和 R-R 間隔、Bland-Altman 圖和平均心率
Figure5.
J-J interval and R-R interval, Bland-Altman plot and mean heart rate of Subject 1
2.2 數據分析
如表 1 所示，表中記錄了改進的峰值檢測方法的結果。假陰率和覆蓋率是相對的，隨著不同受試者疊加的偽影長度變化，錯過的假陰性隨之發生相應變化，覆蓋率也會發生改變。表 1 的性能數據代表不同受試者之間的平均誤差、絕對誤差、覆蓋率、假陰率和假陽率。算法對 BCG 信號疊加偽影的敏感度會影響覆蓋率。當受試者在床上有大量肢體活動時，其測試結果不是很理想，原因是：由于考慮到要保持與參考數據的同步，該方法只適用于處理含短暫偽影的 BCG 數據，對于較長偽影的 BCG 數據，不能正確地進行 BCG 節拍估計。通常情況下，該方法可以準確估計大部分 BCG 信號的心率，但如果由于 BCG 信號失真等原因，導致算法估計的當前心跳信號位置和實際位置相差較大時，估計的 J-J 間期發生偏移，從而產生累加誤差，影響下次 J-J 間期的估計，進而假陽率增加，其結果不可靠。但重要的是，該方法的誤差水平仍然維持在較低的范圍內，而且其絕對誤差也保持在可接受的限度內。
表1
改進的峰值檢測方法的性能
Table1.
The performance of improved peak detection method
此外，本文還使用相同的數據集評估了 Choi 等[5]所提出傳統峰值檢測方法的性能，如表 2 所示，和 Brüser 等[2]所提出層次聚類算法的性能，如表 3 所示。Choi 等[5]提出的方法可以在一段時間內找到 4 個子區間內的局部峰值，并從這些局部峰值中選擇最大峰值作為 J 峰，同時還附加了一些額外的規則；Brüser 等[2]使用層次聚類，先從記錄中推斷出心跳模板，然后根據該模板選擇出信號中最適的位置作為心跳位置，進而確定心跳間隔。從平均假陰率來看，改進的峰值檢測方法為 3.29%，而忽略掉偽影的 Choi 等[5]為 5.92%，表明偽影占比大概在 5.92% 左右，這說明了改進的峰值檢測方法通過重建部分偽影數據用于心跳間隔的估計。從平均假陽率來看，Brüser 等[2]幾乎為 0，但平均覆蓋率只有 89.39%，可以看出這種對波形要求較高的算法，提高了對心跳位置的識別率而降低了覆蓋率。改進的峰值檢測方法的覆蓋率、平均假陽率、平均誤差和絕對誤差分別為 96.71%、1.01%、0.86% 和 6.79 ms，這表明與 Choi 等[5]和 Brüser 等[2]相比，該方法獲得的 J-J 間期偏差更小，具有更高的覆蓋率和更好的整體準確性。總體來說，改進的峰值檢測方法要優于 Choi 等[5]和 Brüser 等[2]的方法。
表2
傳統 J 峰檢測算法的性能
Table2.
Performance of traditional J-peak detection algorithms
表3
層次聚類算法的性能
Table3.
Performance of hierarchical clustering algorithm
3 結束語
本文提出了一種改進 BCG 信號估計心率的方法，該方法是利用 BCG 信號的 J 峰提取心跳節拍，通過尋找局部 BCG 信號中相鄰波峰、波谷的最大差值的位置，將此位置設定為 J 峰，最后獲得 BCG 節拍。對于有偽影疊加的 BCG 數據段，本文設計了一種偽影模型可以近似估計該數據段的 BCG 節拍。
從 13 名受試者獲得 260 min 測試數據，得到約 20 000 次心跳的 BCG 信號樣本，本文對方法進行了驗證，可以看到實驗的數據輸出與黃金標準參考[9]之間具有良好一致性。實驗結果還表明，即使存在復雜的偽影和其他弱噪聲的情況下，該方法也能保證較低的誤差水平和較好的抗干擾能力。同時，該方法也與已有的其他方法進行了比較，實驗結果同樣具有較高的覆蓋率和較好的準確率。
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
