一种分数阶90°移相器及其qrs检波方法
【专利摘要】本发明公开了一种分数阶90°移相器及其QRS检波方法，属于数字信号处理技术领域。基于Grünwald?Letnikov分数阶微积分，定义分数阶中心差分算子，其本质为一维分数阶90°移相器，其相频特性与传统一阶导数相同，恒为90°相移，由于其幅频特性为频率的分数阶次幂函数，可以通过调节分数阶次在检测精度与抗噪性之间取得一定平衡，将心电信号x(n)通过所设计的分数阶90°移相器，进行R波检测，在已确定R波的基础上，位于R波之前的极小值点为Q波，位于R波之后的极小值点为S波。本发明计算简便，对受噪严重、波形异常的心电信号具有较好的检波效果。
【专利说明】
-种分数阶90°移相器及其QRS检波方法
技术领域
[0001 ]本发明公开了一种分数阶90°移相器及其QRS检波方法，属于数字信号处理技术领 域。
【背景技术】
[0002] -个正常的屯、电信号是由一系列的波形组成的，主要包含P波、QRS复合波、T波、ST 段W及听皮。QRS复合波的检测是屯、电信号波形检测中的首要问题。精确可靠的QRS波检测是 诊断屯、电疾患的首要依据，是其他波形诊断的重要参考。QRS复合波的检测主要是基于其形 态、频率特性与其他波形的差异来实现的。近年来，QRS复合波检测手段主要有：阔值法，小 波分析法，模板匹配法，神经网络法等等。
[0003] QRS复合波检测的首要问题是波峰值点的检测，即R波的检测。由于R波具有幅度大 和斜率局的特点，因此最为常见的R波检测方法有幅度法、面积法和斜率法。幅度法和面积 法在遇到非常规的屯、电信号，如高T波、高P波时容易出错;斜率法则优于运两种方法。斜率 法，又称差分阔值法，是将滤波后的信号进行差分运算，再采用预先给定的阔值进行判别。 另一种阔值法则考虑屯、电信号QRS复合波的宽度检测，当检测到QRS时输出一个脉冲信号， 并与某个事先设定的阔值进行比较，从而检测出QRS复合波。总的来说，阔值法设计思路简 单，计算量小，易于工程实现，但同时阔值的设置存在随机性，精度不能得到一致的保证。
[0004] 小波变换可W将信号分解成不同尺度下的子信号，属于局部化分析方法。它研究 信号频谱和小波变换尺度之间的关系，通过将屯、电信号中不同波形和干扰源在频谱上进行 对比，从频谱的分布上检测出QRS波。小波变换通过反复实验，学习并优化阔值的选择，运用 多标准判定，保证了较高的检测精度。小波变换的另一个重要特点是抗噪性强，减少漏检和 误检的概率，其缺点是计算量大，不适宜实时工程应用。
[0005] 模板匹配法将QRS波群、T波等近似认为是单独的模板，通过幅值归一化处理，信号 中不同的幅频分量对应各个模板的能量分布，通过不同幅频分量的模板之间的匹配来检测 出QRS波群，将对QRS波群的检测转化为模板间的匹配鉴别。该方法对幅值判断和频率判断 进行了统一，优点是抗干扰能力强，精度高，缺点是耗时长，计算量大。
[0006] 神经网络法首先对信号进去降噪处理，得到噪声较小段信号的QRS复合波模板作 为参考模板。然后按次序在屯、电信号中找到匹配度较高的屯、电信号QRS复合波。神经网络的 特点是判别速度快，但需要较长的训练时间。
【发明内容】
[0007] 发明目的：为了克服现有技术中存在的不足，本发明提供一种分数阶90°移相器及 其QRS检波方法，该移相器相频特性与传统一阶导数相同，恒为90°相移，由于其幅频特性为 频率的分数阶次幕函数，可W通过调节分数阶次在检测精度与抗噪性之间取得一定平衡。 [000引为实现上述目的，本发明采用的技术方案为:一种分数阶90°移相器，包括中屯、差 分Griinwald-Letnikov微积分算子及；，所述中屯、差分G;riinwald-Letn;Lkov微积分算子口;包 括前向Griinwald-Letnikov微积分算子公/、.后向Griinwald-Letnikov微积分算子公及比 例运算器；所述前向Griinwald-LetnAov微积分算子滤波器巧'、后向Griinwald-LetnAov微 积分算子滤波器相互并联在一起之后再与比例运算器串联;所述输入信号分别通过前向 Griinwa 1 d-Le tn i koV微积分算子滤波器马V、后向Griinwa 1 d-Letn i koV微积分算子滤波器 进行滤波，滤波后的信号通过比例运算器调节后得到输出信号。
[0009] 优选的：所述中屯、差分Griinwald-LetnAov微积分算子滤波器1?与前向Griinwald- LetnAov微积分算子滤波器巧和后向Grilnwald-LetnAov微积分算子滤波器巧:的关系为：
[0010]
[00川其中，V辛化，nez。
[0012] 优选的:所述输出信号为：
[0013] Wei。）=X(ei。）· isgn( ω ) I ω IV
[0014] 其中，V为分数阶微积分滤波器的阶次，取任意非偶实数;X(ei")为输入信号；i为 虚数单位，ω为信号频率。
[0015] 优选的：所述中屯、Griinwald-Letnikov微积分算子滤波器巧、前向Griinwald- LetnAov微积分算子滤波器及后向Grilnwald-Letnikov微积分算子滤波器巧：均通过 差分方程的形式来构建。
[0016] 一种基于分数阶90°移相器的QRS检波方法，将屯、电信号x(n)通过所设计的分数阶 90°移相器，进行R波检测，在已确定R波的基础上，位于R波之前的极小值点为Q波，位于R波 之后的极小值点为S波。
[0017]优选的：将屯吨信号x(n)通过所设计的分数阶90°移相器，利用阀值法滤除屯吨信 号中的高频小幅毛刺，然后进行R波检测。
[0018]优选的:所述R波检测，其步骤为：1、区间法选出备选的正波和负波;2、排除法剔除 不符合要求的极值点，选择相邻两个极值点中正确的极值点；3、配对法进一步排除T波和P 波;4、坐标修正法确定原始信号中R波的正确坐标。
[0019]优选的：所述QRS波形检测，Q波和S波检测:位于R波之前的极小值点为Q波，之后的 极小值点为S波。
[0020] 本发明提供的分数阶90°移相器及其QRS波形检测，相比现有技术，具有W下有益 效果：
[0021] (1)分数阶90°移相器，其相频特性与传统一阶导数相同，恒为90°相移，然而由于 其幅频特性为频率的分数阶次幕函数，可W通过调节分数阶次在检测精度与抗噪性之间取 得一定平衡。
[0022] (2)设计方法简便，通过对信号进行分数阶微积分运算，分别实现对信号的分数阶 前向滤波和后向滤波，通过并联设计最终得到分数阶90°移相器，设计方法简便且算法效率 局。
[0023] (3)设计灵活性大，通过调节分数阶微积分的阶次，可W依据实际需要的设计指标 选择合适的阶次。
[0024] 综上所述:将本发明基于分数阶90°移相器的检波方法应用于屯、电信号QRS波形检 测，可取得较好的检测效果。设计方法简便，设计灵活性大，通过调节分数阶微积分的阶次， 可W依据实际需要的设计指标选择合适的阶次。本发明的分数阶90°移相器的检波方法可 应用于非平稳信号处理，尤其适用于批信号处理领域。
【附图说明】
[0025] 图1为本发明分数阶90°移相器的流程图；
[0026] 图2为卷积模板曲线图，其中V分别取-0.5和0.5两种情况，m= 100;
[0027] 图3为本发明分数阶90°移相器的QRS波形检测流程图；
[0028] 图4为分数阶90°移相器作用于输入信号前后的变化图；
[0029] 图5为阔值处理前后信号对比图；
[0030] 图6为检测R波第一步，采用区间法选出备选的正波和负波后的曲线图；
[0031] 图7为检测R波第二步，采用排除法剔除不符合要求的极值点后的曲线图；
[0032] 图8为检测R波第Ξ步，采用配对法进一步排除T波和P波后的曲线图；
[0033] 图9为检测R波第四步，采用坐标修正法确定原始信号中R波的正确坐标的曲线图；
[0034] 图10为美国麻省理工学院MIT-BIH屯、电数据库第106号数据检测出的QRS波群。
【具体实施方式】
[0035] 本发明所解决的技术问题在于获得一种比现有的基于一阶差分检波方法更理想 的分数阶信号处理算法，所设计的滤波器其本质为一维分数阶90°移相器，其相频特性与传 统一阶导数相同，恒为90°相移，由于其幅频特性为频率的分数阶次幕函数，可W通过调节 分数阶次在检测精度与抗噪性之间取得一定平衡。在具体说明本
【发明内容】
之前，有必要对 本说明书所用符号涵义进行简要说明。
[0036] 根据一阶左导数的定义，可得到：
[0037]
[003引其中，f(t)为输入信号。
[0039] 将巧定义为一阶η重左导数，可得到：
[0040]
[0041 ]其4
为二项式系数。
[0042] 当η取分数时，该二项式系数可替换夫
痒中V为任 意非整数，由此我们定义V阶前向Griinwald-Letnikov微积分算子用符号表示，对函数f (t)进行V阶前向61'化1￥31(1-1^61：]1化0￥微积分，得到下式：
[0043]
[0044] 其中，t和a分别是前向Griinwald-Letnikov微积分算子的上限和下限，该微积分算 子将η阶微分与η阶积分运算统一用一个微积分算子来表示，且将整数域的运算推广到非整 数域。当ν>0时，上式表示微分运算，当ν<0时，上式表示积分运算，微积分算子由此得名， 信号进行V阶前向Griinwald-LetnAov微积分的傅里叶变换表示为：
[0045] 巧巧7脚二快卽/(衍）
[0046] 因此，V阶前向Griinwald-LetnAov微积分算子巧"的频域函数为：
[0047]
[004引对应的幅频和相频响应分别为：
[0049]
[0050] 根据一阶右导数的定义，可得到：
[0化1 ]
[0052]其中f(t)为输入信号。
[0化3]将定义为一阶η重右导数，可得到：
[0化4]
[0化5]其中
为二项式系数。
[0056]当η取分数时，该二项式系数可替换夫
，其中V为任 意非整数。由此我们定义V阶后向Griinwald-Letnikov微积分算子用符号巧'表示。对函数f (t)进行V阶后向Griinwald-LetnAov微积分，得到下式：
[0化7]
[0化引其中，b和t分别是后向Griinwald-Letnikov微积分算子的上限和下限。信号进行V 阶后向Griinwald-LetnAov微积分的傅里叶变换表示为：
[0化9]
[0060] 因此，V阶后向Griinwald-LetnAov微积分算子巧的频域函数为：
[0061]
[0062] 对应的幅频和相频响应分别为：
[0063]
[0064] 前向61'化1*日1(1-1^61：]1化〇￥微积分算子巧'和后向6111]1*日1(1-1^61：]111?)￥微积分算子巧* 可被视为两种特殊的滤波器，其幅值响应随信号频率和分数阶微积分阶次变化而变化，信 号通过前向61'化1*曰1(1-1^61：]114〇￥微积分算子坏相角滞后^端"批），而信号通过后向61'0 nwa 1 d-Letn;?koV微积分算子巧相角超前^sgn(份)。
[0065] 下面结合附图对本发明作更进一步的说明。
[0066] 本发明设计的分数阶90°移相器，如图1所示，包括中屯、Griinwald-Letnikov微积分 算子滤波器巧，所述中屯、61'11脯日1(1-1^61：]114〇￥微积分算子滤波器'￡^包括前向61'11脯日1(1- LetnAov微积分算子滤波器巧'、.后向Griinwald-LetnAov微积分算子滤波器及比例运 算器，所述前向Griinwald-LetnAov微积分算子滤波器〇/'、后向Griinwald-LetnAov微积分 算子滤波器相互并联在一起之后再与比例运算器串联;所述输入信号分别通过前向Gru nwald-Letnikov微积分算子滤波器公巧日后向Griinwald-Letnikov微积分算子滤波器巧进 行滤波，滤波后的信号通过比例运算器
调节后得到输出信号。
[0067] 用数学表达式来说明1?与〇/、〇巧勺关系，如下式 [006引
[0069] 其中V辛化，neZ。
[0070] 因此，对函数f (t)进行V阶中屯、Griinwald-LetnAov微积分可W表示为：
[0073] 假设将中屯、Griinwald-Le化化OV微积分算子应用于一个有限长的输入信号序列 f(t)，相应的傅里叶变换为
[0074]
[007引因此，中屯、Griinwald-LetnAov微积分算子巧的频率函数为
[0076]
[0077] 对应的幅频特性和相频特性分别为：
[007引
[0079]具体实施时，由于数字信号处理是基于离散信号进行的，因此我们需将上述时域 函数进行离散化操作，即用差分方程的形式来构建数字信号滤波器。取h= 1，且当或 h 足够大的时候，可省略公式中的取极限运一步骤。前向GriinwaId-LetnAov微积分算 子0/、后向Griinwald-Letnikov微积分算子巧"和中屯、Griinwald-Letnikov微积分算子巧的 差分近似定义如下：
[0084] 为了计算出前向Grilnwald-Letnikov微积分算子邱的差分方程，取卷积模板 =[.0,a…口2 口 11.].，模板长度为m+1。那么，前向Griinwald-Letnikov微积分算子巧的差分方 程可写为
[0085] Δ;/的=/約 * 始J
[00化]同理，为了计算出后向Griinwald-LetnAov微积分算子公J的差分方程，取卷积模板 MfV = 口口.1 口2 ''··c?。,],模板长度为m+1。那么，后向Griinwald-LetnAov微积分算子.巧·的差分方 程可写为
[0087] Δ?,''/.(/) = ,/ν)*Μ;
[0088] 当中屯、差分算子1?应用在有限长度的信号f(t)上时，差分方程形式为
[0089]
[0090] 其中
卷积模板长度为2m+l。
[0091] 图2所示为卷积模板曲线图，其中V分别取0.5和-0.5两种情况，m= 100,即卷积 模板长度为201，两条卷积模板曲线奇对称且在两端趋于0,表现了90°移相的通性，且保证 了适度的计算复杂度。
[0092] 图3为QRS波形检测流程图，首先，我们将原始屯、电信号x(n)通过所设计的分数阶 90°移相器得到信号y(n)作为后续实验对象，见图4,接着，利用常见的阔值法滤除屯、电信号 中的高频小幅毛躁，使得滤波后的输出信号较平稳的区间更加光滑，便于后续R波检测，图5 是阔值作用后的屯、电信号，两条水平线分别为正阔值与负阔值。阔值法作用后得到的信号， 实验中记为信号X壯。
[0093] 接着，第Ξ个环节通过4个步骤得到R波的正确坐标。
[0094] 图6为检测R波的第一步，用区间法选出备选的正波和负波。将信号等距离划分区 间，逐个区间检索出正、负极值。每个区间包含50个采样点（由于信号采样频率为360Hz，因 此50个采样点对应的区间时长大约为ο . 14s)。如果区间内有非零的极值点（绝对值），则把 运些极值点选出来并记录坐标。R-R峰值间隔医学统计上大约为0.4-1.2秒之间，每个R-R峰 间，大约有3-9个捜索区间不等。图6(b)中，圆圈表示的是区间法确定的正负极值点。可W看 到，通过区间法初步筛选出来的极值点很多，且由于区间起始点的随机性和区间大小的设 定具有先决性，检测出的"假"极值点很多。
[0095] 图7为第二步，用排除法剔除不符合要求的极值点，选择相邻两个的极值点中正确 的极值点。测量相邻两个非零正极值点间的距离。选择出所有相邻距离小于区间长度的幅 值较小的点，并删除。
[0096] 图8为第Ξ步，配对法进行进一步排除T波和P波。为了防止偶尔出现的高幅值T波 和P波，需要使用配对法进行排除，找到正确的正负极值对。测量两个相邻的正负极值点（正 负极值点放在一起考虑）。有两种情况：1)若距离小于等于区间长度，即50个采样点或0.14 秒，且正极值点出现在前，负极值点在后，则两个极值点为一个正确的正-负配对。2)若距离 大于区间长度，两个极值点不构成一对配对，需要删除。图8(b)中显示了将错误归为R波的T 波排除后的效果。
[0097] 图9为第四步，坐标修正法确定原始信号中R波的正确坐标。一旦我们在差分信号 中找到R波，捜索将鉴于此时确定的信号坐标，回到原始信号，将此坐标作为新的捜索中屯、 点，向此坐标前后拓宽一定的长度，得到新的捜索窗，窗口大小大概取在125ms~250ms。通 过运个步骤的处理，弥补移相和滤波过程中造成的信号点R波相位的偏移。图9(b)展示的是 在原始信号x(n)中使用本方法检测出来的R波，图中X指原始屯、电信号，R指检测出来的R波。 [0〇9引图10为美国麻省理工学院MIT-BIH屯吨数据库第106号数据检测出的QRS波群，由 于Q波和S波分别是位于R波前后的负向极值点波，且距离R波较近，QRS波群的长度一般为 0.08~0.12s(即29-44个采样点之间），因此WR波为中屯、，前后各22个采样点作为检测窗 口。位于R波之前的极小值点为Q波，之后的极小值点为S波。
【主权项】
1. 一种分数阶90°移相器，其特征在于:包括中心Griinwald-Letnikov微积分算子滤波 器 1，所述中心6101歷&1(1-1^1：11丨1^〇￥微积分算子代包括前向6101歷31(1-1^1：11丨1^〇￥微积分算 子1)丨和后向61'11簡31(1-1^1：11丨1^〇￥微积分算子样:，:所述输入信号分别通过前向61'1111￥31(1-Letnikov微积分算子Μ和后向Griinwa 1 d-Letni kov微积分算子乃丨进行滤波，滤波后的信号 通过比例运算器调节后得到输出信号。2. 根据权利要求1所述的分数阶90°移相器，其特征在于:所述中心Griinwald-Letnikov 微积分算子滤波器·与前向Griinwald-Letnikov微积分算子滤波器.1>/以及后向Griinwald-Letnikov 微积分算子滤波器 的关 系为： 其中，v关2η，ηΕΖ。3. 根据权利要求2所述的分数阶90°移相器，其特征在于:所述输出信号为： Y(eliJ) = X(eliJ) · isgn( ω ) I ω Iv 其中，v为分数阶微积分滤波器的阶次，取任意非偶实数;X(e1M)为输入信号;i为虚数单 位，ω为信号频率。4. 根据权利要求3所述的分数阶90°移相器，其特征在于:所述中心Griinwald-Letnikov 微积分算子滤波器乃】、前向6101簡&1(1-1^1：11丨1?)￥微积分算子滤波器￡ )，17以及后向61'11簡31(1-Letnikov微积分算子滤波器!>；：均通过差分方程的形式来构建。5. -种基于权利要求1所述的分数阶90°移相器的QRS检波方法，其特征在于:将心电信 号x(n)通过所设计的分数阶90°移相器，进行R波检测，在已确定R波的基础上，位于R波之前 的极小值点为Q波，位于R波之后的极小值点为S波。6. 根据权利要求5所述的QRS检波方法，其特征在于:所述R波检测的步骤为：1、区间法 选出备选的正波和负波;2、排除法剔除不符合要求的极值点，选择相邻两个极值点中正确 的极值点;3、配对法进一步排除T波和P波;4、坐标修正法确定原始信号中R波的正确坐标。7. 根据权利要求5所述的QRS检波方法，其特征在于:所述QRS波形检测、Q波和S波检测： 位于R波之前的极小值点为Q波，之后的极小值点为S波。
【文档编号】A61B5/0472GK105870549SQ201610192289
【公开日】2016年8月17日
【申请日】2016年3月30日
【发明人】王建宏, 方靖淮, 陈玉娟, 于志华, 金晶亮, 殷姝
【申请人】南通大学