网站地图
-
脉搏信号是人体一种生理信号，包含着人体常见的心血管系统病理信息[1]。脉搏波蕴藏着人体大量的健康状态信息状态，可以为诊断疾病或者人体亚健康状态提供重要的科学依据。传统的检测仪器用粘性金属电极和导电凝胶与患者身体接触，给患者带来不便。早在20世纪30年代，HERTZMAN就首次提出了光电容积脉搏波信号描记法 。光电容积脉搏波描记术(photoplethysmography，PPG)是一种众所周知的用于监测生理信息的非侵入性方法, 可以使用光电容积脉搏描记法来监测受测试者的健康状况。生理信息监测还可以作为一种运动监测手段, 有助于人们合理安排运动强度和运动量[3]。
在剧烈的运动状态下，运动伪迹的存在使得脉搏波不能以理想的波形存在，为后续心率提取带来困难。为了减小运动伪迹，参考文献[4]中提出了一种新的方法使得PPG信号中的伪影(motion artifact, MA)去除结合了准确的心率(heart rate, HR)频率估计。参考文献[5]中采用了一种时间序列方法，该方法通过修正个体来修正运动伪影PPG脉冲，所提出的方法将每个脉冲替换为全局平均脉冲和特定脉冲的凸组合局部平均脉冲而不会丢失时间信息的脉冲。参考文献[6]和参考文献中使用小波变换去除运动伪迹，这种方法虽然可以将脉搏信号进行分解，但是需要分析噪声在哪一层，而且对于需要确定小波基及设定阈值。参考文献[8]中利用递归最小二乘(recursive least square, RLS)和归一化最小均方(normalizd least mean square, NLMS)自适应滤波器用于对PPG信号进行降噪。参考文献[9]中使用奇异谱分析将脉搏波分解，将分解后分量的主频与加速度信号的主频进行分析，去除加速度中噪声频率成分，重构后实现去噪。参考文献[10]中提出了一种基于多参考自适应噪声抵消技术的腕式PPG信号的鲁棒HR估计算法，其主要挑战是设计一个合格的参考噪声信号到自适应滤波器。参考文献[11]中提出3种方法去除运动伪迹, 包括3个部分：(1)PPG和加速度计信号的时频谱估计；(2)通过从PPG信号中减去加速度计信号的时频谱来去除运动伪影; (3)用三次样条插值法剔除影响心率的残余运动伪影。
自适应滤波对于参考信号有要求，需要选择合适的参考信号。基于单参考信号，本文中加入了三轴陀螺仪信号，并且运用奇异谱分析和快速横向递归最小二乘(fast transverse recursive least square, FTRLS)算法有效去除了频域上的运动干扰谱峰。
-
奇异谱分析(singular spectrum analysis，SSA)先将1-D时域信号首先构造成高维的轨迹矩阵，然后进行奇异值分解(singular value decomposition, SVD)，将奇异分量归类分组，最后重构矩阵[12]。计算轨迹矩阵和奇异值分解称为时间序列分解，分组和重构步骤称为重建[13]。奇异谱分析具体步骤如下。
(1) 计算信号的轨迹矩阵X。将1维信号YN=(y1, …, yN)进入多维的序列X1, …, XK。L是窗口长度, N为信号长度，2≤L≤N，K=N-L+1。轨迹矩阵X =[X1, …, XK]。
(2) 计算矩阵S = XXT。将矩阵S的奇异值分解，计算矩阵S的特征值和特征向量。λ1, …, λL是S的逐渐递减的特征值，U1, …, UL是这些特征值相应的特征向量。
$ \boldsymbol{V}_{i}=\boldsymbol{X}^{\mathrm{T}} \boldsymbol{U}_{i} / \sqrt{\lambda_{i}},(i=1, \cdots, d) $
(1) 式中，d=max{i}。轨迹X可以写成X = X1+…+ Xd。则：
$ \boldsymbol{X}_{i}=\sqrt{\lambda_{i}} \boldsymbol{U}_{i} \boldsymbol{V}_{i}^{\mathrm{T}} $
(2) 式中, ($\sqrt {{\mathit{\boldsymbol{\lambda }}_i}} $, Ui, Vi)是SVD的特征三元组。
(3) 分组。基本矩阵Xi被分割成m个不相交子集I1, …, Im。将对应于I子集的矩阵XI定义为XI = XI (1)+… XI (p)。轨迹矩阵X现在可以被写成：
$ \boldsymbol{X}=\boldsymbol{X}_{\boldsymbol{I}}(1)+\cdots+\boldsymbol{X}_{\boldsymbol{I}}(m) $
(3) (4) 对角线平均。每个矩阵XI被转换成一个新的长度N的时间序列。设X是元素xij的L×K矩阵，1≤i≤L，1≤j≤K，L'=min(L, K)，K′=max(L, K)，N=L+K-1。若L < K，则xij'=xji。然后，对角线平均将矩阵X转移到g0, …, gN－1，按照下列公式计算：
$ g_{k}=\left\{\begin{array}{l} \frac{1}{k} \sum\limits_{m=1}^{k+1} x_{m, k-m+2}{}^{\prime},\left(0 \leqslant k \leqslant L^{\prime}-1\right) \\ \frac{1}{L^{\prime}} \sum\limits_{m=1}^{L^{\prime}} x_{m, k-m+2}{}^{\prime},\left(L^{\prime}-1 \leqslant k \leqslant K^{\prime}\right) \\ \frac{1}{N-K} \sum\limits_{m=K-K^{\prime}+2}^{N-K^{\prime}+1} x_{m, k-m+2}{}^{\prime},\left(K^{\prime} \leqslant k \leqslant N\right) \end{array}\right. $
(4) -
FTRLS算法[14]可视为两个横向滤波器，联合过程估计器和辅助滤波器并行工作。FTRLS算法大致可分解为前向预测、后向预测和联合估计3个部分, 在收敛过程中实现3个部分之间的参量的互相交换、相互作用，从而实现FTRLS算法[15]。下文中，e(k, N)为先验误差，ε(k, N)为后验误差，下标f和b分别表示前向和后向。3个部分的大致过程如下。
(1) 前向预测。在前向预测关系式中，瞬时前向后验误差为：
$ \begin{aligned} \varepsilon_{\mathrm{f}}(k, N) &=\boldsymbol{x}(k)-\boldsymbol{w}_{\mathrm{f}}{}^{\mathrm{T}}(k, N) \boldsymbol{x}(k-1, N)=\\ &\boldsymbol{x}^{\mathrm{T}}(k, N+1)\left[\begin{array}{c} 1 \\ -\boldsymbol{w}_{\mathrm{f}}(k, N) \end{array}\right] \end{aligned} $
(5) 式中, x(k) 为PPG信号，wf(k, N)为滤波器在k时刻的前向预测权重向量。
前向预测权重向量[16]更新方程为:
$ \begin{gathered} \boldsymbol{w}_{\mathrm{f}}(k, N)=\boldsymbol{w}_{\mathrm{f}}(k-1, N)+ \\ \phi(k-1, N) \varepsilon_{\mathrm{f}}(k, N) \end{gathered} $
(6) 式中, ϕ(k-1, N)为中间变量。
(2) 后向预测。在后向预测关系式中，后向的后验和先验预测误差关系表达式为：
$ e_{\mathrm{b}}(k, N)=\frac{\varepsilon_{\mathrm{b}}(k, N)}{\gamma(k-1, N)} $
(7) 式中，γ(k-1, N)为k-1时刻后验误差与前验误差的转换因子。后向预测权重向量更新方程[17]为：
$ \begin{gathered} \boldsymbol{w}_{\mathrm{b}}(k, N)=\boldsymbol{w}_{\mathrm{b}}(k-1, N)+ \\ \phi(k-1, N) \varepsilon_{\mathrm{b}}(k, N) \end{gathered} $
(8) (3) 联合估计。其先验误差可写成：
$ e(k, N)=\boldsymbol{d}(k)-\boldsymbol{w}^{\mathrm{T}}(k-1, N) \boldsymbol{x}(k, N) $
(9) 式中，d(k)为期望信号。联合过程估计滤波器权重向量为：
$ \boldsymbol{w}(k, N)=\boldsymbol{w}(k-1, N)+\phi(k, N) \varepsilon(k, N) $
(10) -
正常脉搏波形和存在运动伪迹干扰脉搏波波形如图 1所示。在自行车上高速骑行时，运动伪迹已经明显严重破坏了脉搏波波形，在强烈运动条件下脉搏波已不是正常波形，所以去除运动伪迹显得尤为重要。
Figure 1. Motion artifact pulse wave and normal pulse wave
由于加速度计适用于长时间追踪角度变化，加速度计计算角度没有累积误差，三轴陀螺仪由角速度积分得到的角度由于温漂、单次迭代等原因, 往往偏差较大[18]。所以用互补滤波法使得它们的优势融合起来，进行角度矫正[19]。
FTRLS算法是通过快速横向滤波算法代替传统的RLS算法[20]。比起传统的RLS算法, FTRLS降低了计算复杂度，减少了算法收敛时所花费的时间。本文中采用的FTRLS算法框架如图 2所示。利用PPG信号中的运动伪迹分量具有同时加速度信号的频率，首先将校正后的三轴加速度计通过SSA，经过嵌入生成矩阵、SVD、归类分组和重构，分组时将相同频率成分的分量分为一组；其次，将奇异谱分析分组后的信号为参考信号，通过快速横向滤波RLS滤波器，将构造出的新分量依次输入到3级的FTRLS自适应滤波器中, 每个FTRLS自适应滤波器接收前级产生的剩余信号和本级对应的输入的参考信号分量。如图 3所示，PPG信号和x轴加速度参考信号sx输入到FTRLS自适应滤波器，输出信号作为2级输入和y轴加速度信号sy进行自适应滤波，同理输出作为3级输入与z轴加速度信号sz自适应滤波。
Figure 2. Filtering framework
Figure 3. Three-stage adaptive filtering
-
采用mimic数据库，采集了8名参与者(3名男性和5名女性)的腕部PPG数据，年龄为22岁~32岁(平均26.5岁)。参与者进行了一种或多种不同类型的锻炼，数据库是在加速度计和陀螺仪运动测量两种情况下进行体育锻炼时记录的信号。此数据库与以前的公共数据库不同，它包括在胸部同时采集的心电图(electrocardiography, ECG)，三轴陀螺仪信号、±2g和±16g三轴加速度信号、以及ECG采样时间和R峰时间，因为它包括在行走、跑步、轻松骑自行车和高阻力下骑自行车[21]。例如在高阻力下骑自行车时采集的12s到14s脉搏数据，如图 4所示。
Figure 4. The 12s~14s record information
表 1表示每个数据记录的持续时间。“0”表示数据库中参与者没有参与活动。
Table 1. Duration of each data record
namber walk run exercise activity low_bike high_bike 1 9min 48s 0 9min 39s 9min 48s 2 6min 39s 0 5min 41s 6min 54s 3 4min 47s 5min 7s 4min 54s 4min 41s 4 0 4min 52s 0 0 5 0 5min 8s 4min 40s 0 6 5min 36s 5min 2s 4min 40s 0 7 6min 42s 4min 47s 0 0 8 3min 40s 0 0 0 -
本文中采用互补滤波方法，利用三轴加速度计和三轴陀螺仪得到真实角度，利用陀螺仪角度来矫正加速度信号。三轴速度计需要滤掉高频信号，三轴陀螺仪需要滤掉低频信号，所以互补系数正好相加为1，调整它们所占的比重进行滤波。如图 5所示, z轴角度矫正，加速度波形受到的干扰被滤掉，陀螺仪低频信号也被滤去。
Figure 5. z-axis angle correction
矫正后的加速度信号，利用SSA算法分组为不同频率成分的信号重构。SSA算法可以对窗口长度和分组两个参量进行优化。对于窗口长度，最好选择小于将要分解信号长度一半的值，可以平衡误差。在测试中，选择的窗口长度为350, 重建效果最好。将加速度数据矫正后经过SSA分析，得到不同频率的信号以及主要成分。
SSA算法相当于压缩重构，去除冗余信息，减小信息大量存储问题，提高后续计算速度问题[22]。考虑不同特征值进行重建时，原始信号与重构信号互相关结果如表 2所示。当信息被丢失过多时，信号将无法得到正确的重构，结果表明，只考虑前8个特征三元组就可以在不损失太多信息情况下重建信号。
Table 2. Cross-correlation results of original signal and reconstructed signal
number 4 5 6 7 8 9 10 cross conelation/% 89.12 95.87 98.65 99.51 99.51 99.98 99.98 -
采用信号为1000采样点，采样率为256。将分组以后的加速度信号作为FTRLS的参考信号。结果如图 6所示, 自行车高速骑行下s1_high_resistance_bike的原始信号和FTRLS滤波去除运动伪迹以后的信号，以及对应的功率谱密度。FIR滤波器阶数为10，定义遗忘因子为1。去噪前和去噪后的功率谱对比，可以看到在FTRLS滤波以后频率窗口在500左右的第二大波峰去除掉，从而运动干扰带来的波峰去除，为后续心率计算带来便利，从滤波后的脉搏波的图中可以看到还保留了重博波信息。
Figure 6. Pulse wave to remove motion artifacts
-
为了验证提出的基于3级FTRLS滤波方法，将3轴加速度角度矫正后利用奇异谱方法分组重建，采用最小均方算法(least mean square, LMS)和RLS滤波，结果如图 7所示。原始信号是s1_high_resistance_bike的数据，LMS算法、RLS算法都有一定的滤波效果，但是频率窗口在500左右时，只有FTRLS滤除了功率谱中第二大峰。
Figure 7. Three methods of s1_high_resistance_bike to remove motion artifacts
本文中采用均方根误差(root mean square error，RMSE)ERMSE和信噪比(signal-to-noise ratio，SNR)RSNR对评价滤波效果评估。若RSNR值越大，ERMSE值越小，表明滤波效果越好[23]。均方根误差定义为：
$ E_{\mathrm{RMSE}}=\sqrt{\frac{1}{N}\left[\sum\limits_{k=1}^{N} s(k)-y(k)\right]^{2}} $
(11) 式中，s(k)为信号，y(k)为滤波后信号。
信噪比计算公式为：
$ R_{\mathrm{SNR}}=10 \lg \left(P_{\mathrm{s}} / P_{\mathrm{n}}\right) $
(12) 式中，Pn为噪声方差，Ps为信号方差。
从表 3中可知，本文中提出的滤波算法去除运动伪迹比LMS, RLS的信噪比分别提升了12.2%, 6.7%。同时，均方根误差分别减少了30%, 11%。相比LMS, RLS方法，本文中处理脉搏波的速度更快，完全满足脉搏波在实时处理的需求。
Table 3. Comparison of three filtering methods
method RSNR/dB ERMSE running time/s LMS 27.49 25.89 25.99 RLS 28.91 20.26 20.57 FTRLS 30.87 17.89 18.95 -
本文中提出的方法打破了自适应滤波单一的三轴加速度参考信号，采用了陀螺仪信号作为参考信号，结合了奇异谱分析和FTRLS算法，能有效去除噪声, 还能保留重博波信息。仿真实验验证了这种去除运动伪迹算法的有效性，在高速自行车骑行下脉搏波中去除运动伪迹取得了很好的效果。
在运动状态下脉搏信号的级联滤波算法
Cascade filtering algorithm of pulse signal in motion state
-
摘要: 为了减小在剧烈环境中脉搏波信号产生的运动伪迹，采用mimic数据库中由光子器件构成的光电传感器采集的脉搏信号，通过陀螺仪信号和三轴加速度互补滤波来矫正三轴加速度的角度，运用奇异谱分析将矫正后的三轴加速度信号分组为不同频率成分的信号，作为三级快速横向递归最小二乘(FTRLS)算法的参考信号自适应消除运动伪迹。结果表明，与最小均方算法、递归最小二乘法相比，信噪比分别提升了12.2%、6.7%，均方根误差分别减少了30%，11%，并且提高了计算速度。该研究对去除频域上脉搏波的运动干扰谱峰，保留重博波信息是有帮助的。
-
关键词:
- 光子器件 /
- 运动伪迹 /
- 互补滤波 /
- 奇异谱分析 /
- 快速横向递归最小二乘算法
Abstract: In order to effectively reduce the motion artifacts generated by the pulse wave signal in the severe environment, the pulse signal collected by the photoelectric sensor composed of photonic devices in the mimic database was used, and the gyro signal and the three-axis acceleration complementary filtering were used to correct the three-axis acceleration. Singular spectrum analysis was used to group the corrected three-axis acceleration signals into signals with different frequency components, which were used as the reference signal of the three-stage fast transverse recursive least square (FTRLS) algorithm to adaptively eliminate motion artifacts. Compared and analyzed with the least mean square algorithm and recursive least square method, the signal-to-noise ratio increased by 12.2% and 6.7%, the root mean square error was reduced by 30% and 11%, respectively, and the calculation speed was increased. The experimental results show that this research is helpful to remove the pulse wave motion interference spectrum peak in the frequency domain and retain the information of the heavy wave. -
Table 1. Duration of each data record
namber walk run exercise activity low_bike high_bike 1 9min 48s 0 9min 39s 9min 48s 2 6min 39s 0 5min 41s 6min 54s 3 4min 47s 5min 7s 4min 54s 4min 41s 4 0 4min 52s 0 0 5 0 5min 8s 4min 40s 0 6 5min 36s 5min 2s 4min 40s 0 7 6min 42s 4min 47s 0 0 8 3min 40s 0 0 0 
下载: 导出CSV
Table 2. Cross-correlation results of original signal and reconstructed signal
number 4 5 6 7 8 9 10 cross conelation/% 89.12 95.87 98.65 99.51 99.51 99.98 99.98
下载: 导出CSV
Table 3. Comparison of three filtering methods
method RSNR/dB ERMSE running time/s LMS 27.49 25.89 25.99 RLS 28.91 20.26 20.57 FTRLS 30.87 17.89 18.95
下载: 导出CSV
-
[1] LIN B. Research on chronic cardiovascular disease detection technology based on ECG singnal and PPG signal[D]. Xi'an: Xidian University, 2019: 1-40(in Chinese). [2] HERTZMAN A B. The blood supply of various skin areas as estimated by the photoelectric plethysmograph[J]. American Journal of Physio-logy, 1938, 124(2): 328-340. doi: 10.1152/ajplegacy.1938.124.2.328 [3] MA M Y. Research on noise cancellatlon in heart rate monitoring system[D]. Beijing: Beijing University of Posts and Telecommunications, 2018: 11-62(in Chinese). [4] ZHANG Q, XIE Q, WANG M, et al. Motion artifact removal for PPG signals based on accurate fundamental frequency estimation and notch filtering[C]// 2018 40th Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC). NewYork, USA: IEEE, 2018: 2965-2968. [5] MAJEED A M, SUJIT S, AHUL A, et al. Motion artifact removal of photoplethysmogram (PPG) signal[C]// Annual International Conference of the IEEE Engineering in Medicine and Biology Society 2019. New York, USA: IEEE, 2019: 5576-5580. [6] ZHANG Y H, GAO Ch Q, WANG R, et al. Methods of pulse wave analysis and its application[J]. Beijing Biomedical Engineering, 2019, 38(3): 319-326(in Chinese). [7] WANG Q, XU Y J, ZENG G Q, et al. Pulse wave extraction and dual domain analysis[J]. Application of Electronic Technology, 2019, 45(2): 53-57 (in Chinese). [8] ARUNKUMAR K R, BHASKAR M. CASINOR: Combination of adaptive filters using single noise reference signal for heart rate estimation from PPG signals[J]. Signal Image and Video Processing, 2020(3): 1-7. [9] ZHANG Z, PI Z, LIU B. TROIKA: A general framework for heart rate monitoring using wrist-type photoplethysmographic signals during intensive physical exercise[J]. IEEE Transactions on Biomedical Engineering, 2015, 62(2): 522-531. doi: 10.1109/TBME.2014.2359372 [10] CHOWHURY S S, HYDER R, HAFIZ M B, et al. Real-time robust heart rate estimation from wrist-type PPG signals using multiple reference adaptive noise cancellation[J]. IEEE Journal of Biomedical and Health Informatics, 2018, 22(2): 450-459. doi: 10.1109/JBHI.2016.2632201 [11] KONG Y S, CHON K. Heart rate estimation using PPG signal du-ring treadmill Exercise[C]// Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC) 2019. New York, USA: IEEE, 2019: 3253-3256. [12] HUANG H Ch. Heart rate monitoring from photo plethysmographic signals during motion state[D]. Nanjing: Southeast University, 2016: 10-50(in Chinese). [13] WEIDE L, JUAN Zh. An innovated integrated model using singular spectrum analysis and support vector regression optimized by intelligent algorithm for rainfall forecasting[J]. Autonomous Intelligence, 2019, 2(1): 46-55. doi: 10.32629/jai.v2i1.37 [14] HE H C. Research on adaptive time delay estimation method for pulsar navigation[D]. Changsha: University of Hunan, 2018: 10-66(in Chinese). [15] JUN Z, JING J Z, QIANG Ch. An improved FTRLS filtering algorithm and its simulation analysis[C]//Proceedings of 2017 17th IEEE International Conference on Communication Technology(ICCT 2017). New York, USA: IEEE, 2017: 1750-1755. [16] DINIZ P S R. Adaptive filtering: Algorithms and practical implement[M]. Berlin, Germany: Springer, 2014: 1-355. [17] CIODDI J, KAILATH T. Fast recursive-least-squares transversal filters for adaptive filtering[J]. IEEE Transactions on Acoustics, Speech and Signal Processing, 1984, 32(2): 304-337. doi: 10.1109/TASSP.1984.1164334 [18] LI Q D, ZHANG Y L. Improved attitude solving method based on Kalman filter and complementary filtering[J]. Journal of Qingdao University of Science and Technology (Natural Science Edition), 2019, 40(6): 85-89(in Chinese). [19] SHENG G R, GAO G W, ZHANG B Y. Research on quadrotor attitude calculation based on quaternion[J]. Modern Electronic Technology, 2020, 43(14): 8-12(in Chinese). [20] LENNART L, MARTIN M, DAVID F. Fast calculation of gain matrices for recursive estimation schemes[J]. International Journal of Control, 1978, 27(1): 1-19. doi: 10.1080/00207177808922343 [21] DELARAM J, ALEAANDER J C. Description of a database containing wrist PPG signals recorded during physical exercise with both accelerometer and gyroscope measures of motion[EB/OL]. (2020-09)[2020-09]. http://www.voiceofthebody.io/. [22] YANG Zh L, SHI W, CHEN H X. Adaptive compression sensing of optical fiber perimeter alarm signal[J]. Laser Technology, 2020, 44(1): 74-80(in Chinese). [23] JI Ch P, ZHAO L N. Research of wave filter method of human sphygmic signal[J]. Laser Technology, 2016, 40(1): 42-46(in Ch-inese). -
点击查看大图
图(7) / 表(3)
计量
- 文章访问数: 3919
- HTML全文浏览量: 2606
- PDF下载量: 17
- 被引次数: 0
