一种心磁信号背景噪声S变换域去除方法及系统

摘要

本发明公开了一种心磁信号背景噪声S变换域去除方法及系统，属于心脏磁信号处理技术领域，包括以下步骤：S1：采集噪声数据点；S2：对噪声进行S变换；S3：绘制最大值二维时频曲面；S4：确定区分临界界面；S5：比较心磁信号二维时频矩阵与区分临界界面的关系；S6：进行反S变换。本发明采用S变换方法在时频域内实现背景噪声的消除，具有实现简单，背景噪声去除效果好，运行速度快，无需购置外加硬件设施等特点；该方法通用性强，适合于无屏蔽无遮掩情况下受高信噪比干扰的心磁信号数据的背景噪声去除去及心磁信号参数精确提取，值得被推广使用。

说明书

技术领域

本发明涉及心脏磁信号处理技术领域，具体涉及一种心磁信号背景噪声S变换域去除方法及系统。

背景技术

心脏疾病是危害人类健康的主要疾病之一。心脏功能的早期检查是预防和治疗心脏疾病的关键。目前，心脏功能的检测和诊断方法主要采用心电图(Electrocardiograph，ECG)来检测和分析心脏电活动。心电图通常由P波、QRS波群和ST波段组成，其中P波代表左、右心房除极化过程，QRS波群代表左、右心室除极化等过程，S-T波段反映心肌细胞复极化状态，QRS波群蕴涵的瞬态信息对识别心脏系统的电磁状态变化极其重要，也能提供心脏系统某些病态的先兆特征。

随着磁传感器件的发展，超导量子干涉器作为无损电磁探测的探测器，具有高灵敏度，宽带宽，良好线性程度等各方面的优势，广泛地应用于心磁信号测量与材料无损检测中。心磁图(Magnetocardiograms，MCG)的波形和心电图一样，即由P波、QRS波和ST波段组成，当人体某部位发生变化时，体内的电流及磁场会发生变化，通过对心脏磁场的变化而探测心脏电生理是MCG在临床使用的理论基础。MCG与ECG相比有它的独到之处，ECG与心肌细胞的跨膜电流强度有关，主要依赖于心肌组织的电导率和靠近体表的体电流强度；MCG依赖于磁场分布，它源于整体电流，特别是心肌组织的涡旋电流源和具有高电流密度的体电流部分。胸前磁场的实际测量结果表明，MCG可以提供心脏系统电磁场源性质的更多重要细节，特别是涡旋电流源信息，而这些细节不可能被ECG测到。然而，在无遮掩无屏蔽的普通的检测室或实验室情况下，高灵敏度的高温射频超导量子干涉仪采集到的心磁信号中存在大量的背景噪声，其主要成份为50Hz工频干扰信号及其它外部环境噪声。低频磁噪声干扰约为几个纳特斯拉(10

发明内容

本发明所要解决的技术问题在于：如何从强背景噪声中有效提取心磁信号及相关波段特征，提供了一种心磁信号背景噪声S变换域去除方法。

本发明是通过以下技术方案解决上述技术问题的，本发明包括以下步骤：

S1：采集噪声数据点

实时采集仪器所在环境下的背景噪声N个数据点；

S2：对噪声进行S变换

对采集的N个数据点背景噪声采用S变换方法变换成二维时频矩阵；

S3：绘制最大值二维时频曲面

连续采集多次，计算每次采集信号对应的二维时频矩阵，得到多组二维时频矩阵，计算统计过程中多组每一个时频点时频矩阵的最大值，得到统计情况下的背景噪声S变换二维时频矩阵，绘制矩阵得到最大值二维时频曲面；

S4：确定区分临界界面

取最大值曲面的设定置信度曲面作为噪声与心磁信号的区分临界界面；

S5：比较心磁信号二维时频矩阵与区分临界界面的关系

对采集的心磁信号采用S变换方法变换成二维时频矩阵，当采集的心磁信号二维时频矩阵大于区分临界界面时，保留心磁信号的二维时频矩阵；当心磁信号的二维时频矩阵小于或等于区分临界界面时，将心磁信号的二维时频矩阵置为0，对心磁信号的二维时频矩阵进行修改；

S6：进行反S变换

对修改后的心磁信号二维时频矩阵进行反S变换，得到去除背景噪声后的心磁信号。

更进一步地，在所述步骤S1中，所述仪器所在环境为心磁信号测量前的仪器所在环境。

更进一步地，在所述步骤S2中，S变换为采用高斯窗函数的短时傅里叶变换，其表达式为：

更进一步地，在所述步骤S4中，设定置信度为95％。

更进一步地，在所述步骤S6中，反S变换为采用高斯窗函数的短时傅里叶逆变换，其表达式为：

本发明还提供了一种心磁信号背景噪声S变换域去除系统，采用上述的背景噪声去除方法对心磁信号中的背景噪声进行去除，包括：

噪声采集模块，用于实时采集仪器所在环境下的背景噪声N个数据点；

噪声变换模块，用于对采集的N个数据点背景噪声采用S变换方法变换成二维时频矩阵；

绘制模块，用于连续采集多次，计算每次采集信号对应的二维时频矩阵，得到多组二维时频矩阵，计算统计过程中多组每一个时频点时频矩阵的最大值，得到统计情况下的背景噪声S变换二维时频矩阵，绘制矩阵得到最大值二维时频曲面；

临界界面确定模块，用于取最大值曲面的设定置信度曲面作为噪声与心磁信号的区分临界界面；

比较模块，用于对采集的心磁信号采用S变换方法变换成二维时频矩阵，当采集的心磁信号二维时频矩阵大于区分临界界面时，保留心磁信号的二维时频矩阵；当心磁信号的二维时频矩阵小于区分临界界面时，将心磁信号的二维时频矩阵置为0，对心磁信号的二维时频矩阵进行修改；

反变换模块，用于对修改后的心磁信号二维时频矩阵进行反S变换，得到去除背景噪声后的心磁信号；

中央处理模块，用于向其他模块发出指令，完成相关动作；

所述噪声采集模块、噪声变换模块、绘制模块、临界界面确定模块、比较模块、反变换模块均与中央处理模块电连接。

本发明相比现有技术具有以下优点：该心磁信号背景噪声S变换域去除方法，采用S变换方法在时频域内实现背景噪声的消除，全过程可以软件自己运行，无需人为多余操作；该方法实现简单，背景噪声去除效果好，运行速度快，无需购置外加硬件设施；该方法通用性强，适合于无屏蔽无遮掩情况下受高信噪比干扰的心磁信号数据的背景噪声去除去及心磁信号参数精确提取，值得被推广使用。

附图说明

图1是本发明实施例中心磁信号背景噪声去除方法的流程示意图；

图2a是本发明实施例中噪声S域内的时频特性图；

图2b是本发明实施例中的平均功率谱(K＝1000)；

图2c是本发明实施例中平均功率谱与95％置信线；

图2d是本发明实施例中比值随频率变化关系图；

图3是本发明实施例中含背景噪声的心磁信号；

图4是本发明实施例中含背景噪声心磁信号S变换后三维时频图；

图5是本发明实施例中噪声与心磁信号区分临界面示意图；

图6是本发明实施例中去除背景噪声后的心磁信号。

具体实施方式

下面对本发明的实施例作详细说明，本实施例在以本发明技术方案为前提下进行实施，给出了详细的实施方式和具体的操作过程，但本发明的保护范围不限于下述的实施例。

本实施例做了如下实验内容：

1、噪声S变换域统计特性

设均值为0，方差为σ

E[n(t)n(u)]＝σ

其中，E表示期望，u为时间，δ(t)为冲激函数。

设高斯噪声n(t)经S变换后的时频矩阵为S(τ,f)，则高斯白噪声的平均功率谱为：

其中，*表示取复共轭。

当噪声为有色噪声时，采用一阶马尔可夫模型可得有色噪声x(t)的时域叠代关系为：

x(t)＝wx(t-1)+n(t) (3)

其中，w为滞后因子。

相应的有色噪声傅立叶变换下的平均功率谱为：

其中，m为离散频点。

为了研究有色噪声在S变换时频域内的功率谱特性，首先对有色噪声在S变换二维时频域下的特性分布作定性的分析，猜测其总体分布规律，然后利用Monte Carlo方法对其规律进行验证。

由于噪声服从高斯分布，采用快速S变换算法对其进行变换后，其时频特性仍然服从高斯分布。具体推导过程如下：

第一步，计算输入高斯噪声的傅立叶变换，变换结果的实部和虚部均服从高斯分布；

第二步，对变换结果进行频移，频移结果的实部和虚部仍服从高斯分布；

第三步，计算频移结果与定值高斯窗的乘积，乘积结果不影响实部和虚部分布；

第四步，对乘积结果的反傅立叶变换，变换结果实部和虚部仍服从高斯分布。

由于平均功率谱|S(τ,f)|

采用Monte Carlo方法验证的主要步骤如下：

1)在给定显著性水平α下，求置信区间1-α；

2)构造均值为μ，方差为σ

3)对n(t)或x(t)进行S变换，可得瞬时功率谱|S(τ,f)|

4)重复步骤2)、步骤3)K次，计算平均功率谱

5)计算位于置信区间1-α所在的曲面与平均功率谱曲面的比值

2、基于噪声统计特性的背景噪声去除方法

设含噪信号为：x(t)＝u(t)+n(t)，其中u(t)为信号源，n(t)为噪声。为了利用噪声在S变换时频域内的统计特性，首先做如下的原假设和被假设：

H

其中，H

拒绝域为：

将公式(5)代入公式(6)，可得1-α置信区间内的原假设条件为：

不等式左边为满足原假设在1-α置信区间的时频域，即x(t)为噪声所在的时频区域：

同理，x(t)满足被假设为信号源所在的时频区域为：

根据上述实验步骤，在显著性水平α＝0.05，置信区间为95％的情况下，实验统计采用均值为0，方差为0.1的随机高斯白噪声，采样点数N＝128，重复实验统计次数K＝1000。单次随机噪声所对应的S变换三维时频等高图如图2(a)所示，由图2(a)可知，噪声在S域内幅值的大小随着频率的增大而增大，这一特征与式(2)所示的理论结果一致。随着统计次数的增加，噪声的平均功率谱曲面逐渐平滑，如图2(b)所示。

当τ＝32时，计算统计所得的平均功率谱与频率的关系及95％置信线处的功率谱与频率的关系如图2(c)所示。95％所在置信线所处的功率谱与平均功率谱之间的比值如图2(d)所示。由图2(d)可知，两者的比值在3左右波动。这一实验结果与

本实施例提供了一种心磁信号背景噪声去除方法，包括以下步骤：

(1)：实时采集仪器所在环境下的背景噪声N个数据点；

(2)：对采集的N个数据点背景噪声采用S变换方法变换成二维时频矩阵S[N，M]，M＝N/2，取整数；

(3)：连续采集1000次，计算每次采集信号对应的二维时频矩阵，可得1000组二维时频矩阵S1[N，M]、S2[N，M]、…、S1000[N，M]；

(4)：计算1000统计过程中每一个时频点时频矩阵的最大值{max{S1[1，1]、S2[1，1]、…、S1000[1，1]}，max{S1[2，1]、S2[2，1]、…、S1000[2，1]}，…，max{S1[N，M]、S2[N，M]、…、S1000[N，M]}可得1000次统计情况下的噪声S变换二维时频矩阵，绘制矩阵可得最大值二维时频曲面；

(5)：在置信度为95％的情况下，取最大值曲面95％置信度曲面作为噪声与心磁信号的区分临界界面S[Nmax，Mmax]；

(6)：对采集的心磁信号采用S变换方法变换成二维时频矩阵S[N1，M1]，当采集的心磁信号二维时频矩阵S[N1，M1]大于区分临界界面S[Nmax，Mmax]，保留S[N1，M1]；当S[N1，M1]小于或等于S[Nmax，Mmax]时，将S[N1，M1]置为0，对S[N1，M1]进行修改；

(7)：对修改后二维时频矩阵S[N1，M1]进行反S变换，即可得去除背景噪声后的心磁信号。

所述的心磁信号背景噪声S变换域去除设计方法，其特征在于，所述步骤(1)中，所述仪器所在环境为心磁信号测量所在同一环境。

在本实施例中，所述步骤(2)中S变换为采用高斯窗函数的短时傅里叶变换，其表达式为：

在本实施例中，所述步骤(5)中置信度为95％是指1000次统计过程中，每一个时频点有950次的值都落在区分临界界面S[Nmax，Mmax]以下。

在本实施例中，所述步骤(7)中反S变换为采用高斯窗函数的短时傅里叶逆变换，其表达式为：

本发明方法适合任意无磁屏蔽室条件下利用高温射频超导量子干涉仪系统所测量到的心磁信号。为了方便说明，下面以理想心磁信号受10dB高斯白噪声背景噪声干扰为例进行去除背景噪声进行相关设计说明。

设计时，主要包括以下几点：首先，实时采集仪器在未接触人体时所在环境下的背景噪声，经S变换后可得二维时频矩阵，如图2(a)所示；重复采集1000次，可得1000次背景噪声经S变换后的二维时频矩阵，计算每个时频点在1000次的最大值，可得如图2(b)所示曲面。在95％置信度下，可得最大值对应曲线与95％置信线如图2(c)所示，计算最大值对应曲线与95％置信线比值可得图2(d)所示曲线。当心磁信号在无磁屏蔽室条件下实时采集时，会受到电力线干扰及环境噪声干扰，实时采集的心磁信号如图3所示，其对应的S变换时频图如图4所示，根据假设检验方法，对比心磁信号三维时频图与1000次噪声统计置信曲面，如图5所示，当心磁信号S变换幅值小于临界面是，则认为是背景噪声；反之则为心磁信号，将小于临界曲面的心磁信号S变换时频矩阵数据置0；最后采用反S变换，即可恢复去除背景噪声后的心磁信号，如图6所示。心磁信号背景噪声99％能有效去除掉。

综合上述仿真结果及分析表明，该心磁信号背景噪声去除方法是可行的。

综上所述，上述实施例的心磁信号背景噪声S变换域去除方法，采用S变换方法在时频域内实现背景噪声的消除，全过程可以软件自己运行，无需人为多余操作；该方法实现简单，背景噪声去除效果好，运行速度快，无需购置外加硬件设施；该方法通用性强，适合于无屏蔽无遮掩情况下受高信噪比干扰的心磁信号数据的背景噪声去除去及心磁信号参数精确提取，值得被推广使用。

尽管上面已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本发明的限制，本领域的普通技术人员在本发明的范围内可以对上述实施例进行变化、修改、替换和变型。