在波浪中检测图案
我试图从心电图读取图像,并检测其中的每一个主波(P波,QRS波群和T波)。 现在我可以读取图像,并得到像(4.2; 4.4; 4.9; 4.7; …)这样的代表心电图数值的vector,问题的一半是什么。 我需要一个能够遍历这个向量的algorithm,并且能够检测到这些波的每一个何时开始和结束。
这是一个图表的例子:
如果他们总是有相同的大小,那将是容易的,但它不是像它的工作,或者如果我知道心电图会有多less波,但也可能会有所不同。 有没有人有一些想法?
谢谢!
更新
我想要实现的例子:
鉴于波
我可以提取vector
[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;]
例如,我想检测一下
P波在[19 – 37]
QRS波群[51 – 64]
等等…
我要做的第一件事就是看看那里已经有了什么 。 的确,这个具体问题已经被深入研究。 以下是一些非常简单的方法的简要概述: 链接 。
我也必须回答另一个答案。 我研究信号处理和音乐信息检索。 从表面上看,这个问题的确出现了类似于发病检测,但是问题的上下文并不相同。 这种types的生物信号处理,即P,QRS和T相的检测,可以利用这些波形中的每一个的特定时域特性的知识。 在MIR中的发现检测不是,真的。 (至less不可靠)
一种适用于QRS检测(但不一定用于音符起始检测)的方法是dynamic时间翘曲。 当时域特性保持不变时,DTW可以非常好地工作。 这里是一个简短的IEEE论文,使用DTW解决这个问题: 链接 。
这是一个很好的IEEE杂志文章,比较了许多方法: 链接 。 你会看到许多常见的信号处理模型已经被尝试过了。 浏览论文,尝试一个你在基础层面上理解的东西。
编辑:浏览这些文章后,基于小波的方法似乎对我来说是最直观的。 DTW也能很好地工作,那里有DTW模块,但小波方法对我来说似乎是最好的。 其他人通过利用信号的衍生物来回答。 我的第一个链接检查了1990年以前的方法,但是我怀疑它们不像更现代的方法那么强大。
编辑:我会尝试给我一个简单的解决办法,当我有机会,但我认为小波适合在这里的原因是因为它们是参数化的各种形状很有用,无论时间或幅度缩放 。 换句话说,如果信号具有相同的重复时间形状,但是在不同的时间尺度和振幅下,小波分析仍然可以将这些形状识别为相似(粗略地说)。 另外请注意,我有点过滤银行到这个类别。 类似的事情。
这个难题的一部分是“ 发现检测 ”,并且已经编写了许多复杂的algorithm来解决这个问题。 这里有更多关于onset的信息。
下一个是海明距离 。 这种algorithm允许您进行模糊比较,input是2个数组,输出是一个整数“距离”或两个数据集之间的差异。 数字越小,2就越相似。 这是非常接近你所需要的,但它不是确切的。 我继续前进,并对Hamming Distancealgorithm做了一些修改来计算一个新的距离,它可能有一个名字,但我不知道它是什么。 基本上,它将数组中每个元素之间的绝对距离相加,并返回总和。 这是在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,这是这两个数组之间的距离。
现在把这些东西放在一起。 您可以使用起始点检测来查找数据集中所有波的起点。 然后,您可以循环这些位置,将每个波与样本P波进行比较。 如果你打了一个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波和t波紧挨着,但是因为这个序列是可预测的,所以它会更容易阅读。
另一方面,这个问题可能有一个基于微积分的解决scheme。 然而,在我的心中,曲线拟合和积分使这个问题更加混乱。 我写的距离函数将find非常相似的区域差异 ,减去两条曲线的积分。
也许可以牺牲开始计算而有利于一次迭代1点,从而执行O(n)距离计算,其中n是图中点的数量。 如果你有一个所有这些距离计算的清单,并且知道那里有50个pQt序列,那么你就可以知道50个最短距离,在p波的所有位置都不重叠 。 答对了! 这是为了简单? 然而,由于距离计算的数量增加,折衷是效率的损失。
你可以使用互相关 。 拿出每个模式的模型样本,并将它们与信号相关联。 你会得到相关性高的高峰。 我希望用这种技术提取qrs和t波,效果很好。 之后,您可以通过在qrs之前的相关信号上查找峰值来提取p波。
互相关是一个相当容易实现的algorithm。 基本上:
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中查找峰值(例如,超过阈值的值)
我要做的第一件事就是简化数据。
分析绝对数据,而不是分析从一个数据点到下一个数据点的变化量。
这是一个快速的class轮, 分离的数据作为input,并输出该数据的增量。
perl -0x3b -ple'( $last, $_ ) = ( $_, $_-$last )' < test.in > test.out
根据您提供的数据运行它,这是输出:
0; 0; 20; 0; 0; 1; -1; -1; 0; 0; 0; 0; 1; 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; -1; 0; 1; 2; 4; 6; 9; 7; 7; 6; -4; -6 -8; -7; -5; -4; 0; 1; 0; – 1; 1; 1; 0; 1; 0; 1; 1; 0; 0; 0; 0; 0; 0; 1; 0; 1; 1; -1; 0; 1; 0; 0; 0 ; 1; 0; 1; 1; 2; 2; 0; 1; 1; 1; 1; 1; 1; 1; 0; 0; 1; 0; 0; -1 -2 -1 -2 -2 -2 -2 ; 0; 1; 1; 0; 1; 0; 1; 0; 1; 0; 1; 1; 0; 0; 0; 0; 0; 0; 0; 0; 1; 1; 1; 0; 0; 0; 0; 0; 0; 0; 0; 1; 1; 1; 0; 0; 1; 0; 0; 0; 0; 0; 0; 0; 1; 1; 0; 0; 0; 0 -1; 0; 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; 1; 5; 5; 7; 8; 9; 4; -7 -5 -8 ; -7 -6; -2; -1; 0; 0; 0; 0; 0; 1; 0; 0; 1; -1; 0; 1; 0; 1; 1; 0; 0; 0 ; 1; 0; 0; 0; 1; 0; 1; 0; 0; 0; 1; 1; 0; 2; 1; 1; 1; 1; 1; 1; 1; 1; 1; -1; 1; 0; 0; 1; -2 -2 -2 -2 -1 0 -1 -2 -1 0 -1 -1; 0; 1; -1; 1; 0; 1; 1; -1; 1; 0; 1; 0; 0; 0; 1; 1; 0; 0; 1; 0; 1; 0; 1; 0; 0; 1; 1; 0; 1; 1; 0; 1; 0; 0 ; 0; 0; 1; -1; 0; 1; 1; 0; 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; 0; 0; 0; 1; 0; 0; 4; 3; 9; 8; 11; 4; 1-5; -6 -8 -8 -4 -2 -2; 0; 0; 0; 1; 1; 0; 0; 1; 0; 0; 1; -1; 0; 1; 0; 0; 0; 1; -1; 0; 1; 1; 0; 0; 0; 0; 1; 0; 1; 0; 1; 2; 1; 1; 2; 0; 1 ; 1; 1; 1; 0; 0; 1; 1; 0; 0; -35; 0; 0; 0;
在上面的文本中插入的新行并不是最初出现在输出中的。
完成之后,findqrs复合体是微不足道的。
perl -F';' -ane'@F = map { abs($_) > 2 and $_ } @F; print join ";", @F'< test.out
;; 20 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; 4; 6; 9; 7; 7; 6; -4; -6 -8; -7; -5; -4;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;; 5; 5; 7; 8; 9; 4; -7 -5 -8; -7 -6
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;; 4; 3; 9; 8; 11; 4; 1-5; -6 -8 -8 -4;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; – 35 ;;;
20
和-35
数据点是从原始数据开始和结束0
。
要find其他数据点,您将不得不依赖模式匹配。
如果你看第一个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 __________/
你会发现类似于p波的t波。 主要的区别是什么时候发生。
这应该是足够的信息,让你开始。
这两个单线只是例子,不build议日常使用。
另外两个尖峰和山谷是否也是复杂的?
在我头顶,我认为你需要做的是计算每个点的这个图的斜率。 那么你还需要看看斜率有多快(二阶导数)。 如果你有一个突然的变化,那么你知道你已经达到了某种尖锐的高峰。 当然,你想要限制对变化的检测,所以你可能想要做一些事情,比如“如果斜率在时间间隔T内改变了X”,那么你就不会拿起图中的小颠簸。
自从我做了任何math以来,这已经有一段时间了……这似乎是一个math问题;)噢,我还没有做过任何types的信号分析:)。
只是补充一点。 你也可以尝试我认为的信号平均。 例如,平均最后3或4个数据点。 我想你也可以检测到突然的变化。
我不是这个问题的专家,但是从更一般的知识来看,我只能说我是头一回:让我们假设你知道QRS波群(或者其他的特征之一,但是在这个例子中我将使用QRS波群)大概发生在一个长度为L的固定时间段。我想知道你是否可以把它看作是一个分类问题,如下所示:
- 将信号分成长度为L的重叠窗口。每个窗口都有或没有完整的QRS波群。
- 傅里叶变换每个窗口。 您的function是每个频率的信号强度。
- 训练决策树,支持向量机等在一些手被注释的数据。
一种很可能会产生好结果的方法是曲线拟合:
- 将连续的波分成多个区间(可能最好是在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)
是一个简单的函数,可以用来模拟三个波的每一个。 -
使用拟合algorithm(例如,Levenberg-Marquardtalgorithmhttp://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为每个intervall的数据集。
参数t_p,t_qrs和t_p是你感兴趣的。
我认识一个在这个领域工作的人 。 在这里你可以find他的出版物清单 。 如果我没有记错的话,他会使用隐马尔可夫模型来可靠地检测已知形状的训练集中的波浪,但是您可以在论文中find更多细节。
你已经有了很多好的答案。 我只是感到惊讶没有人build议从PhysioToolkit的 ' WFDB软件包 ',特别是ecgpuwave
这是一个奇妙的问题! 我有几个想法:
dynamic时间扭曲可能是一个有趣的工具。 您可以为您的三个类build立“模板”,然后使用DTW可以看到您的模板和信号的“块”之间的相关性(将信号分解成例如0.5秒的比特,即0-5。 1-.6 .2-.7 …)。 我使用类似于加速度计数据进行步态分析的工作,效果相当好。
另一种select是组合的信号处理/机器学习algorithm。 再次将信号分解为“块”。 再次制作“模板”(每个类需要十几个),对每个块/模板进行FFT ,然后使用朴素贝叶斯分类器 (或另一个ML分类器,但NB应该将其切割)分类为你的三个class。 我也对步态数据进行了尝试,能够获得高达98%的精度,并且可以调用相对复杂的信号。 让我知道这是如何工作,这是一个非常令人兴奋的问题。
“ 小波变换 ”可能是一个相关的关键字。 我曾经参加过一个使用这种技术在嘈杂的心电图中检测不同心跳相位的人的演示。
就我的理解而言,它有点像傅里叶变换,但是使用(缩放的)a(在你的情况下是心跳形)脉冲的拷贝。
小波已被certificate是定位这种types的数据峰值的最佳工具,其中峰值是“不同的大小” – 小波的缩放属性使其成为这种types的多尺度峰值检测的理想工具。 这看起来像一个非平稳的信号,所以使用DFT并不是一些人所build议的正确的工具,但是如果这是一个探索性的项目,你可以看看使用信号的频谱(估计主要是使用自相关的FFT信号。)
这里有一篇很好的论文来回顾几种高峰检测方法 – 这将是一个很好的开始。
– 保罗
我的答案是关于检测时间序列数据中的模式的类似问题,在这里 – https://stackoverflow.com/a/11903770/1149913 – 并包括Python代码。
我的方法是“切换自回归隐马尔可夫模型”(谷歌短语的一些相关的出版物)。
我还没有彻底地读过对方的答案,但我已经扫描了他们,我注意到没有人build议看傅立叶变换来分割这些波。
对我来说,这似乎是math中谐波分析的一个明确的应用。 可能有几个微妙的地方,我可能会失踪。
离散傅立叶变换系数为您提供组成您的离散时间信号的不同正弦分量的幅度和相位,这实际上是您想要查找的问题状态。
虽然我可能会错过这里的东西
首先,标准心电图波形的各个组成部分可能会从任何给定的图上缺失。 这样的阴谋通常是不正常的,通常表示某种问题,但不能保证他们在那里。
第二,认识它们就像科学一样,尤其是在出现问题的地方。
我的方法可能是尝试训练一个neural network来识别组件。 你会给它的前30秒的数据,归一化,所以最低点是在0和最高点在1.0,它将有11个输出。 非exception评级的输出将是最后10秒的权重。 0.0现在是-10秒,现在是1.0。 产出将是:
- 最近的P波开始的地方
- 最近的P波结束了
- 最近一次P波exception评级,一个极端“缺席”。
- 最近的QRS综合症开始的地方
- 最近一次QRS波群的Q部分变成R部分。
- 最近一次QRS波群的R部分变成S部分。
- 最近一次QRS复杂结束的地方。
- 最近一次QRS综合征的exception评分,其中一个极端是“缺席”。
- 最近的T波开始的地方
- 最近的T波结束了。
- 最近一次T波exception评级,一个极端“缺席”。
我可以仔细检查一下其他types的分析,或者使用其他types的分析以及neural network的输出来给出答案。
当然,这个neural network的详细描述不应该被视为规定性的。 我敢肯定,我不一定会select最理想的产出,我只是想提一些他们可能的想法。
Alaor,这是一个非常好的问题。 你有没有想过neural network? 你有足够的信息,你可以教他们识别波(训练他们)或使用模糊逻辑…我会去neural network的方法。
希望这可以帮助!
你尝试过biosppy吗? 如果你使用Python,你可以直接使用它。 如果没有,您可以阅读他们的源代码并修改为您的语言。
使用BioSPPY
目前还不可能实施T波分析,因为目前它只包含R波分析。 例如Tstart Tpeak Tend
不会自动实现
你将不得不使用你自己的分析。
我的build议是尝试和实施下面的方法
http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3201026/
这是我最近发现的,发现很有趣
另一个值得关注的t波分析方法是ECGlib团队
http://ieeexplore.ieee.org/document/6713536/
希望这可以帮助