基于MUSIC的按需心率估计算法

想象未来几十年后的世界，您的孙子们可能不知道医院这个词，所有健康信息都是通过传感器远程记录和监测。想象您的家里配备了不同的传感器来测量空气质量、温度、噪声、光照和气压，并且根据您的个人健康信息，系统调整相关环境参数以优化您的家居环境。在实现美好未来的道路上，ADI 公司处于一个独特的有利位置，通过提供相互补充的传感器、软件和算法来增加其在数字健康市场的份额。

心率 (HR) 监测是许多现有可穿戴产品和临床设备的关键特性。这些设备一般测量光电容积脉搏波 (PPG) 信号，为获得该信号，须利用 LED 照射人体皮肤，然后用光电二极管测量血流引起的反射光强度变化。PPG 信号形态与动脉血压 (ABP) 波形相似，这使得该信号成为受科学界欢迎的非侵入性心率监测工具。PPG 信号的周期性与心脏节律相对应。因此，可以根据 PPG 信号估算心率。然而，受血液灌流不良、环境光线以及最重要的运动伪像 (MA) 1 的影响，心率估算性能会降低。业界已提出许多信号处理技术来消除 MA 噪声，包括 ADI 公司的运动抑制和频率跟踪算法，通过使用一个靠近 PPG 传感器放置的三轴加速度传感器来实现。当没有运动时，最好能有一个按需算法来向跟踪算法提供快速且更精确的心率估算。本文改造了多信号分类 (MUSIC) 频率估计算法，以利用 ADI 医疗健康手表平台，根据手腕上的 PPG 信号实现高精度按需心率估算，图 1 所示为其框图。该图的细节将在后面的内容中说明。

ADI 医疗健康手表提供的 PPG 信号概述

当 LED 发光时，血液和组织会吸收不同数量的光子，导致光电检测器检测到不同的结果。光电检测器测量血液脉动的变化并输出一个电流，该电流随后经放大和滤波以供进一步分析。

图 1.利用腕上 PPG 信号的基于 MUSIC 的按需心率估计算法。

图 2a 显示了一个由交流 (ac) 和直流 (dc) 分量组成的一般 PPG 信号。PPG 波形的直流分量检测组织、骨骼和肌肉反射的光信号，以及动脉和静脉血液的平均血容量。交流分量则表示心动周期的收缩期和舒张期之间发生的血容量变化，交流分量的基频取决于心率。图 2b 是来自 ADPD107 手表的 PPG 信号，这在之前的《模拟对话》文章中已介绍过。ADI 多感知手表的目标是测量人体手腕上的多项生命体征。ADI 手表有 PPG、心电图 (ECG)、皮肤电活动 (EDA)、加速度 (ACC) 和温度传感器。本文仅关注 PPG 和 ACC 传感器。

现在我们仔细看看 PPG 和 ABP 波形的相似之处。ABP 波形是由于左心室射出血液造成的。主压力沿全身血管网流动并到达多个部位，动脉阻力和顺应性的显著变化引起反射。第一个部位是胸主动脉和腹主动脉之间的接合处，其引起第一次反射，通常称为收缩晚期波。第二个反射部位是腹主动脉和髂总动脉之间的接合处。主波被再次反射回来，产生一个很小的下降，称为重搏切迹，这可以在第一次和第二次反射之间观察到。还有其他较小的反射，这些反射在 PPG 信号中被平滑掉2。本文的重点是心率估计，其仅取决于 PPG 信号的周期性，此算法不考虑 PPG 的确切形态。

图 2a.含交流和直流部分的典型 PPG 信号

图 2b.ADI 医疗保健手表PPG信号

PPG 信号预处理

PPG 信号易受周边组织的不良血液灌流和运动伪像的影响是众所周知的 1。为将这些因素的影响降至最小，以免干扰随后的 PPG 分析和心率估计，须有一个预处理阶段。需要一个带通滤波器来消除 PPG 信号的高频成分(如电源)和低频成分(如毛细血管密度和静脉血容量的变化、温度变化等等)。图 3a 显示了滤波之后的 PPG 信号。使用一组信号质量指标来找到适合于按需算法的 PPG 信号第一个窗口。第一次检查涉及 ACC 数据和 PPG 信号，以确定是否能检测到一段无运动的数据，然后衡量其他信号质量指标。如果三个方向上存在高于 ACC 数据绝对值的特定阈值的运动，则按需算法会拒绝根据这样的数据窗口进行估计。下一信号质量检查是基于数据段特征的某种自相关。图 3b 显示了经滤波的 PPG 信号的一个自相关例子。可接受信号段的自相关表现出如下特性：具有至少一个局部峰值，并且对应于最高可能心率的峰值不超过某一数量;局部峰值从高到低递减，间隔时间递增;以及其他一些特性。仅计算与有意义的心率(位于 30 bpm 到 220 bpm 范围内)相对应的间隔时间的自相关。

当有足够的数据段连续通过质量检查时，算法的第二阶段就会使用基于 MUSIC 的算法算出准确的心率。

图 3a.经过带通滤波的图 1b 中 PPG 信号

图 3b.图 2a 中信号图的自相关

基于 MUSIC 的按需心率估计算法

MUSIC 是一种基于子空间的方法，使用谐波信号模型，可以高精度地估算频率3。对于受到噪声破坏的 PPG 信号，傅立叶变换 (FT) 可能表现不佳，因为我们需要的是高分辨率心率估计算法。此外，FT将时域噪声均匀分布到整个频域中，限制了估算的确定性。使用 FT 很难在较大峰值附近观察到较小峰值 4。因此，在本研究中，我们使用基于 MUSIC 的算法进行心率的频率估计。MUSIC 背后的关键思想是噪声子空间与信号子空间正交，所以噪声子空间的零点会指示信号频率。下面的步骤是这种心率估计算法的总结：

• 从数据中删除平均和线性趋势

• 计算数据的协方差矩阵

• 对协方差矩阵应用奇异值分解 (SVD)

• 计算信号子空间阶数

• 形成信号或噪声子空间的伪谱

• 找出 MUSIC 伪谱的峰值作为心率估计值

MUSIC 必须应用奇异值分解，并且必须在整个频率范围内搜索频谱峰值。我们来看一些数学算式，以使上述步骤更清晰。假设经滤波的PPG信号有一个长度为 m 的窗口，表示为 xm 且 m≤L (其中 L 为给定窗口中经滤波 PPG 信号的总样本数)。那么，第一步是形成样本协方差矩阵，如下所示：

然后对样本协方差矩阵应用 SVD，如下所示：

其中，U 为协方差矩阵的左特征向量，Λ 为特征值的对角矩阵，V 为右特征向量。下标 s 和 n 分别代表信号和噪声子空间。正如之前提到的，使用信号已经通过信号质量检查阶段的先备知识，对基于 MUSIC 的算法进行修改以用于心率估计，因此预处理步骤之后，信号中唯一存在的频率成分是心率频率。接下来形成信号和噪声子空间，假设模型阶数只包含一个单音，如下所示：

Us = U(1: p, :); Un = U(p + 1: end, :)

Us = U(1: p, :); Un = U(p+ 1: end, :)

其中 p = 2 为模型数。仅考虑有意义心率限值内的频率。这会大大减少计算量，使嵌入式算法的实时实现成为可能。搜索频率向量定义为：

其中，k 为心率目标频率范围内的频点，L 为 xm (t) 中数据的窗口长度。然后，下面的伪谱使用噪声子空间特征向量找出 MUSIC 的峰值，如下所示。

4, which shows a sharp peak at 1.96 Hz, and which translates to 117.6 bpm HR.

这里使用伪谱一词，是因为它表明所研究信号中存在正弦分量，但它不是一个真正的功率谱密度。图 4 显示了基于 MUSIC 的算法处理 5 秒数据窗口得到的示例结果，在 1.96 Hz 处有一个很陡的峰值，换算为心率是 117.6 bpm。

图 4.使用 PPG 数据的基于 MUSIC 估计的一个示例

基于 MUSIC 的按需心率估计算法的结果

我们已经在一个包含 1289 个测试案例 (data1) 的数据集上测试了该算法的性能，并且在数据开始时，测试对象被要求静止。表 1 给出了基于 MUSIC 算法的结果，并指出估计的心率是否在参考 (ECG) 的 2 bpm 和 5 bpm 精度范围内，以及估计时间的第 50 百分位数(中值)和第 75 百分位数。表 1 中的第二行显示了对于一个包含 298 个测试案例 (data2) 的数据集，存在周期性运动(如步行、慢跑、跑步)时该算法的性能。通过检测运动，如果任一数据被视为不可靠而遭到拒绝，或者是认为不受运动影响而精确估算得到心率，则认为该算法是成功的。在内存使用方面，假设缓冲区大小为 500 (即 100 Hz 时为 5 秒)，对于目标频率范围 (30 bpm 至 220 bpm)，所需总内存约为 3.4 kB，每次调用花费 2.83 周期。

表 1.基于 MUSIC 的按需心率估计算法的性能数值

结语

基于 MUSIC 的按需算法是 ADI 公司医疗保健业务部门生命体征监测小组提出的众多算法之一。在我们医疗健康手表中使用的按需算法与这里讨论的基于 MUSIC 的方法不同，前者的计算成本较低。ADI 公司为传感器(嵌入式)和边缘节点提供软件和算法功能，使其从数据中获取有价值的信息，仅将最重要的数据发送到云端，让我们的客户和合作伙伴可以在本地做出决策。我们选择应用的标准是，其成果对于我们的客户来说非常重要，并且我们拥有独特的测量专业技术。本文只是对 ADI 公司研发的算法的简单介绍。凭借我们在传感器设计方面的现有专业知识，以及我们在生物医学算法开发(包括嵌入式和云计算)方面的努力，ADI 公司将拥有独特的优势来为全球医疗健康市场提供最先进的算法和软件。