头部旋转运动下自适应非接触鲁棒性心率检测方法*

巴图巴雅尔·欧赟 赵跃进3)† 孔令琴3)‡ 董立泉3) 刘明3) 惠梅

1) (北京理工大学光电学院,北京 100081)

2) (北京理工大学,精密光电测试仪器及技术北京市重点实验室,北京 100081)

3) (北京理工大学长三角研究院(嘉兴),嘉兴 314019)

基于人脸视频的生理信号检测面临的主要挑战是运动伪影噪声.针对受试者头部刚性旋转运动引起的伪影噪声,本文提出利用头部运动信息构建自适应滤波器的非接触式心率检测方法.该方法利用人脸二维和三维的特征点计算受试者运动中头部的偏航和俯仰欧拉角度,并将其作为调控过程噪声协方差的信号质量指数,进而构建了自适应Kalman 滤波器,实现了稳健的心率估计.实验结果表明:本文提出的方法可有效抑制头部刚性旋转运动引起的噪声,平均绝对误差为2.22 beat/min,均方根误差为2.76 beat/min,与现有方法相比准确度分别提升9%与24.6%,具有统计显著性.本文提出的头部旋转角度自适应非接触鲁棒性心率检测方法在自发运动的真实场景下能有效提升检测的准确性,扩大了成像式光电容积描记技术在视频健康监测领域的使用场景.

1 引言

心脑血管疾病已经成为人类健康的主要威胁,根据心脑血管疾病早发现、早治疗,积极干预的防治需求[1],寻求一种操作简单、结果准确、可重复性好、并可应用于大规模人群日常检测的测量技术具有重要意义.成像式光电容积描记(imaging photoplethysmography,IPPG)技术是在光电容积描记(photoplethysmography,PPG)技术基础上发展的一种基于成像设备的非接触生理参数检测技术[2,3],克服了传统PPG 技术存在的缺陷,正逐渐应用于远程医疗监控等领域,成为近年来研究的热点[4,5].

基于IPPG 的心率测量技术研究已经取得了重要进展,可以在受试者静止及轻微非刚性运动状态下实现心率测量[6,7],运动伪影是目前限制该技术的主要瓶颈.解决运动伪影的常用方法主要包括基于数据的颜色通道线性组合法、基于先验知识的颜色通道组合法及基于图像预处理的方法.基于数据的颜色通道线性组合法是依据血红蛋白的光吸收在整个光谱范围内的相对强度差异性,利用数据统计线性组合人体皮肤在三色通道的信号实现运动鲁棒性的方法,通常绿色通道具有较强的脉搏信号强度[8].2008 年,Hülsbusch[9]构建了双色通道线性组合法,首次利用该强度的差异性将PPG 信号和运动伪影分成两个独立的信号;2010 年,Poh 等[10,11]扩展了这项工作,提出了通过独立分量分析(independent component analysis,ICA)方法,并将非高斯性作为独立标准线性组合三通道信号;Lewandowska 等[12]在2011 年提出主成分分析(principal component analysis,PCA)法重新定义了三通道的独立线性组合,使用盲源分离技术(blind source separation,BSS)分离了PPG 信号及噪声.基于数据的颜色通道线性组合法认为携带脉搏信号的成分是先验未知的,同时显示出最强的周期性,但由于运动伪影同样具有强周期性使其不适用于真实运动场景.由于基于数据的线性组合法需要较长观测时间,有研究提出了基于先验知识的颜色通道组合法,根据不同的经验推理加权组合三色通道信号克服运动伪影.其中具有代表性的研究如下:2013 年Haan 和Jeanne[13]提出的基于皮肤反射物理模型组合三通道信号获得色度(chrominance,CHROM)信号的方法,首次实现利用IPPG技术测量单人受试者在健身器材上运动时的心率;2017 年,Wang 等[14]提出了颜色失真滤波算法(color-distortion filtering,CDF),选取不同的颜色特征变化作为噪声和信号的分离依据,从而实现对RGB 信号滤波;2017 年,Wang 等[15,16]提出的基于正交肤色平面(plane orthogonal to skin,POS)的方法,在RGB 空间中寻找一个与皮肤正交的平面以提取脉冲信号.上述研究均是针对信号处理部分,实际上对图像进行跟踪配准等预处理也可以实现在一定程度减小运动伪差影响,如图像配准法[17]、感兴趣区域(region of interest,ROI)空间平均法[6]、ROI 跟踪法[18]、依据时频分析法的抗运动频谱峰值跟踪[19],及适用于长距离移动的基于自适应变焦系统心率检测法[20]等.

然而,现有方法均只适用于呼吸、头部轻微移动等非刚性运动场景,当受试者进行旋转头部和行走等自发性运动时无法准确估计心率.针对上述问题,本文提出了一种利用角度信息构造自适应滤波器抑制头部刚性旋转运动伪差的方法,提供单一稳健的心率估计.该方法在目前基于IPPG 心率测量研究的基础上扩大了在自发运动状况下的使用场景,满足了在头部刚性旋转运动场景下非接触的实时准确测量心率,从而有效监控运动中受试者的健康状态,进一步推动了IPPG 在生理信号检测和视频健康监测领域的发展.

基于人脸视频[21]的IPPG 心率测量技术,利用成像设备对包含被测部位的信息进行视频采集[22],将脉搏信号即由血液容积变化引起的光强变化用视频图像的方式记录下来,再通过对视频图像处理提取出脉搏波信号[23-25],最后通过脉搏波特征量的分析实现心率信号的提取[4,26,27].事实上,在普通的生活场景中,如在健身、长期远程连接医疗监控技术时,难以避免受试者相对于相机的运动.下面将针对受试部位(脸部)刚性旋转运动引起的运动伪差,分析其产生原因,并对本文提出的噪声抑制方法进行阐述.

2 方 法

基于人脸视频的IPPG 心率测量技术,利用成像设备对包含被测部位的信息进行视频采集,将脉搏信号即由血液容积变化引起的光强变化用视频图像的方式记录下来,再通过对视频图像处理提取出脉搏波信号,最后通过脉搏波特征量的分析实现心率信号的提取.

2.1 刚性运动伪差

在运动的场景中,将人脸视作一个刚体,即其在三维空间中的旋转运动可以分为3 种类型:俯仰(pitch)、偏航(yaw)、翻转(roll),如图1 所示.

图1 人脸的3 个旋转运动自由度Fig.1.Three rotational freedom degrees of human head.

当生物皮肤组织在各个方向具有均匀照度的理想环境(积分球)中,且观察角度小于60°时,可近似视作扩展的朗伯辐射源[28].假设人脸在理想光照环境中运动,即将人脸视作扩展朗伯辐射源,其辐亮度与方向无关.

当人脸保持静止时,成像系统像平面的辐照度E可表示为

其中,τ表示成像系统透过率;L表示ROI 辐亮度;D/f′表示光学系统数值孔径.

CCD 相机接收到的总功率P可表示为

其中,AEP表示入射光瞳面积.

原始IPPG 信号正比于CCD 接受总功率与ROI 的像面积之比

人脸作为扩展朗伯辐射源,所以当人脸以角度α旋转时辐亮度L不随角度α变化,只有ROI 的像面积随角度α变化,因此(3)式应改写为

由于当人脸做翻转运动时CCD 相机获取的像面大小实际上没有改变,并且在真实场景中头部的上下、左右方向的旋转运动更频繁,故在本文中仅考虑俯仰和偏航运动对心率信号的影响.

人脸做俯仰和偏航运动会引起 cos2α的变化,从而对IPPG 信号引入误差.然而真实场景中的环境与均匀照度的积分球的差异较大,人脸ROI 作为接收面,入射光强会随着接收面法线和辐射方向的夹角改变而变化,导致CCD 相机采集视频时人脸ROI 作为辐射面的辐亮度L也会随着之改变.

综上所述,当人脸做俯仰或偏航运动时,原始IPPG 信号会受旋转角度α变化及ROI 入射光强变化的调制,继而引入了运动伪影.

2.2 旋转运动自适应滤波器

依据贝叶斯算法,在测量心率前对真值的预测为该时刻的先验概率,在得到测量值后对真实值的估计为后验概率.由于后验概率中包含测量过程中产生的噪声,例如运动伪影等,故利用先验概率修正信号实现对错误心率估计的滤波.本文仅考虑测量过程中旋转头部运动的影响,如无特别说明下文中提到的测量噪声均指旋转头部运动引入的伪影.

依据2.1 节的分析,本文利用贝叶斯框架,设计了一种头部旋转运动自适应滤波器,其实质是根据人脸旋转的角度信息自适应地改变过程噪声协方差估计,在新数据到达时准确估计心率动态变化.

2.2.1 自适应过程噪声协方差

为从人脸视频中获取头部旋转角度信息,即头部姿态欧拉角,首先利用 CLNF 算法定位面部特征点坐标,并将二维点反向投影到三维人脸模型,得到三维点坐标,再通过标定获得相机标定参数,最后利用最小二乘法求解透视位姿得出头部姿态欧拉角[29-32].本文仅考虑俯仰角度(pitch)和偏航角度(yaw):

对所获角度的绝对值进行归一化并取反操作,使数值范围在0—1 之间.并满足旋转角度较大时数值趋近于0,角度较小时数值趋近于1:

令信号质量指数为

利用信号质量指数θSQI调节测量噪声协方差估计,构成头部旋转自适应协方差RSQI:

其中,R0取经验值0.1[33].由此将获得的自适应测量噪声协方差与Kalman 滤波器结合,使得滤波器在头部旋转欧拉角较大即数据置信度较低时,测量噪声协方差较大.

2.2.2 基于RSQI调控Kalman 增益构造自适应滤波器

本文提出的贝叶斯决策是基于张玉燕等[34]和Nemati 等[35]提出的Kalman 滤波器基础上实现的.Kalman 滤波器提供了一种有效的递归解决方案,它通过组合之前状态值和当前观测值,推断或重构最符合系统状态真实值的最优估计.

假设第k个的心率真实状态值为xk,观测值为yk,分别满足如下的离散时间差分方程:

其中,A为状态转移矩阵;xk-1为前一时刻的心率状态值;uk为该时刻系统的控制输入变量;H为系统的测量参数;wk-1为过程激励噪声,vk为测量噪声,二者均为高斯噪声且相互独立,他们的正态分布均值为0,协方差分别为.由于人体在稳定状态下心率不会突变,因此假定A为单位矩阵;忽略控制输入变量、将测量参数H设为1;忽略过程激励噪声的影响,将Q设置为经验值1.(10)式和(11)式改写为

假设第k-1个心率为均值为、协方差为的正态分布,即第k-1 个时刻的后验估计的协方差为

预测步骤

根据贝叶斯规则,第k个时刻的先验估计及协方差为

即预测步骤为:当前时刻的先验概率由上一时刻的后验概率决定.

更新步骤

可以由预测步骤得到第k个时刻的状态更新,即后验估计及协方差为

(16)式和(17)式中,Kalman 增益是Kalman 滤波过程中的核心步骤.利用(9)式得到的自适应测量噪声协方差来调节Kalman 增益系数,定义如下:

当出现伪影时,估计测量噪声协方差RSQI会增大,显著降低Kalman 增益,此时降低测量数据的置信度,更信任该时刻的先验概率分布.

图2 为本文提出自适应滤波器的算法流程图.首先,使用头部姿势欧拉角信息构建心率信号质量指数,用于估计自适应测量噪声协方差RSQI.然后,利用RSQI计算Kalman 增益系数;根据头部姿势角度值的动态变化,结合Kalman 增益系数、心率测量值和当前时刻的先验概率,实现对新到达数据的自适应滤波.最终,将当前状态的后验概率更新为后续状态的先验概率.当某一时刻后不再有新的心率测量值传入时,退出该更新迭代循环.

图2 头部旋转角度信息作为自适应滤波器质量控制的流程图Fig.2.Flow chart of head rotation angles as a quality control of our proposed adaptive filter.

3 实验系统

3.1 实验方案

实验方案如图3 所示.实验中采集的视频帧率为30 frame/s,分辨率为1024×768,并保存为AVI格式.所有测量均采用环境光稳定的房间,避免了环境光强突变对实验结果的干扰.在所有实验中,CCD 相机采集人脸视频,然后通过数据线传输到电脑端进行存储和分析.同时,利用指夹式血氧仪(德州仪器公司,AFE4400SPO2EVM)进行心率检测,并将其测量结果作为参考真值.

图3 实验方案图 (a) 实验装置示意图;(b) 采集到的视频片段((i)头部姿势的欧拉角为0°;(ii) 头部俯仰运动;(iii) 头部偏航)Fig.3.Experimental plan diagram:(a) Experiment set-up;(b) video clips ((i) the Euler angle of the head pose is 0°;(ii) the head pitch segment;(iii) the head yaw segment).

参加实验的志愿者为21 名年龄在23—28 周岁之间的在校研究生,志愿者在参加实验前被告知了实验目的且均是自愿参加实验,所有参加实验的志愿者均身体健康无心脏病等疾病.实验中,要求受试者距离CCD 相机0.5 m 保持静止,并正向面对相机(视线与相机法线平行)随机做出头部俯仰及偏航的动作.

3.2 感兴趣区域提取

本文设计了如下由面部关键点构建的动态脸颊与鼻翼ROI 区域提取方法:选取面部68 个关键点中第1,17,42 和47 个关键点构建人脸ROI 区域,如图4(a)所示.由p1和p2关键点决定的矩形区域作为ROI 区域,如图4(b)所示.

图4 中,p1和p2关键点的坐标由下式确定:

式中pxn和pyn(n=1,16,41,46) 代表 点n的x坐 标和y坐标;加(减)50 代表向右(左)偏移50 个像素值,如图4(b)中蓝色圆点所示;h代表当前帧的人脸检测高度.采用两侧关键点分别向内偏移50 个像素值决定矩形框左右边界,保证ROI 区域内不会随着头部的旋转运动(如图4(c))而采集到环境背景,从而很好地避免了引入唇色噪声和背景噪声.

图4 矩形ROI 区域选取示意图 (a)面部68 个关键点检测;(b)脸颊矩形ROI 区域;(c)偏航运动时ROI 区域Fig.4.Flow chart of rectangular ROI area selection:(a) 68 feature points detection;(b) rectangular ROI area;(c) ROI area during yaw motion.

3.3 数据处理

本文提出的心率检测方法算法流程如图5 所示.首先采用开源工具OpenFace2.0 对采集到的视频进行人脸检测,基于特定人脸关键点选取毛细血管分布丰富的脸颊和鼻翼区域作为ROI,对每一帧的人脸图像提取ROI 区域可以有效克服因受试者头部移动等运动引入的噪声伪影.然后将每一帧ROI 的像素均值按视频时间序列连接生成原始IPPG 信号.然后将每一帧ROI 的像素均值按视频时间序列连接生成原始IPPG 信号.再采用先验平滑滤波(prior smoothing filter,PSF)[36]对IPPG信号进行去趋势处理,去除信号的低频漂移趋势噪声.通过小波滤波算法消除高频噪声,具体为首先选择db12 小波基对原始信号进行6 层分解,由于人体脉搏信号频率一般在0.7—4.0 Hz 之间,因此选择第4,5,6 层的细节信号重构脉搏波信号(对应的频率范围为0.469—3.75 Hz).最后对去噪后的信号建立10 s 的滑动窗口进行快速傅里叶变换[37],设置滑动步长为0.2 s,并对窗口内的有效数据进行补零操作避免频谱泄漏,获得其对应功率谱,并利用功率谱峰峰值所在频率估算心率.并对指夹式血氧仪获取的一维时域脉搏信号同样以10 s 的窗口大小及0.2 s 的滑动步长进行快速傅里叶变换,通过其功率谱计算得到参考真实心率值.

图5 本文提出方法的流程图Fig.5.Flow chart of method proposed in this paper.

在上述基础上,结合人脸二维及三维特征点计算头部姿势欧拉角度,进而根据角度绝对值获得信号质量指数θSQI,利用θSQI估计测量噪声协方差RSQI从而调控Kalman 增益,以此构建一个头部旋转运动的自适应滤波器.该滤波器能够依据头部旋转角度的变化动态滤除旋转运动引入的噪声,进而实现单一稳健的心率检测.

4 结果与讨论

为验证本文提出算法的可行性,将所提出的滤波器与目前性能较佳的CHROM 方法[13]结合并进行对比实验.

图6 为某一受试者的测试结果.图中浅灰色和深灰色虚线分别为头部偏航和俯仰角度变化曲线,蓝线为参考设备指夹式血氧仪测得的心率值.绿线为仅采用CHROM 算法测得的心率值,红线为加入本文提出的自适应滤波器算法后测得的心率值.从图6 可以看出,当局部出现运动伪影时,CHROM不能有效抑制运动伪差,而结合本文提出的算法后,能稳健地抑制运动伪差,实现心率准确测量.

图6 本文提出的滤波器与CHROM 结合的对比结果Fig.6.Comparison results of the filter proposed in this paper combined with CHROM.

为定量描述上述两种方法的21 组实验结果与真实心率值之间的一致性,用相关性图进行了结果分析,结果如图7 所示.

图7 给出了CHROM 方法和本文提出的滤波器与CHROM 结合的21 组实验结果相关性图.图中任意一点的横轴表示实验的心率估计值,纵轴表示指夹式血氧仪获得的心率参考真值.从图7 的相关性结果中可以看出,与传统的CHROM 方法相比,利用本文所提方法获得的心率估计值与心率真实值更趋于一致,表明头部旋转自适应滤波器显著提高了心率估计结果与真实心率结果的一致性.

为了进一步证实这一结论说法,将本文设计的滤波器与不同的传统方法(GREEN,CDF,POS,CHROM)结合前后的性能(误差)进行定量分析,结果如表1 所列.

表1 为本文设计的滤波器与不同的传统方法(GREEN,CDF,POS,CHROM)结合前后测得心率值的MAE 和RMSE 误差结果,从表1 可以看出,4 种方法增加了本文的自适应滤波器后MAE及RMSE 均有明显提升效果;其中在CHROM 的基础上组合本文的自适应滤波器后的效果最优,其MEA 及RMSE 值分别为2.22 beat/min 和2.74 beat/min,较CHROM 的结果提升了9%和25%,R2相关系数提升了3%.尤其应注意,CHROM叠加CDF 再与本文的自适应滤波器结合后,R2相关系数最优,高达0.8901,说明CDF 方法与本文的滤波器结合能有效提升心率信号与真实值的相关性.上述结果分析表明,本文提出的头部刚性旋转自适应滤波器能有效抑制旋转运动引入的伪影并且显著提升估计心率值的准确性.

图721 组实验结果相关性图(图中红线表示y=x 的线性关系) (a) CHROM 方法的相关性图;(b)本文提出的滤波器与CHROM 结合的相关性图Fig.7.Correlation plots of 21 groups of experimental results (the red lines in the plots indicate linear relationship of y=x):(a) Correlation plot of CHROM method;(b) correlation plot of combining CHROM with our proposed adaptive filter.

表1 本文提出的自适应滤波器与不同的传统方法结合前后实验结果对比Table 1.Comparison of experimental results before and after combining the filter proposed in this paper with different traditional methods.

需要说明的是,本文在处理头部姿态欧拉角度时应用的归一化操作使用了视频片段中的最大欧拉角度用于确定信号质量指数,在实时监控受试者的心率值时,会有一定的处理时间延迟,因此我们会在未来的工作中继续进行相关的实验研究.Kalman 滤波器参数的初始化值是影响该技术性能的重要因素,本文并没有确定心率检测的最优滤波器参数,相关工作将在日后的研究中进行探讨.此外,本文提出的由头部姿态信息估计自适应测量噪声协方差,并以此构建自适应Kalman 滤波器抑制运动伪差的方法,适用于任何由运动伪影引入噪声的时间序列信号,可以扩展到其他非接触式生理参数检测及跟踪技术领域.

5 结论

本文针对非接触式心率检测时受试者刚性运动引入伪影的问题,提出了一种根据头部姿态角度作为主要失真指标提高心率测量准确率的方法.该方法利用人脸检测技术获得受试者头部姿态欧拉角度并以此计算信号质量指数θSQI,根据θSQI估计测量噪声协方差来调控Kalman 增益,从而构建一个头部刚性旋转运动的自适应滤波器,利用该滤波器有效提高了心率检测准确性.与传统的适应于运动场景的心率检测算法的对比分析表明,本文提出并设计的自适应滤波器有效解决了头部刚性旋转运动引起伪差的问题.该方法使基于人脸视频的心率检测技术在现实场景中有效监控受试者自发运动时的健康状态成为可能,进一步推动了IPPG 技术在生理信号检测和视频健康监测领域的发展.