如何探测(心电扫描仪)波的模式?

我试着从一个心电扫描仪上读出一张图片,然后检测其中的每一个主波(P 波、 QRS 波和 T 波)。我可以读取图像并得到一个矢量(如 (4.2; 4.4; 4.9; 4.7; ...))。我需要一个算法,可以通过这个矢量检测每个波的开始和结束时间。举个例子:

alt text

如果他们总是有相同的大小,或者如果我提前知道心电图有多少波,这将是很容易的。考虑到这波浪潮:

alt text

我提取矢量:

[0; 0; 20; 20; 20; 19; 18; 17; 17; 17; 17; 17; 16; 16; 16; 16; 16; 16; 16; 17; 17; 18; 19; 20; 21; 22; 23; 23; 23; 25; 25; 23; 22; 20; 19; 17; 16; 16; 14; 13; 14; 13; 13; 12; 12; 12; 12; 12; 11; 11; 10; 12; 16; 22; 31; 38; 45; 51; 47; 41; 33; 26; 21; 17; 17; 16; 16; 15; 16; 17; 17; 18; 18; 17; 18; 18; 18; 18; 18; 18; 18; 17; 17; 18; 19; 18; 18; 19; 19; 19; 19; 20; 20; 19; 20; 22; 24; 24; 25; 26; 27; 28; 29; 30; 31; 31; 31; 32; 32; 32; 31; 29; 28; 26; 24; 22; 20; 20; 19; 18; 18; 17; 17; 16; 16; 15; 15; 16; 15; 15; 15; 15; 15; 15; 15; 15; 15; 14; 15; 16; 16; 16; 16; 16; 16; 16; 16; 16; 15; 16; 15; 15; 15; 16; 16; 16; 16; 16; 16; 16; 16; 15; 16; 16; 16; 16; 16; 15; 15; 15; 15; 15; 16; 16; 17; 18; 18; 19; 19; 19; 20; 21; 22; 22; 22; 22; 21; 20; 18; 17; 17; 15; 15; 14; 14; 13; 13; 14; 13; 13; 13; 12; 12; 12; 12; 13; 18; 23; 30; 38; 47; 51; 44; 39; 31; 24; 18; 16; 15; 15; 15; 15; 15; 15; 16; 16; 16; 17; 16; 16; 17; 17; 16; 17; 17; 17; 17; 18; 18; 18; 18; 19; 19; 20; 20; 20; 20; 21; 22; 22; 24; 25; 26; 27; 28; 29; 30; 31; 32; 33; 32; 33; 33; 33; 32; 30; 28; 26; 24; 23; 23; 22; 20; 19; 19; 18; 17; 17; 18; 17; 18; 18; 17; 18; 17; 18; 18; 17; 17; 17; 17; 16; 17; 17; 17; 18; 18; 17; 17; 18; 18; 18; 19; 18; 18; 17; 18; 18; 17; 17; 17; 17; 17; 18; 17; 17; 18; 17; 17; 17; 17; 17; 17; 17; 18; 17; 17; 18; 18; 18; 20; 20; 21; 21; 22; 23; 24; 23; 23; 21; 21; 20; 18; 18; 17; 16; 14; 13; 13; 13; 13; 13; 13; 13; 13; 13; 12; 12; 12; 16; 19; 28; 36; 47; 51; 46; 40; 32; 24; 20; 18; 16; 16; 16; 16; 15; 16; 16; 16; 17; 17; 17; 18; 17; 17; 18; 18; 18; 18; 19; 18; 18; 19; 20; 20; 20; 20; 20; 21; 21; 22; 22; 23; 25; 26; 27; 29; 29; 30; 31; 32; 33; 33; 33; 34; 35; 35; 35; 0; 0; 0; 0;]

例如,我想探测:

  • [19 - 37]中的 P 波。
  • [51 - 64]中的 QRS 复合体。
  • 等等。
19651 次浏览

那另外两个尖峰和山谷也是 qrs 复合体吗?

Off the top of my head, I think what you need to do is calculate the slope of this graph at each point. Then you also need to see how quickly the slope changes (2nd derivative???). If you have an abrupt change then you know you've hit some kind of sharp peak. Of course, you want to limit the detection of the change, so you might want to do something like "if the slope changes by X over time interval T", so that you don't pick up the tiny bumps in the graph.

我已经有一段时间没有做过数学题了... ... 这看起来像是一道数学题;)哦,我也没有做过任何类型的信号分析:)。

再加一点。我想你也可以尝试信号平均。例如,平均最后3或4个数据点。你也可以通过这种方式检测到突然的变化。

这个难题的一部分是“ 开始检测”,已经编写了许多复杂的算法来解决这个问题。这里有更多关于 发病的信息。

下一件是 汉明距离。该算法允许进行模糊比较,输入为2个数组,输出为2个数据集之间的整数“距离”或差。数字越小,两者越相似。这非常接近你需要的,但它不精确。我继续对汉明距离算法做了一些修改来计算一个新的距离,它可能有一个名字,但我不知道它是什么。基本上,它将数组中每个元素之间的绝对距离相加,然后返回总距离。下面是 Python 中的代码。

import math


def absolute_distance(a1, a2, length):
total_distance=0
for x in range(0,length):
total_distance+=math.fabs(a1[x]-a2[x])
return total_distance


print(absolute_distance([1,3,9,10],[1,3,8,11],4))

这个脚本输出2,即这两个数组之间的距离。

现在我要把这些碎片拼起来。您可以使用 Onset 检测来查找数据集中所有波的开始。然后你可以循环通过这些位置比较每一个波与样本 P-Wave。如果你撞上了 QRS 波群,距离将是最大的。如果你击中另一个 P 波,数字不会是零,但会小得多。任何 P 波和 T 波之间的距离都会非常小,但是如果你做出以下假设,这就不是问题了:

The distance between any p-wave and any other p-wave will be smaller than the distance between any p-wave and any t-wave.

这个序列看起来像这样: pQtpQtpQt... p-wave 和 t-wave 紧挨着,但是因为这个序列是可预测的,所以更容易阅读。

另一方面,这个问题可能有一个基于微积分的解决方案。然而在我看来,曲线拟合和积分使这个问题更加混乱。我写的距离函数会找到 面积差,它是非常相似的减去两条曲线的积分。

也许可以牺牲起始计算,以便每次迭代1个点,从而执行 O (n)距离计算,其中 n 是图中的点数。如果你有一个所有这些距离计算的列表,并且知道那里有50个 pQt 序列,那么你就会知道50个最短距离,即 不要重叠中所有 p 波的位置。这样的简单性怎么样?然而,由于距离计算次数的增加,这种权衡就是效率的损失。

我不是这个特定问题的专家,但是从我的头脑中更多的一般知识: 让我们假设你知道 QRS 复合体(或者其他特征之一,但是我将在这个例子中使用 QRS 复合体)发生在一个大致固定的时间段 L。我想知道你是否可以把这当作一个分类问题,如下所示:

  1. 将你的信号分割成长度为 L 的重叠窗口。每个窗口都有或没有完整的 QRS 波群。
  2. 你的傅里叶变换是每个频率的信号强度。
  3. 在一些手工注释的数据上训练决策树、支持向量机等。

曲线拟合是一种很可能产生良好结果的方法:

  • 将连续波划分为区间(可能最好在 qrs 复合体的尖峰之间有一半左右的区间边界)。一次只考虑一个间隔。
  • 定义一个模型函数,它可以用来近似所有可能的心电图曲线的变化。这并不像一开始看起来那么困难。 模型函数可以构造成三个函数之和,每个波的起点(t _)、振幅(a _)和宽度(w _)都有参数。

       f_model(t) = a_p   *  f_p  ((t-t_p  )/w_p) +
    a_qrs *  f_qrs((t-t_qrs)/w_qrs) +
    a_t   *  f_t  ((t-t_t  )/w_t)
    

    函数 f_p(t)f_qrs(t)f_t(t)是一些简单的函数,可以用来模拟三个波的每一个。

  • 使用拟合算法(如 Levenberg-MarQuardt 算法 http://en.wikipedia.org/wiki/Levenberg%E2%80%93Marquardt_algorithm)确定每个时间间隔的数据集的拟合参数 a _ p,t _ p,w _ p,a _ qrs,t _ qrs,w _ qrs,a _ t,t _ t,w _ t。

    您感兴趣的参数是 t _ p、 t _ qrs 和 t _ p。

你可以使用 互相关。取每个模式的模型样本,并将它们与信号相关联。相关性高的地方会出现峰值。我希望这种提取 qrs 和 t 波的技术能够得到很好的结果。然后,你可以通过寻找 qrs 之前相关信号的峰值来提取 p 波。

互相关是一个非常容易实现的算法,基本上:

x is array with your signal of length Lx
y is an array containing a sample of the signal you want to recognize of length Ly
r is the resulting correlation


for (i=0; i<Lx - Ly; i++){
r[i] = 0;
for (j=0; j<Ly ; j++){
r[i] += x[i+j]*y[j];
}
}

寻找 r 中的峰值(例如,超过阈值的值)

我要做的第一件事是简化数据。

不要分析绝对数据,而是分析从一个数据点到下一个数据点的变化量。

下面是一个快速的一行程序,它将使用 ;分离的数据作为输入,并输出该数据的增量。

perl -0x3b -ple'( $last, $_ ) = ( $_, $_-$last )' < test.in > test.out

在您提供的数据上运行它,输出如下:

0; 0; 20; 0; 0; -1; -1; -1; -1; -1; -0; 0; 0; 0; -1; 0; 0; 0; 0; 1; 1; 1; 1; 1; 1; 1; 1; 0; 0; 0; -2; -1; -2; -2; -1; -1; -1; -1; 1; -1; -1; 0; -1; 0; -1; 2; 4; 6; 9; 7; 7; 6; -4; -6; -8; -7; -5; -4; 0; -1; 0; -1; 1; 0; -1; 0; 2; 2; 0; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 0; -1; -2; -2; -2; -2; -2; -0; -1; -1; 0; -1; -1; 0; -1; 0; -1; 0; 1; -1; 0; -1; 0; 0; 0; 0; 0; 0; 0; -1; 0; 0; 0; 0; 0; 0; -1; 1; -1; 0; 0; 0; 0; 0; 0; -1; 1; 0; 0; 0; 0; 0; -1; 0; 0; 0; 0; 0; -1; 0; -1; 0; 1; -1; 0; 0; -1; 0; 0; -1; 1; 5; 5; 7; -8; -7; -6; -2; -1; 0; 0; 0; 0; 1; 0; 0; 1; -1; 0; 1; 0; 1; 0; 0; 0; 1; 1; 1; 1; 1; 1; 1; -1; 1; 1; 0; 0; -1; -2; -2; -1; 0; -1; -2; -1; -1; 0; 0; 0; 0; 0; 0; 0; 1; 0; 0; 1; -1; 0; -1; 0; -1; 0; 0; 0; 0; 0; 0; 1; -1; 0; 1; -1; 0; 0; 0; 0; 0; 1; -1; 0; 0; 1; -1; 0; 1; 0; -2; 0; -1; -2; 0; -1; -2; -1; -1; -1; 0; 0; 0; 0; 0; 0; 0; 0; -1; 0; 4; 3; 9; 8; 11; 4; -5; -6; -8; -4; -2; -2; 0; 0; -1; 1; 0; 0; 1; 0; 0; 0; 1; 0; 0; 0; 1; -1; 0; 1; 1; 0; 0; 0; 0; 1; 2; 1; 1; 1; 1; 1; 0; 0; 0; 0;

在上面的文本中插入了一些新行,这些新行最初并没有出现在输出中。


完成这些操作之后,找到 qrs 复合体就变得很简单了。

perl -F';' -ane'@F = map { abs($_) > 2 and $_ } @F; print join ";", @F'< test.out

; ; 20; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; 4; 6; 9; 7; 6; -4; -6; -8; -7; -5; -4;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ;
; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ;

20-35数据点来自以 0开始和结束的原始数据。

要找到其他数据点,你必须依靠模式匹配。


如果你看第一个 p 波,你可以清楚地看到一个模式。

0;0;0;0;0;0;1;0;1;1;1;1;1;1;0;0;2;0;-2;-1;-2;-1;-2;-1;0;-2;-1;1;-1;0;-1;0;0;0;0;
#           \________ up _______/   \________ down _________/

但是要在第二个 p 波上看到图案就不那么容易了。这是因为第二个分布得更广

0;0;0;1;0;1;1;0;1;0;0;1;1;1;0;0;0;-1;-1;-2;-1;0;-2;0;-1;0;-1;0;1;-1;0;0;-1;0;0;0;
#     \________ up _______/       \________________ down ________________/

第三个 p 波比另外两个更不稳定。

0;0;0;0;0;1;-1;0;1;0;0;2;0;1;0;1;1;1;-1;0;-2;0;-1;-2;0;-1;-1;-2;-1;0;0;0;0;0;
#                \_______ up ______/  \__________ down __________/

你会发现 t 波的方式和 p 波相似。主要区别在于它们何时发生。


这些信息应该足够让你开始了。

这两个一行程序仅供例子使用,不推荐日常使用。

要做的第一件事是 看看外面已经有什么了。事实上,这个具体的问题已经得到了大量的研究。下面是一些真正简单的方法的简要概述: 链接

我还要回答另一个问题。我研究信号处理和音乐信息检索。从表面上看,这个问题确实与起始检测相似,但是问题的上下文并不相同。这种类型的生物信号处理,即检测 P、 QRS 和 T 相,可以利用这些波形的 特定时域特性知识。在 MIR 中的起始检测不能,真的。(至少不可靠。)

有一种方法可以很好地检测 QRS 波(但不一定适用于音符起始检测) ,那就是动态时间规整。当时域特性保持不变时,DTW 能够很好地工作。下面是一篇使用 DTW 解决这个问题的 IEEE 简短论文: 链接

这是一篇不错的 IEEE 杂志文章,比较了许多方法: 链接。您将看到已经尝试了许多常见的信号处理模型。浏览一下这篇文章,然后尝试一个你能基本理解的。

编辑: 在浏览了这些文章之后,基于小波的方法对我来说似乎是最直观的。DTW 也可以很好地工作,并且存在 DTW 模块,但是小波方法对我来说似乎是最好的。另一个人则通过利用信号的导数来回答这个问题。我的第一个链接检查了1990年以前的那些方法,但我怀疑它们不如更现代的方法那么健壮。

编辑: 当我有机会的时候,我会尝试给出一个简单的解决方案,但是我认为小波适合这里的原因是因为它们在参数化各种形状时非常有用,不管 时间或振幅缩放是什么。换句话说,如果你有一个信号具有相同的重复时间形状,但在不同的时间尺度和振幅,小波分析仍然可以识别这些形状是相似的(粗略地说)。还要注意的是,我有点把过滤器组归为这一类。差不多。

我没有仔细阅读对方的回答,但我扫描了它们,我注意到没有人建议看傅里叶变换来分割这些波。

在我看来,这似乎是 傅里叶分析在数学中的明确应用。我可能忽略了一些细微之处。

离散傅里叶变换系数给出了构成离散时间信号的不同正弦分量的振幅和相位,这基本上就是你的问题状态,你想要找到的。

不过我可能漏掉了什么。

小波变换”可能是一个相关的关键词。我曾经参加过一个人的演讲,他使用这种技术来检测噪声心电图中不同的心跳阶段。

就我有限的理解而言,它有点像一个傅里叶变换,但是使用(按比例缩放)的副本,在你的例子中是心跳形状的脉冲。

首先,标准心电图波的各个组成部分可以从任何给定的图中丢失。这样的情节通常是不正常的,通常表明一些问题,但你不能保证他们在那里。

其次,识别它们既是科学也是艺术,尤其是在出现问题的情况下。

我的方法可能是训练一个神经网络来识别组件。你会给它前30秒的数据,规范化,所以最低点是0,最高点是1.0,它将有11个输出。不是异常评级的输出将在最后10秒加权。0.0表示距离现在 -10秒,1.0表示现在。产出将是:

  1. 最近的 P 波是从哪里开始的
  2. 最近的 P 波在哪里结束
  3. 最新 P 波异常评级,其中一个极端“缺失”。
  4. 最近的 QRS 综合征是从哪里开始的
  5. 最近 QRS 波群的 Q 部分变成了 R 部分。
  6. 最近 QRS 波群的 R 部分变成了 S 部分。
  7. 最近的 QRS 波群结束的地方。
  8. 最近 QRS 波群的异常评分,其中一个极端是“缺失”的。
  9. 最近的 T 波开始的地方。
  10. 最近一次 T 波结束的地方。
  11. 最近一次 T 波异常评级,其中一个极端是“缺失”的。

我可能会用人们建议的其他类型的分析再次检查这一点,或者使用那些其他类型的分析以及神经网络的输出来给出你的答案。

当然,这种对神经网络的详细描述不应被视为规范性的。例如,我确信我没有必要选择最优的输出,我只是抛出了一些关于它们可能是什么的想法。

小波已被证明是定位这类数据中峰值的最佳工具,因为这类数据的峰值是“不同大小”的——小波的缩放特性使其成为这类多尺度峰值检测的理想工具。这看起来像一个非平稳信号,所以使用 DFT 不会是正确的工具,因为一些人建议,但如果这是一个探索性的项目,你可以看看使用的信号频谱(估计基本上使用 FFT 的自相关信号)

这里 是一篇很好的论文,回顾了几种峰值检测方法——这将是一个很好的起点。

保罗

这是一个很棒的问题! 我有一些想法:

动态时间规整 可能是一个有趣的工具。您可以为您的三个类建立“模板”,然后使用 DTW 可以看到您的模板和信号的“块”之间的相关性(将信号分解为,比如,.5秒,即。0.5.1.6.2.7).我曾经使用过类似的加速度计数据进行步态分析,效果相当不错。

另一种选择是结合信号处理/机器学习算法。再把你的信号分成几块。再次制作“模板”(每个类需要12个左右) ,取每个块/模板的 FFT,然后使用 朴素贝叶斯分类器(或另一个 ML 分类器,但 NB 应该删除它)对你的三个类进行分类。我也在步态数据上做过尝试,并且能够得到98% 以上相对复杂的信号准确率召回率。让我知道这是如何工作的,这是一个非常令人兴奋的问题。