在运动状态下脉搏信号的级联滤波算法

2021-06-24孙再敏

激光技术 2021年4期

文武，孙再敏

(1.重庆邮电大学通信与信息工程学院，重庆 400065；2.重庆邮电大学通信新技术应用研究中心，重庆 400065；3.重庆信科设计有限公司，重庆 401121)

引言

脉搏信号是人体一种生理信号，包含着人体常见的心血管系统病理信息[1]。脉搏波蕴藏着人体大量的健康状态信息状态，可以为诊断疾病或者人体亚健康状态提供重要的科学依据。传统的检测仪器用粘性金属电极和导电凝胶与患者身体接触，给患者带来不便。早在20世纪30年代，HERTZMAN就首次提出了光电容积脉搏波信号描记法[2]。光电容积脉搏波描记术(photoplethysmography，PPG)是一种众所周知的用于监测生理信息的非侵入性方法,可以使用光电容积脉搏描记法来监测受测试者的健康状况。生理信息监测还可以作为一种运动监测手段,有助于人们合理安排运动强度和运动量[3]。

在剧烈的运动状态下，运动伪迹的存在使得脉搏波不能以理想的波形存在，为后续心率提取带来困难。为了减小运动伪迹，参考文献[4]中提出了一种新的方法使得PPG信号中的伪影(motion artifact,MA)去除结合了准确的心率(heart rate,HR)频率估计。参考文献[5]中采用了一种时间序列方法，该方法通过修正个体来修正运动伪影PPG脉冲，所提出的方法将每个脉冲替换为全局平均脉冲和特定脉冲的凸组合局部平均脉冲而不会丢失时间信息的脉冲。参考文献[6]和参考文献[7]中使用小波变换去除运动伪迹，这种方法虽然可以将脉搏信号进行分解，但是需要分析噪声在哪一层，而且对于需要确定小波基及设定阈值。参考文献[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)算法有效去除了频域上的运动干扰谱峰。

1 理论研究

1.1 奇异谱分析

奇异谱分析(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是这些特征值相应的特征向量。

(1)

式中，d=max{i}。轨迹X可以写成X=X1+…+Xd。则：

(2)

(3)分组。基本矩阵Xi被分割成m个不相交子集I1,…,Im。将对应于I子集的矩阵XI定义为XI=XI(1)+…XI(p)。轨迹矩阵X现在可以被写成：

X=XI(1)+…+XI(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

(4)

1.2 FTRLS

FTRLS算法[14]可视为两个横向滤波器，联合过程估计器和辅助滤波器并行工作。FTRLS算法大致可分解为前向预测、后向预测和联合估计3个部分,在收敛过程中实现3个部分之间的参量的互相交换、相互作用，从而实现FTRLS算法[15]。下文中，e(k,N)为先验误差，ε(k,N)为后验误差，下标f和b分别表示前向和后向。3个部分的大致过程如下。

(1)前向预测。在前向预测关系式中，瞬时前向后验误差为：

εf(k,N)=x(k)-wfT(k,N)x(k-1,N)=

(5)

式中,x(k)为PPG信号，wf(k,N)为滤波器在k时刻的前向预测权重向量。

前向预测权重向量[16]更新方程为:

wf(k,N)=wf(k-1,N)+

φ(k-1,N)εf(k,N)

(6)

式中,φ(k-1,N)为中间变量。

(2)后向预测。在后向预测关系式中，后向的后验和先验预测误差关系表达式为：

(7)

式中，γ(k-1,N)为k-1时刻后验误差与前验误差的转换因子。后向预测权重向量更新方程[17]为：

wb(k,N)=wb(k-1,N)+

φ(k-1,N)εb(k,N)

(8)

(3)联合估计。其先验误差可写成：

e(k,N)=d(k)-wT(k-1,N)x(k,N)

(9)

式中，d(k)为期望信号。联合过程估计滤波器权重向量为：

w(k,N)=w(k-1,N)+φ(k,N)ε(k,N)

(10)

2 滤波算法

正常脉搏波形和存在运动伪迹干扰脉搏波波形如图1所示。在自行车上高速骑行时，运动伪迹已经明显严重破坏了脉搏波波形，在强烈运动条件下脉搏波已不是正常波形，所以去除运动伪迹显得尤为重要。

Fig.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自适应滤波。

Fig.2 Filtering framework

Fig.3 Three-stage adaptive filtering

3 数据库

采用mimic数据库，采集了8名参与者(3名男性和5名女性)的腕部PPG数据，年龄为22岁～32岁(平均26.5岁)。参与者进行了一种或多种不同类型的锻炼，数据库是在加速度计和陀螺仪运动测量两种情况下进行体育锻炼时记录的信号。此数据库与以前的公共数据库不同，它包括在胸部同时采集的心电图(electrocardiography,ECG)，三轴陀螺仪信号、±2g和±16g三轴加速度信号、以及ECG采样时间和R峰时间，因为它包括在行走、跑步、轻松骑自行车和高阻力下骑自行车[21]。例如在高阻力下骑自行车时采集的12s到14s脉搏数据，如图4所示。

表1表示每个数据记录的持续时间。“0”表示数据库中参与者没有参与活动。

Table 1 Duration of each data record

Fig.4 The 12s～14s record information

4 实验

4.1 奇异谱分析

本文中采用互补滤波方法，利用三轴加速度计和三轴陀螺仪得到真实角度，利用陀螺仪角度来矫正加速度信号。三轴速度计需要滤掉高频信号，三轴陀螺仪需要滤掉低频信号，所以互补系数正好相加为1，调整它们所占的比重进行滤波。如图5所示,z轴角度矫正，加速度波形受到的干扰被滤掉，陀螺仪低频信号也被滤去。

Fig.5 z-axis angle correction

矫正后的加速度信号，利用SSA算法分组为不同频率成分的信号重构。SSA算法可以对窗口长度和分组两个参量进行优化。对于窗口长度，最好选择小于将要分解信号长度一半的值，可以平衡误差。在测试中，选择的窗口长度为350,重建效果最好。将加速度数据矫正后经过SSA分析，得到不同频率的信号以及主要成分。

SSA算法相当于压缩重构，去除冗余信息，减小信息大量存储问题，提高后续计算速度问题[22]。考虑不同特征值进行重建时，原始信号与重构信号互相关结果如表2所示。当信息被丢失过多时，信号将无法得

Table 2 Cross-correlation results of original signal and reconstructed signal

到正确的重构，结果表明，只考虑前8个特征三元组就可以在不损失太多信息情况下重建信号。

4.2 FTRLS

采用信号为1000采样点，采样率为256。将分组以后的加速度信号作为FTRLS的参考信号。结果如图6所示,自行车高速骑行下s1_high_resistance_bike的原始信号和FTRLS滤波去除运动伪迹以后的信号，以及对应的功率谱密度。FIR滤波器阶数为10，定义遗忘因子为1。去噪前和去噪后的功率谱对比，可以看到在FTRLS滤波以后频率窗口在500左右的第二大波峰去除掉，从而运动干扰带来的波峰去除，为后续心率计算带来便利，从滤波后的脉搏波的图中可以看到还保留了重博波信息。

Fig.6 Pulse wave to remove motion artifacts

4.3 滤波效果评价

为了验证提出的基于3级FTRLS滤波方法，将3轴加速度角度矫正后利用奇异谱方法分组重建，采用最小均方算法(least mean square,LMS)和RLS滤波，结果如图7所示。原始信号是s1_high_resistance_bike的数据，LMS算法、RLS算法都有一定的滤波效果，但是频率窗口在500左右时，只有FTRLS滤除了功率谱中第二大峰。