CN103040459A - 一种多通道微弱生理信息记录系统中工频干扰的高保真滤除方法 - Google Patents

Abstract

一种多通道微弱生理信息记录系统中工频干扰的高保真滤除方法，根据多通道测量系统的这个特点，利用信噪比较高的通道估算工频电干扰正弦波的参数。计算出的频率就是所有通道的工频电干扰频率，而计算出来的相位，取决于具体的记录系统中各个通道的空间分布来决定是否反相180度，这可以由本身通道信息估算出的工频电干扰正弦参数来判断。使所有的通道尤其是那些信噪比较低的通道，都得到了和高信噪比通道一样的三个参数中的两个参数的准确估计。对于三个参数中的最后一个参数幅度，可以利用这样一个关系来优化低信噪比的通道中干扰的幅度的估算。所有通道中的工频干扰的幅度的相对比例，是保持不变的。从而提高了低信噪比通道的估算准确性。

Description

一种多通道微弱生理信息记录系统中工频干扰的高保真滤除方法

技术领域

本发明涉及一种工频电干扰的高保真滤除方法，具体涉及一种多通道微弱生理信息记录系统中工频干扰的高保真滤除方法。

背景技术

工频电基波和谐波对测量信息的干扰，是电子时代以来所有在与工频相同频段信号测量中绕不过的一个基本问题。例如我国的50赫兹交流电（北美60Hz）正好落在脑电信息中的Gamma波段。而其几次低次谐波也正好包含在许多生理电活动的范围内，例如，心室晚电位和脑电中的Fast Ripple信息。为了消除工频干扰，已经发明了很多有效的滤波方法。经典的电磁屏蔽室是一种在空间衰减电磁干扰的基本方法，但其高成本、不便移动以及有限的衰减率限制了其应用。在检测输入端利用共模抑制来对抗工频电干扰，例如生物电测量中的差分输入放大以及生物磁检测中的梯度线圈设计，都是简单有效的共模抑制方法。但是输入端的任何不对称，都会使得共模变成差动成分而进入输入电路。于是，在测量电路中应用最广泛的就是陷波器-一种窄带的带阻滤波方法。但是这样经典的滤波技术，不仅滤除掉工频电干扰，同样也滤除掉频带内的所有信息。更严重的是，如果突遇大的信号起伏变动窄带往往会引起自激震荡，即便是自适应滤波技术在此种情况下也会遭遇有限响应时间的问题。另外一大类信号分解方法应用于工频干扰滤除，例如独立源分析（Independent Component Analysis），信号空间投射（signal space projection），信号空间分离法（Signal SpaceSeparation）等，要么对信号分布有严格的要求例如信号与噪声的相对严格的正交分布，要么对测量系统有严格的限制例如SSS仅适用于脑磁测量的系统中，并且对测量线圈的空间位置有很高的精度要求。而且，从实际应用来看，这些方法都有很大的工频干扰残余在信号中。

通常，工频电干扰的频率和幅度，是随着时间不断变化着的。但是，在相对较短的时间内，例如一、两秒钟内，这样的变化非常之小，以至于我们可以将工频干扰看作是恒定的。因此，工频电干扰的基波和任意次谐波，在这样的时间段内都可以看作是一个稳定的正弦波。如果我们在测量的数据中将此正弦波看作一个信号，而将其它测量到的其它信息作为噪声，从而利用傅立叶变换计算出正弦波的三个基本参数频率、幅度和相位，那就能唯一的确定工频电干扰，进而在测量数据中减去该正弦波从而实现工频滤波。显然，这样的滤波技术具有窄带且无自激的理想特性。

有很多公开发表的利用傅立叶变换以及所加的窗函数来准确估算平稳信号中正弦波的基本方法。发明人早先也发展了一种在外加干扰突变情况准确计算正弦波并成功地应用于生理信息检测中的工频电干扰滤除，在该方法中，发明人利用仿真技术研究了干扰宽度、信号长度、采样率、以及信噪比对正确估算正弦波三个基本参数的影响。当其它因素保持不变的条件下，信噪比（注意这里所谓的信号是指要提取的工频干扰正弦波）越高，越有利于正确估算正弦波。如果工频干扰的正弦波估算不准确，就会显然在利用该方法减除工频干扰后仍有干扰残留。尤其是，当频率估算不准确时，相减后就会产生一个对应于这个相差的差频正弦波残留在真正的信号中。使信号产生误差。这对于解释信号和后续的分析，都会导致新的误差。

在多道生理信息测量记录系统中，例如脑电图，脑磁图等，总会有一些通道受到的电干扰小于其它通道。这本来并不是一件坏事，但是，在利于估算正弦波减除工频电干扰时，却容易造成比受到大干扰的通道更大的误差。

发明内容

本发明的目的在于为临床电生理诊断和实验室研究提供一种有效、高保真的多通道微弱生理信息记录系统中工频干扰的高保真滤除方法，该滤波方法对工频电干扰的影响减到最小，同时不影响真正的信号分量，即便对于和工频电干扰同样频率的信号分量也是如此。

为达到上述目的，本发明采用的技术方案是：

1）首先对N个通道逐一进行如下分析：

对每一个通道的数据细分为M个数据单元，每一个数据单元长度取决于工频的变化快慢，在每一个数据单元里，工频电干扰对应的正弦波的三个基本参数即频率、幅度和相位保持不变，对每一个数据单元的数据进行包络检查，查验是否有冲击干扰；

如果没有冲击干扰，利用传统的STEM方法估算正弦波的三个基本参数即频率、幅度和相位；

如果有冲击干扰，利用RAW-STEM方法来提取正弦波的三个参数即频率、幅度和相位；

同时，计算出这一数据单元的标准方差SD，

再对下一个通道进行以上操作，得到所有通道在所有时刻即所有数据单元的工频电干扰的数据库，这个数据库就是正弦波的参数估计值，其中的幅度估计值的准确可靠性是用幅度估计值与标准方差的比值来衡量的；

2）利用提取的幅度与标准方差之比值，作为加权系数，对信噪比最高的数个通道中估算的频率进行加权平均，如公式（1）所示，从而进一步优化参数估算的准确性

$f_0=\frac{1}{{\mathrm{\Σ}}_{n=1}^L\frac{a_n}{{\mathrm{SD}}_n}}\·{\mathrm{\Σ}}_{n=1}^L\frac{a_n}{{\mathrm{SD}}_n}{f}_n---\left(1\right)$

式中，L是所选出的具有最高信噪比的通道数，a_n是利用传统STEM或者RAW-STEM估算出的第n个通道的工频电干扰的幅度，SD_n是第n个通道在当前时刻的标准方差，f_n是当前通道在当前时刻的工频干扰频率估算值，而f₀则是整个测量系统所有通道工频电干扰的当前估算频率；

3）把公式（1）中的f_n代之以趋于一致的各个通道的相位估计值，估算出准确的相位；

4）选择所有通道信噪比最高的那些时刻即数据单元，来求得各个通道在这些时刻的幅度值，利用公式（1）的加权方法，f_n代之以各个通道相应的幅度估计值，来计算各个通道幅度的比例值，这样，在任意时刻，利用该比例关系和该时刻的一个或数个最高信噪比通道提供的幅度值，求出其余较低信噪比通道的在该时刻的幅度。

本发明滤除干净而不损失信号分量，即便信号包含有工频基波和谐波同样的频率成分。经滤除工频干扰后的信号，可以继续进行诸如Coherence的后续分析而不受滤波器本身特性的影响。

本发明的基本原理根据是，在多通道测量记录系统中，所有的通道应该观察到相同的工频电干扰源，这样的干扰源可以是一个或者多个干扰源在空间分布，尽管干扰的幅度通常会因通道而异。换句话说，所有通道受到的工频电干扰，幅度可以不同，但是频率应该是一致的。而其相位，也可能一致，也可能有些通道和另外的通道刚好反相180度，这取决于测量系统输入端的设置。根据多通道测量系统的这个特点，利用信噪比（这里的信号是电干扰正弦波，而噪声是背景记录信号等，下同，除非另作特别说明）较高的通道估算工频电干扰正弦波的参数。其计算出的频率就是所有通道的工频电干扰频率，而计算出来的相位，取决于具体的记录系统中各个通道的空间分布来决定是否反相180度，这可以由本身通道信息估算出的工频电干扰正弦参数来判断。这样一来，所有的通道尤其是那些信噪比较低的通道，都得到了和高信噪比通道一样的三个参数中的两个参数的准确估计。对于三个参数中的最后一个参数幅度，可以利用这样一个关系来优化低信噪比的通道中干扰的幅度的估算。这个关系就是，所有通道中的工频干扰的幅度的相对比例，是保持不变的。这样一来，所有通道中的工频正弦波的基本参数，都可以按照该多道测量系统中信噪比最高的通道来估算，从而提高了低信噪比通道的估算准确性。

附图说明

图1是本发明的流程图；

图2是对一段NetStation脑电记录的脑电图的工频干扰（60Hz）进行滤波处理。（A）为处理前的频谱图，（B）为对干扰频率附近展开以详细对比滤波技术的影响。（C）为用该发明滤除60Hz工频干扰后的频谱图，而（D）为滤除干扰后的频谱在59～61Hz展开。

具体实施方式

参见图1，首先对N个通道逐一分析。而每一个通道的数据，又细分为M个的数据单元。只要这样的数据单元足够的短例如数秒钟，就可以合理的假设，其中的工频电干扰的正弦波是平稳的从而三个基本参数即频率、幅度和相位保持不变。进一步对这个数据单元的数据进行包络检查，查验是否有冲击干扰。如果没有冲击干扰，就利用传统的STEM方法（Grandke T.Interpolationalgorithms for discrete Fourier transforms of weighted signals.IEEETrans Instrum Meas1983;IM-32:350–5.）估算正弦波的三个基本参数，如果有冲击干扰，利用RAW-STEM方法（Tongsheng Zhang and Yoshio Okada,Recursive artifact windowed–single tone extraction method(RAW–STEM)as periodic noise filter for electrophysiological signals withinterfering transients,Journal of Neuroscience Methods155(2006)308–318）来提取正弦波的三个参数。同时，计算出这一单元数据的标准方差SD。由这一段数据提取的正弦波幅度与这个标准方差之比，可以判断所估算的正弦波的三个基本参数的准确性。该比值越高，所估算的参数越准确。当这个流程图所展示的程序运行完结时，就输出了反映所有通道在所有时刻的工频电干扰的数据库，这个数据库就是正弦波的参数估计值和对这个估计值准确性进行量度的幅度与标准方差的比值。

有了这样一个数据库，可以根据信噪比最高的数个通道，例如占总通道数20％的通道的工频干扰正弦波参数，来确定整个多通道系统的工频干扰在某个时刻的频率和相位。而幅度信息，可以选择信噪比最高的时刻所提供的各通道幅度比例，结合本时刻最高信噪比通道的幅度来确定。这件，就使得所有通道在每个时刻的正弦波参数估计达到最优。

事实上在具体应用中，即便是信噪比最高的通道，其信噪比也是有限的。而发明者早先的研究表明，背景噪声和生理信息对工频干扰正弦波的参数的估算，总是在正确值两侧变化。因此，本发明进一步利用上述步骤中提取的幅度与标准方差之比值，作为加权系数，对信噪比最高的数个通道，例如从总通道数20％的通道中估算的频率进行加权平均，如公式（1）所示。从而进一步优化参数估算的准确性。

$f_0=\frac{1}{{\mathrm{\Σ}}_{n=1}^L\frac{a_n}{{\mathrm{SD}}_n}}\·{\mathrm{\Σ}}_{n=1}^L\frac{a_n}{{\mathrm{SD}}_n}{f}_n---\left(1\right)$

式中，L是所选出的具有最高信噪比的通道数，a_n是利用传统STEM或者RAW-STEM估算出的第n个通道的工频电干扰的幅度，SD_n是第n个通道在当前时刻的标准方差，f_n是当前通道在当前时刻的工频干扰频率估算值。而f₀则是整个测量系统所有通道工频电干扰的当前估算频率。

而对前述各个通道的相位估计值，有两种可能，一是所有通道趋于一个值，二是所有通道的相位趋于两个相差180度的两个值。不论是哪种情况，把公式（1）中的f_n代之以趋于一致的各个通道的相位估计值，估算出准确的相位；

各个通道在每个时刻的幅度，是相对复杂的一个物理量。在本发明中对其优化估计的方法也相对复杂。假设干扰源相对于多通道测量系统的空间位置不变-这在大多数情况下是一个合理的假设，那么，一个符合逻辑的推理就是各个通道所观察到的工频电干扰的基波和任意次谐波的幅度，其比例保持不变。因此，我们可以选择所有通道信噪比最高的那些时刻（例如刺激前的记录时段），来求得各个通道在这些时刻的幅度值，进而利用公式（1）的加权方法，f_n代之以各个通道相应的幅度估计值，来计算各个通道幅度的比例值。这样，在任意时刻，利用该比例关系和该时刻的一个或数个最高信噪比通道提供的幅度值，就可以求出其余较低信噪比通道的在该时刻的幅度。

下面是一些对常见多通道电生理系统利用该方法滤除工频电干扰的结果。

图2是对一段NetStation脑电记录的脑电图的工频干扰（60Hz）进行滤波处理。（A）为处理前的频谱图，（B）为对干扰频率附近展开以详细对比滤波技术的影响。（C）为用该发明滤除60Hz工频干扰后的频谱图，而（D）为滤除干扰后的频谱在59～61Hz展开。注意幅度轴在0.1处截断以突显滤波前后的区别。很明显，工频干扰的基波（60Hz）得到了完全的滤除而没有出现经常在陷波器技术中出现的频谱洼陷。

Claims (1)

1.一种多通道微弱生理信息记录系统中工频干扰的高保真滤除方法，其特征在于包括以下步骤：

1）首先对N个通道逐一进行如下分析：

对每一个通道的数据细分为M个数据单元，每一个数据单元长度取决于工频的变化快慢，在每一个数据单元里，工频电干扰对应的正弦波的三个基本参数即频率、幅度和相位保持不变，对每一个数据单元的数据进行包络检查，查验是否有冲击干扰；

如果没有冲击干扰，利用传统的STEM方法估算正弦波的三个基本参数即频率、幅度和相位；

如果有冲击干扰，利用RAW-STEM方法来提取正弦波的三个参数即频率、幅度和相位；

同时，计算出这一数据单元的标准方差SD，

再对下一个通道进行以上操作，得到所有通道在所有时刻即所有数据单元的工频电干扰的数据库，这个数据库就是正弦波的参数估计值，其中的幅度估计值的准确可靠性是用幅度估计值与标准方差的比值来衡量的；

2）利用提取的幅度与标准方差之比值，作为加权系数，对信噪比最高的数个通道中估算的频率进行加权平均，如公式（1）所示，从而进一步优化参数估算的准确性

$f_0=\frac{1}{{\mathrm{\Σ}}_{n=1}^L\frac{a_n}{{\mathrm{SD}}_n}}\·{\mathrm{\Σ}}_{n=1}^L\frac{a_n}{{\mathrm{SD}}_n}{f}_n---\left(1\right)$

式中，L是所选出的具有最高信噪比的通道数，a_n是利用传统STEM或者RAW-STEM估算出的第n个通道的工频电干扰的幅度，SD_n是第n个通道在当前时刻的标准方差，f_n是当前通道在当前时刻的工频干扰频率估算值，而f₀则是整个测量系统所有通道工频电干扰的当前估算频率；

3）把公式（1）中的f_n代之以趋于一致的各个通道的相位估计值，估算出准确的相位；

4）选择所有通道信噪比最高的那些时刻即数据单元，来求得各个通道在这些时刻的幅度值，利用公式（1）的加权方法，f_n代之以各个通道相应的幅度估计值，来计算各个通道幅度的比例值，这样，在任意时刻，利用该比例关系和该时刻的一个或数个最高信噪比通道提供的幅度值，求出其余较低信噪比通道的在该时刻的幅度。

Images