一种非接触式心率变异性测量方法与流程

本发明涉及一种非接触式心率变异性测量方法，应用在人体心理压力的评价技术领域。

背景技术

目前心率变异性的测量主要基于心电信号，即通过心电计算心动周期，并根据逐拍心动周期计算心率变异性。由于心电测量需要与人体接触，其应用受到很多局限，例如被试者可能会抗拒心电测量。

目前非接触式的心率测量，通过摄像机拍摄人面部视频获得；如申请号：201310172275.9公开的非接触式心率测量系统及测量方法，提出当拍摄视频后，通过对面部图像的多个颜色分量进行盲源分离，筛选出与心率相关的分量，然后通过频谱分析得到心率。但是在对人体的疲劳程度、心理状态、压力大小等进行评价和分析时，仅有心率是不够的，还需要进行心率变异性的测量。

上述方法在测量心率变异性时存在缺陷，主要体现在：1、心率分析需要依赖频谱分析，但是频谱分析只能反映较长时间内的平均心率，而不能反映这一段时间内心率的变化，即心率变异性。2、使用盲源分离后筛选出的与心率相关的信号分量来计算心率时，存在较大的噪声和干扰；由于计算心率变异性比心率更易受噪声的干扰，盲源分离后分量中的噪声虽然对计算心率带来的误差很小，但是会对计算心率变异性引入巨大的误差，甚至错误。3、普通的摄像机采集图像的速率为每秒钟25-30帧，即采样率为25-30hz；而心率变异性测量要求的时间精度为20ms，根据奈奎斯特采样定律，采样率不能低于100hz；因此，普通摄像机难以满足心率变异性测量的精度要求。

技术实现要素：

本发明针对现有的非接触式测量方法难以满足心率变异性测量要求的缺陷，提出一种非接触式心率变异性测量方法。克服了上述三方面的缺陷：在心率分析时不完全依赖于频谱分析，实现了心率变异性的计算；对盲源分离后的分量进行了进一步的滤波，提高了心率变异性计算的准确性；使用高速摄像机，采样频率可以满足心率变异性测量精度要求。

具体步骤如下：

步骤一、针对某检测者，使用高速彩色摄像机非接触式的采集该检测者面部的区域视频图像；

采集过程在明亮的环境中进行，摄像机采集速率不低于每秒钟100帧，采集过程不短于3分钟。

步骤二、根据每帧视频图像，分别提取人体面部特征点，并根据特征点选取各图像的目标区域；

对每一帧图像进行人脸识别，并提取特征点；

提取过程为：用视频图像中的每个点分别与人脸特征dlib跨平台计算机视觉库中的特征点进行对比，衡量这两个点周围的领域图像的相似度，最后，提取视频图像中相似度最高的点，即为人脸图像中的特征点。

提取的特征点包括每侧脸颊边缘各2个、鼻翼每侧各1个以及鼻梁顶端1个；将每侧面部鼻梁顶端、鼻翼一侧及脸颊边缘的特征点依次相连，得到一个四边形；人脸两侧的四边形区域定义为目标区域。

步骤三、对每帧视频图像的目标区域内的所有像素点的红色、绿色和蓝色三个通道的信号进行独立成分分析，并选取心率相关分量k(t)；

具体步骤为：

步骤301，针对每帧视频图像，对目标区域内所有像素点的红色、绿色和蓝色的通道信号值分别求取平均值。

即：

r(i)为目标区域内的第i个像素点的红色通道信号值；g(i)为目标区域内的第i个像素点的绿色通道信号值；b(i)为目标区域内的第i个像素点的蓝色通道信号值。为目标区域内的所有像素点的红色通道信号值的平均值；为目标区域内的所有像素点的绿色通道信号值的平均值；为目标区域内的所有像素点的蓝色通道信号值的平均值；

步骤302，将所有视频图像计算所得的三个通道信号强度平均值分别按视频采集的时间排列，得到三个时间序列r(t),g(t),b(t)。

r(t)是所有视频图像的红色通道信号的平均值按视频采集的时间排列得到的时间序列；g(t)是所有视频图像的绿色通道信号的平均值按视频采集的时间排列得到的时间序列；b(t)是所有视频图像的蓝色通道信号的平均值按视频采集的时间排列得到的时间序列；

步骤303，对三个时间序列进行独立成分分析ica，分解得到三个独立成分分量，记为i(t),j(t),k(t)；

计算公式如下：

其中a是一个3*3的矩阵，通过使i(t),j(t),k(t)相关性最小的最优化算法获得。

步骤304，将三个独立成分分量分别与g(t)分量进行相关性分析，把相关性最大的独立成分分量选出，并命名为心率相关分量k(t)。

相关性最大的分量包含脉搏波信号。

步骤四、对心率相关分量k(t)进行频谱分析和滤波处理，提取脉搏波特征点，得到逐拍心动周期；

首先，对心率相关分量k(t)的频谱分析，提取0.5-3hz的频段中能量最大点所对应的频率，得到心率频率fc；

然后，根据心率频率fc，选取截止频率为(fc-0.3)～(fc+0.3)hz的带通滤波器，对心率相关分量k(t)进行频域滤波，结果记为kc(t)。

进而，对滤波结果kc(t)进行脉搏波特征点提取，并依次记录脉搏波峰值点出现的时刻p1,p2,……，pn。

提取出现频率约为fc的局部极大值点作为脉搏波峰值点。

最后，根据相邻极大值点出现的时刻计算逐拍心动周期y＝{y1,y2,……,yn-1}。

计算公式为：y1＝p2-p1,y2＝p3-p2,……,yn-1＝pn-pn-1；

步骤五、根据逐拍心动周期计算心率变异性的时域和频域参数。

首先，对逐拍心动周期包含的随机噪声进行平滑滤波；

针对当前时刻的逐拍心动周期yn，平滑滤波公式如下：

j表示采样间隔的数量；

然后，将极大值出现的时刻p1,p2,……，pn作为横坐标，将对应的平滑滤波后的逐拍心动周期y1,y2,……,yn-1作为纵坐标，生成新的时间序列y(p)。

进而，对时间序列y(p)进行三次样条曲线插值，使其采样频率等于视频采集的频率fs，插值后的时间序列记为y(t)，称为心率变异性信号。

最后，对心率变异性信号y(t)进行时域和频域分析，得到心率变异性参数。

对心率变异性信号进行时域分析，得到参数：meannn、sdnn、rmssd和sdsd；

其中，n是参与心率变异性计算的心动周期的个数；μ为相邻y(n)差值的均值。y(n)是心率变异性信号y(t)中n时刻的心动周期值，y(n+1)表示心率变异性信号y(t)中n+1时刻的心动周期值。

对心率变异性信号进行频域分析，得到plf、phf、lfhfratio、vlf、lf和hf。

其中，对心率变异性信号y(t)进行傅里叶变换，频谱图当中的0.003-0.04hz的能量综合为vlf；0.04-0.15hz的能量总和为lf；0.15-0.4hz的能量总和为hf；

其它参数计算方法如下：

plf＝lf/(lf+hf)

phf＝hf/(lf+hf)

lfhfratio＝lf/hf

步骤六、根据心率变异性的时域和频域参数对心理压力进行评价和分析。

本发明的优点在于，可以无接触式的测量心率变异性，其测量精度可与心电测量所得的心率变异性相当；由于非接触式的特点，可以减少对检测对象的影响，获得较为客观的心理压力的评价。

附图说明

图1为本发明一种非接触式心率变异性测量方法的整体流程示意图；

图2为本发明提取的人脸特征点和目标区域示意图；

图3为本发明的心率相关分量及其频谱示意图；

图4为本发明的心率相关分量滤波后脉搏波特征点提取结果。

具体实施方式

下面将结合附图和实施例对本发明作进一步的详细说明。

本发明提供了一种非接触式心率变异性测量方法，可在受试者无察觉的情况下，对其心率变异性进行测量，从而较为真实的反映其心理压力水平，用于心理健康状态的评价。包括6个环节：采集视频图像；选取脸颊作为目标区域；独立成分分析得到心率相关分量；频域滤波并提取逐拍心动周期；计算心率变异性参数以及评价心理压力大小。

所述选取脸颊作为目标区域是指：对每一帧采集所得的人脸图像进行人脸识别，提取眼、鼻、口等人脸特征点，根据这些特征点选取脸颊区域。所述的独立成分分析得到心率相关分量是指：对脸颊区域图像的每个像素，统计红、绿、蓝三个通道信号强度的均值，对三个通道的信号进行独立成分分析，得到三个独立分量，其中与绿通道相关性最高的分量被认为是心率相关分量。所述频域滤波并提取逐拍心动周期是指：对心率相关分量进行频域滤波，并在滤波结果上进行脉搏波特征点的提取，根据相邻特征点之间的时间差，计算出逐拍心动周期。所述的计算心率变异性参数是指：根据逐拍心动周期计算可得心率变异性时域参数4个，频域参数6个。所述的评价心理压力大小是指：根据心率变异性的时域参数和频域参数，可以得到心理压力分析结果。

如图1所示，具体步骤如下：

步骤一、针对某检测者，使用高速彩色摄像机非接触式的采集该检测者面部的区域视频图像；

受试者以正面面向摄像机，若摄像机可以连续采集受试者正面图像3分钟以上，该方法即可计算出受试者的心率变异性，从而用于受试者心理压力大小的评价。在采集过程中，只要受试者面部不移出摄像机采集范围，就可以完成计算。期间受试者可以轻微晃动面部，改变表情或说话，上述面部运动都不会影响计算结果。

通过高速彩色摄像机采集面部的视频图像，采用红色、绿色、蓝色三原色制式，采集过程应在明亮的环境中进行，如果在较暗的环境中，应辅助以稳定光源照射受试者的面部，来保证采集到足够亮度的面部图像。摄像机采集速率应不低于每秒钟100帧，采集过程应不短于3分钟；所采集的视频，以图像格式存储于上位机，并进行后续的分析处理；图像数量应不少于1800张。

步骤二、根据每帧视频图像，分别提取人体面部特征点，并根据特征点选取各图像的目标区域；

对每一帧图像进行人脸识别，并提取人脸的特征点；

提取过程为：基于人脸特征dlib跨平台计算机视觉库，对每一帧图像进行匹配运算。用视频图像中的每个点分别与人脸特征dlib跨平台计算机视觉库中的特征点进行对比，衡量这两个点周围的领域图像的相似度，最后，提取视频图像中相似度最高的点，即为人脸图像中的特征点。

按此方法提取的特征点如图2所示，包括每侧脸颊边缘各2个、鼻翼每侧各1个以及鼻梁顶端1个，共计7个特征点；将每侧面部鼻梁顶端、鼻翼一侧及脸颊边缘的特征点依次相连，得到一个四边形；这7个特征点可围成2个区域；人脸两侧的四边形区域定义为目标区域。

步骤三、对每帧视频图像的目标区域内的所有像素点的红色、绿色和蓝色三个通道的信号进行独立成分分析，并选取心率相关分量k(t)；

具体步骤为：

步骤301，针对每帧视频图像，对目标区域内所有像素点的红色、绿色和蓝色的通道信号值分别求取平均值。

目标区域内的每一个像素点，都有红色、绿色和蓝色三个通道信号，记为r(i),g(i),b(i)。将目标区域内的所有像素点的r(i),g(i),b(i)三个值分别求取平均值，即：

r(i)为目标区域内的第i个像素点的红色通道信号值；g(i)为目标区域内的第i个像素点的绿色通道信号值；b(i)为目标区域内的第i个像素点的蓝色通道信号值。为目标区域内的所有像素点的红色通道信号值的平均值；为目标区域内的所有像素点的绿色通道信号值的平均值；为目标区域内的所有像素点的蓝色通道信号值的平均值；

步骤302，将所有视频图像计算所得的三个通道信号强度平均值分别按视频采集的时间排列，得到三个时间序列r(t),g(t),b(t)。

r(t)是所有视频图像的红色通道信号的平均值按视频采集的时间排列得到的时间序列；g(t)是所有视频图像的绿色通道信号的平均值按视频采集的时间排列得到的时间序列；b(t)是所有视频图像的蓝色通道信号的平均值按视频采集的时间排列得到的时间序列；

步骤303，对三个时间序列进行独立成分分析ica，分解得到三个独立成分分量，记为i(t),j(t),k(t)。

如图3所示，将r(t),g(t),b(t)视为由脉搏波信号和其它噪声线性混叠的结果，可根据独立成分分析ica，将它们分解为三个独立分量，记为i(t),j(t),k(t)；

计算公式如下：

其中a是一个3*3的矩阵，可通过使i(t)，j(t)，k(t)相关性最小的最优化算法获得。

步骤304，将三个独立成分分量分别与g(t)分量进行相关性分析，把相关性最大的独立成分分量选出，并命名为心率相关分量k(t)。

按照上述假设，这三个分量中必有一个为脉搏波信号；而已知脉搏波与绿光的相关性最大，因此将上述三个分量分别与g(t)进行相关性计算，观察其结果，其中相关性最大的分量，即包含脉搏波信号，假设为k(t)，可视为心率相关分量。

步骤四、对心率相关分量k(t)进行频谱分析和滤波处理，提取脉搏波特征点，得到逐拍心动周期；

首先，对心率相关分量k(t)的频谱分析，提取0.5-3hz的频段中能量最大点所对应的频率，得到心率频率fc；

然后，根据心率频率fc，选取截止频率为(fc-0.3)～(fc+0.3)hz的带通滤波器，对心率相关分量k(t)进行频域滤波，结果记为kc(t)。

进而，对滤波结果kc(t)进行脉搏波特征点提取，提取出现频率约为fc的局部极大值点作为脉搏波峰值点，并依次记录脉搏波峰值点出现的时刻p1,p2,……，pn。

脉搏波峰值点提取结果如图4所示。

最后，根据极大值点出现的时刻计算逐拍心动周期y＝{y1,y2,……,yn-1}。

将相邻两个脉搏波峰值点出现的时刻相减，得到：y1＝p2-p1,y2＝p3-p2,……,yn-1＝pn-pn-1；

步骤五、根据逐拍心动周期计算心率变异性的时域和频域参数。

首先，由于得到的逐拍心动周期包含随机噪声，可对其进行平滑滤波，方法是使用邻域为5的均值滤波器，针对当前时刻的逐拍心动周期yn，平滑滤波公式如下：

j表示采样间隔的数量；

然后，将极大值出现的时刻p1,p2,……，pn作为横坐标，将对应的平滑滤波后的逐拍心动周期y1,y2,……,yn-1作为纵坐标，生成新的时间序列y(p)。

进而，对时间序列y(p)进行重采样，使其采样频率等于视频采集的频率fs，插值后的时间序列记为y(t)，称为心率变异性信号。

插值方法可选用三次样条插值法。

最后，对心率变异性信号y(t)进行时域和频域分析，得到心率变异性参数，用于心理压力大小的评价和分析。

对心率变异性信号进行时域分析，得到参数：meannn、sdnn、rmssd和sdsd等参数；

其中，n是参与心率变异性计算的心动周期的个数；μ为相邻y(n)差值的均值。y(n)是心率变异性信号y(t)中n时刻的心动周期，y(n+1)表示心率变异性信号y(t)中n+1时刻的心动周期。

对心率变异性信号进行频域分析，得到plf、phf、lfhfratio、vlf、lf和hf等参数。

其中，对y(t)进行傅里叶变换，频谱图当中的0.003-0.04hz的能量综合为vlf；0.04-0.15hz的能量总和为lf；0.15-0.4hz的能量总和为hf；

其它参数计算方法如下：

plf＝lf/(lf+hf)

phf＝hf/(lf+hf)

lfhfratio＝lf/hf

步骤六、根据心率变异性的时域和频域参数构造svm分类器，对心理压力进行评价和分析。

综上所述，本发明无接触式心率变异性测量方法采用了全新的心率变异性测量技术，摆脱了传统心率变异性测量依赖于接触式心电测量的方法，较少了测量过程中对人的干扰，从而可以更真实的反映人的心理压力状态。且本发明无需受试者配合测量，只需要其面向摄像头3分钟以上即可完成测量。如果通过微型摄像机等方式进行拍摄，即可在受试者不知情的条件下进行心率变异性的测量和心理压力的评价，因此大大拓展了心率变异性测量的应用。