基于多区域分析的非接触式热红外视频心率检测方法*

陆磊，成娟

(合肥工业大学生物医学工程系， 合肥 230009)

1 引 言

心率是人体一项重要的生命体征参数，心率检测已经广泛应用于健康监护、情绪感知、疲劳驾驶和睡眠监测等领域。目前，心率检测方法分为接触式和非接触式两种。传统接触式心率检测方法的测量结果准确可靠，但是长时间佩戴易引起不适，且并不适合一些特殊场合的使用。例如，身体大面积烧伤、传染病人、新生儿监护等[1]。

基于图像光电容积描技术(imaging photo plethysmography，iPPG)的非接触式心率检测技术一般采用可见光摄像头，由于其低成本、使用广泛等特点逐渐成为研究热点[2-3]。但是，当环境昏暗或黑暗时，该技术检测效果不佳甚至不可用。由于热红外视频可以实时反映面部像素值由于心脏周期性搏动时血液循环引起的温度变化而产生相应的变化，可实现非接触式心率的检测，故可以作为基于视频的非接触式心率连续检测的补充方案，特别是昏暗或黑暗的情况。

热红外视频心率检测最早是由Garbey等[4]提出的，采用自适应滤波和快速傅里叶变换的方法来实现人体心率的测量。Gault等[5]利用热成像中温度的变化，检测在正常、轻度疼痛、轻度运动三种状态下的热视频心率，采用前额区域的60%的血管段结合小波分析的方法测试了人在正常、轻度疼痛、轻度运动状态下的心率，得到85%测量精度的结果。

Hamedani等[6]将欧拉放大方法应用到热红外视频心率提取，从而实现了心率的检测。Hu等[7]提出利用热红外成像仪和近红外摄像机来搭建双摄像头成像系统。对面部感兴趣区域(region of interest，ROI)进行像素均值时间序列提取，之后用时域信号分析方法检测心率。梁智敏等[8]对静态热红外视频检测，将热红外图片转化为灰度图片，通过对ROI灰度均值变化提取，然后经过小波分析和带通滤波处理后，最后利用小波包重构实现心率的实时检测。Hessler 等[9]使用小波变换过滤信号，计算传感器和热红外视频信号波形之间的相似性，实现了对热红外视频的心率和呼吸率的提取，达到了很好的相关性。以上研究中，基于热红外视频的心率提取研究也有与iPPG类似的技术。如欧拉放大的方法已经运用到热红外的心率检测中[6,10]。两者从根本上来说，均是心脏周期性收缩和舒张，引起面部像素值的相应变化。

近年来，基于iPPG的非接触式心率检测技术研究表明，通过同步分析多个不同的面部感兴趣区域的信号来提升心率检测的性能。Kumar等[11]通过剔除有噪声的面部区域，并使用剩余区域的信噪比作为权重，挑选出信噪比高的面部区域的绿色通道信号作为原始信号的来源，可以获得准确的心率估计值。Favilla等[12]选取了单色人脸视频中的前额，左右两侧脸颊这三个不同的感兴趣区域作为研究对象，得到相对应的三个单通道信号，将三个单通道信号经预处理后进行独立成分分析，得到相应的独立分量，最后用快速傅里叶变换提取频率范围内振幅最大的信号为心率信号，达到了心率检测的目的。Guo等[13]对多个面部感兴趣区域使用联合盲源分离的方法，应用独立向量分析联合分析多个面部感兴趣区域产生的多个颜色通道信号，得到测量的心率更加可靠更加准确。Wei等[14]提出通过对两个面部感兴趣区域中获得的六通道信号，应用二阶盲源分离技术来测量心率的方法，具有较强的鲁棒性。Qi等[15]采用联合盲源分离的方法对多个面部区域中颜色通道信号形成的多维数据集，筛选出接近真实心率的信号，从而实现了对心率的准确的检测。

对于热红外视频，单个ROI仅能生成单通道像素均值时间序列，在相对静止的环境下，多个面部ROI区域共同包含的信息主要为心率信息，故本研究拟借鉴iPPG技术中的多区域分析技术，形成多通道像素均值时间序列。一方面，提出多区域独立成分分析方法，同时，将单通道集合经验模态分解算法扩展为多通道的多变量经验模态分解，从而验证这两种多区域分析方法进行热红外视频心率检测的有效性。另一方面，与传统的单区域分析方法进行性能比对。

2 方法

2.1 多区域独立成分分析

盲源分离是指在没有关于信号混合过程先验信息的前提下，从一组观测到的混合信号中提取未观测到的信号即潜在的源信号的过程。通常情况下，观测信号由一系列传感器获取，每个传感器接受的信号由不同的源信号所构成。设x(t)=[x1(t),x2(t),...,xm(t)]是在t时刻的m维的源信号向量，其分量是相互独立的，设y(t)=[y1(t),y2(t),...yn(t)]T是n维观测信号向量，每个观测信号分量都是m个独立的源信号的线性组合，即：

y(t)=Ax(t)

(1)

这里A为n×m的混合系数矩阵。盲信号分离是研究从给定的混合信号中y(t)重构出源信号x(t)，独立成分分析(independent component analysis，ICA)是一种最常用的盲源分离技术。ICA的目标是找到一个对观测信号y(t)做线性变换的m×n矩阵W,使得y(t)经变换后的x′(t)各个分量之间尽可能独立，即:

x′(t)=Wy(t)

(2)

2.2 多变量经验模态分解

经验模态分解(empirical mode decomposition，EMD)是由Huang等提出的一种信号分析方法，分解的原理是它根据原始信号本身的时间尺度特征将其分解为若干本征模态函数(intrinsic mode function，IMF)。EMD方法尤其适用于处理非平稳和非线性数据方面具有明显优势。但是实际中，原始信号通常包含了大量复杂的噪声信息，这就可能出现模态混叠现象。为此，Wu等[16]提出一种噪声辅助分析方法，即集合经验模态分解(ensemble empirical mode decomposition，EEMD)，有效减少了模态混叠现象的发生。但是，EMD和EEMD都是对一维信号进行分析处理。

多变量经验模态分解法(multivariate empirical mode decomposition，MEMD)是EMD的一种扩展，实现多元信号的多通道同步联合分析[17-18]。利用MEMD来同时提取多个单通道之间的共同模式，可避免处理单个通道的信息而忽视了多个通道之间相互依赖的信息，将多个通道的信号分解为多个本征模式分量。对于输入信号x(t)=[x1(t),x2(t),...,xm(t)]，其算法具体步骤如下：

(1)在一个球面上，选择合适的m-1个采样点，可以得到m维空间的方向向量集nQk;

(6)通过d(t)=x(t)-m(t)，求得d(t)，如果d(t)满足IMF的停止条件，则为一个IMF，否则从步骤(2)开始执行。则MEMD最终会分解出M个IMF，即

(3)

式中，Cm(t)为第m个IMF，r(t)为残余量。

IMF停止条件为极大值点与极限值点个数之和与过零点的个数的差不超过1；极大值点与极小值点构成的包络的平均值应处处接近于0。

3 实验与结果

热红外视频心率检测实验流程，见图1。首先，利用热红外摄像机采集人体面部视频；其次，对视频图像帧确定N个ROI区域，分别计算每个ROI区域内的像素均值，构成时间序列；然后，分别采用MRICA和MEMD方法提取包含心率信息的独立分量和本征模态函数；最后，通过功率谱分析、筛选得到最佳的独立分量和本征模态函数，其主频即为心率频率。同时，为了验证多区域分析方法的优越性，还对比了单区域分析方法的性能，包括单通道滤波和集合验模态分解算法。

图1 热红外视频心率检测实验流程图Fig.1 The flow chart of heart rate detection with thermal infared video

3.1 视频数据采集

数据采集见图2。研究中使用高德C400热红外摄像机拍摄受试者的面部区域，尺寸为140 mm×206 mm×114 mm，波长为8～14 μm，图片分辨率为640×480，拍摄的热红外视频帧率为25 fps，温度灵敏度是0.1℃。使用康泰cms50d指夹式脉搏血压仪记录受试者的心率数据，脉率显示范围为25～250 bpm，误差为±1 bpm，帧率为60 fps。

图2 数据采集示意图Fig.2 Data acquisition schematic diagram

在房间内，每位受试者均坐在椅子上，尽量保持静止。固定好高德摄像头正对受试者，人脸距离摄像头约为0.5 m，每位受试者只拍摄一次，每次拍摄时间不低于1 min。脉搏血压仪夹在受试者手指上，提取心率数据，与高德摄像头同步采集热红外视频。有30位受试者自愿参与本研究，24名男性和6名女性，年龄在20～28岁。心率的检测时间为20 s，研究中读取前60 s的视频，分为不重叠的3段，共得到90个视频片段，其心率分布直方图见图3，心率的分布范围主要为60～90 次/min。

图3 心率分布直方图Fig.3 Histogram of heart rate distribution

3.2 面部感兴趣区域确定

选择鼻子作为感兴趣区域时会引入呼吸带来的运动干扰；选择嘴唇、眼睛作为感兴趣区域时会由于说话或眨眼造成的噪声而难以提取心率信号[19]。通常，研究人员会选择毛细血管分布丰富，且皮肤区域分布平整的脸颊和额头作为感兴趣区域，有助于高质量心率信号的获取。

首先，手动选取第一帧图像，确定3个ROIs，见图4。脸颊两侧、额头的ROIs是以鼻子所在的两条直线相互对称的，其中x=50，d=160，脸颊两侧的ROIs的像素大小为50×50，到对称中心距离为80；额头ROI到对称中心距离为240，像素大小为160×50。之后，计算每帧图像、每个面部感兴趣区域内所有像素点的均值，构成三通道像素均值时间序列。采用单区域分析方法进行对比时，单通道像素均值时间序列为这三个感兴趣区域内所有像素的平均。

图4 感兴趣区域的确定Fig.4 Region of interest

3.3 实验结果和分析

3.3.1基于多区域独立成分分析的视频心率检测 对多通道像素均值时间序列进行预处理，包括归一化和滤波[20]。本研究使用4阶巴特沃斯滤波，频率范围为0.75～2.5 Hz，即对应于心率值45～150次/min。

MRICA是通过假设三通道信号为包含心率成分的脉冲信号，以及运动和其他信号等多种成分的线性组合。采用独立成分分析处理后可以得到三个独立分量，选择幅值占比最大的信号作为心率信号。最后由最佳信号的功率谱图最大幅值对应的频率可以估算出检测的心率。选择第7个受试者20～40 s数据做分析。三个独立分量时域图及频谱图见图5。最后通过三个独立分量的频谱图选取能量占比最大的信号为最佳候选信号。由计算得能量占比P1=0.045，P2=0.077，P3=0.095，第三个信号能量占比最大即为最佳候选信号，其主频为f=1.2 Hz，得到HR=f×t=1.2×60 =72 bpm，参考信号参考值也为72 bpm。

图5 独立分量Fig.5 Independent component

3.3.2基于多变量经验模态分解的视频心率检测 将归一化的多通道像素均值时间序列，经MEMD分解后，分解为多个本征模式分量，同时提取3个信号通道之间的共同模式。IMF的个数和信号的长度有关，由于研究中的20 s数据有500帧图像，可以得到9个IMF，3个通道共有27个IMF，我们将其记为IMF1，IMF2，...，IMF27。当IMF主频落在感兴趣的心率频率范围0.75～2.5 Hz之间，即对应于心率值45～150次/min，则我们便将其选定为候选的IMF。每个通道均至少有一个候选IMF，在这些候选IMF中选择主频能量占比最大的IMF作为最佳本征模式分量，得到对应的频率和心率检测值。选择受试者4的第41～60 s数据进行分析，得到IMF4、IMF5、IMF1、IMF22为候选IMF，图6是确定的候选IMF时域图和频谱图，计算4个候选IMF的能量占比，分别为P1=0.090，P2=0.072，P3=0.063，P3=0.087。我们发现IMF4的能量占比最大，即为最佳候选IMF，其主频为1.51 Hz，最终得到HR=f×t=1.51×60=90.6 bpm，而相应时间段内的参考心率为91.6 bpm。

图6 候选IMFFig.6 IMF candidates

3.3.3基于多区域分析的视频心率检测结果 为了验证基于多区域分析方法的有效性，将基于单区域分析的单通道滤波和EEMD算法进行对比。四种方法的结果见表1，以4个评价指标来展现，分别为平均绝对误差MAE、标准差SD、均方根误差RMSE、皮尔逊相关系数r[21]。

表1 不同方法的心率检测性能对比Table 1 Comparison of heart rate detection performance of different method

由表1可知，MRICA方法的心率检测结果最好，平均绝对误差为3.17 bpm，标准差为2.93 bpm，均方根误差为4.3 bpm，相关系数为0.87。对比MRICA和MEMD，以及EEMD和单通道滤波方法的结果，可知多区域分析方法得到心率检测的结果优于单区域分析方法。其中单通道滤波效果最差，主要原因为信号本身微弱，单纯做滤波处理，心率信号易淹没在噪声中。EEMD的结果要优于单通道滤波，EEMD方法有一定的去噪效果，但仍有一些包含心率信号的IMF淹没在噪声之中。MEMD是将三个通道数据分解为多个IMF，同时提取三个信号通道之间的共同模式，减少了心率信号的淹没，得到更多包含心率信息的IMF。MRICA是通过假设三通道信号为包含心率成分的脉冲信号等多种成分的线性组合，可以达到较好的心率检测结果。

为了更好的展现算法性能，本研究采用Bland-Altman和散点图来分别进行一致性和相关性分析。图7为不同方法在数据库上的心率检测结果Bland-Altman图。由于单通道滤波算法性能较差，后续仅与采用EEMD算法的单区域分析方法进行对比。

图7 不同方法的Bland-Altman图Fig.7 The Bland-Altman plots of different methods

由图7可知，MRICA方法的心率差值均值为0.2 bpm，一致性界限为-8.3～8.7 bpm。MEMD方法的心率差值均值为-1.6 bpm，一致性界限为-13.2～10 bpm。EEMD方法的心率差值均值为-2.8 bpm，一致性界限为-35.9～30.3 bpm。相较于EEMD，多区域分析方法的一致性界限范围明显缩小，落在一致性界限范围外的点也大幅减少，具有更好的一致性。

图8为不同方法在数据库上心率检测结果的散点图。其中横轴表示从参考信号中获取的心率真值，纵轴表示不同算法获得的心率测量值。可以看出，EEMD算法的散点图中，点的分布较为分散，这说明心率估计误差相对较大。而MRICA和MEMD方法的散点图中点的分布更加集中于基准线y=x附近，这说明多区域分析方法所获得的的心率测量值与心率真值之间具有很强的相关性，性能更好。

图8 不同方法的散点图Fig.8 The scatter plots of different methods

4 总结

本研究利用相对静止环境下，面部多个感兴趣区域共同包含的主要为心率信息这一特点，分别采用多区域独立成分分析和多变量经验模态分解算法进行心率检测，实验结果表明，多区域分析方法相比于单区域分析方法，更能表征潜在的心率信息成分，且多区域独立成分分析方法性能最佳，平均绝对误差3.17 bpm，均方根误差2.93 bpm，标准差4.30 bpm，相关性达到0.87，为基于视频的非接触式心率的24 h连续检测提供解决方案。同时，本研究尚存在不足，目前只是探究了相对静止环境下热红外视频的心率检测研究，未进行实时检测心率，如何减少实际环境中的运动干扰及实时检测心率是亟需探索的方向。而一些基于神经网络的方法已经可以提取普通摄像头的视频心率。未来，一些基于神经网络的方法在热红外视频心率的提取也是探索的方向[22-23]。