CN103761424B - 基于二代小波和独立分量分析肌电信号降噪与去混迭方法 - Google Patents
基于二代小波和独立分量分析肌电信号降噪与去混迭方法 Download PDFInfo
- Publication number
- CN103761424B CN103761424B CN201310755427.8A CN201310755427A CN103761424B CN 103761424 B CN103761424 B CN 103761424B CN 201310755427 A CN201310755427 A CN 201310755427A CN 103761424 B CN103761424 B CN 103761424B
- Authority
- CN
- China
- Prior art keywords
- signal
- noise reduction
- aliasing
- electromyographic signal
- wavelet
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Landscapes
- Measurement And Recording Of Electrical Phenomena And Electrical Characteristics Of The Living Body (AREA)
Abstract
本发明涉及一种基于二代小波变换和独立分量分析相结合的肌电信号降噪与去混迭方法。目前有的信号去混迭方法运用于肌电信号时显现出一些不足。本发明首先采用二代小波降噪算法先滤除信号中的噪声。二代小波降噪算法运用于肌电信号,使降噪后的信号变得平滑,波形当中不必要的震荡得到抑制,肌电信号特征更加明显。再对降噪后的含混信号进行ICA分离,可快速有效的去除信号中的混迭成分。利用本发明对肌电信号做以上预处理,可大大去除信号中的干扰,为后续的肌电信号的特征提取,动作识别等研究提供了便捷。
Description
技术领域
本发明属于多通道信号去混迭领域,涉及一种应用于肌电信号的消噪与去混迭方法。
背景技术
肌电信号包含了丰富的肌肉运动的信息,它能够反映肢体的运动模式。目前,肌电信号已经被广泛应用于临床诊断、康复工程、运动医学等领域。表面肌电信号(sEMG)是一种在皮肤表面募集到的,伴随着肢体运动导致肌肉收缩而产生的生物电信号。基于肌电信号的人机接口通过对表面肌电信号的处理和模式识别可识别出人体不同的动作模式,控制外部环境设备工作,如基于肌电信号的仿生假手的控制,基于肌电信号的跌倒辨识。
基于肌电信号的人机接口研究中,要通过肌电信号对多个动作实现模式识别,往往需要用多个电极同时进行多路肌电信号的采集。由于人体本身存在天然电流,各神经元之间又相互连接,以及采样电极相邻通道之间存在耦合串扰等因素,使得多通道肌电信号采集到的信号存在相邻肌电信号的混迭成分。为了得到纯粹的肌电信号,进行进一步研究分析,我们需要一种有效的方法对肌电信号进行去混迭处理。
罗志增等针对脑电信号中的噪声及夹杂的眼电、心电等伪迹,采用小波分解重构消噪和基于熵估计的RADICAL算法进行消除。对3路脑电观测信号进行小波消噪和白化处理后,通过雅可比旋转矩阵分别对其中的两两组合信号用RADICAL算法进行分离,得出最优分离矩阵完成盲源分离。张洪渊等利用随机变量概率密度函数非参数估计的核函数法,对混合信号的概率密度函数及其导数进行估计,可以直接应用于所有以非线性函数代替评价函数的信号源盲分离算法。以上方法都是不错的信号盲分离算法,但是基于sEMG是一种非线性非平稳信号,以及人体活动时各肌肉的协同工作的特点,使得以上方法运用于sEMG时显现出一些不足。
sEMG的有用能量主要分布在10-500Hz之间,非常微弱(其幅值仅为μV级),极易受到噪声的干扰,在对sEMG作去混迭处理前必须滤除其中的噪声。以改善肌电信号的盲源分离效果。
传统的小波降噪方法对信号进行多尺度分解,采用一个小波基函数去逼近不同尺度上的信号,在小波基函数与逼近信号的差异处会产生细节信息,这些细节信号在阈值处理时会被当作噪声被滤除,使降噪后的信号丢失部分有用信息。这样,使得常规的小波降噪方法无法准确的刻画非平稳信号。
综上所述,已有的信号消噪与去混迭方法在运用于肌电信号时都存在不足之处。本发明提出了一种基于二代小波变换和独立分量分析(ICA)相结合的肌电信号降噪与去混迭方法。
发明内容
本发明的目的就是针对现有技术的不足,提出一种基于二代小波变换和独立分量分析(ICA)相结合的肌电信号降噪与去混迭方法。第二代小波变换继承了传统小波的多分辨率的特性;不依赖傅立叶变换;小波变换后的系数是整数;基于多项式内插的思想,所有的运算都在空间域进行,从而摆脱了对频域的依赖。利用二代小波进行信号降噪能滤除高频噪声并保留有用信息,还能较好地保留sEMG的边缘特性。
独立分量分析(ICA)是根据从传感器监测到的观测信号来估计未知的源信号和混合矩阵的。它从信号高阶统计特性出发,主要目的按是为了提取信号的独立成分,与已有的主分量分析(PCA)方法相比,它不仅可去除分量间的低阶相关性,而且还能够去除分量间的高阶相关性,使之相互独立。因此,ICA是一种有效的盲源分离方法,相对其他方法能更全面地体现数据间的本质结构。
所以,二代小波变换和独立分量分析(ICA)相结合的方法可以有效的降低肌电信号中的噪声并去除相邻通道信号的混迭成分。
为了实现以上目的,本发明方法主要包括以下步骤:
步骤一,对仪器采集到的含噪声和线性混迭的M路观测信号,X=[x1,x2,…,xM]T进行二代小波分解。
对每一路信号x[n]进行二代小波分解的具体步骤如下:
①分裂:将信号序列x[n]分成两个互不相交的子集,通常按奇样本xodd[n]、偶样本xeven[n]分成两个等长的序列。
②预测:采用一个与数据结构无关的预测算子P,保持偶样本不变,由偶样本估计奇样本xodd[n]=P(xeven[n])。如果信号具有局部相关结构,偶样本集合和奇样本集合一定是高度相关的,可以以预测算子P用一个子集估计另一个子集。原值与预测值之差作为小波系数d[n],表示为d[n]=xodd[n]-P(xeven[n])。
③修正:为了保持存在于原始信号中的低频尺度特性,引入修正算子U,以得到下一尺度上的尺度系数:c[n]=xeven[n]-U(d[n])。
上面这三步构成一个提升过程,对输出c[n]重复以上提升过程,可构成一个完整的离散小波变换。得到表示原始数据的低频成分的尺度系数cj[n]和表示原始数据的高频成分的小波系数dj[n]。j为二代小波分解的层级,最高层数记为L。
步骤二,对步骤一所得到的小波系数dj[n]进行阈值处理得到消噪后的高频系数d'j[n],阈值函数如下:
式中α∈[0,1]为一可变参数。当α=0时,阈值函数等效于原始信号。阈值σ=1/0.6745med(|x|),只要选取合适的参数α就可以得到平滑的降噪信号,抑制波形当中不必要的振荡。
步骤三,将步骤二所得的d'j[n]与步骤一所得的cj[n]进行小波重构得到消噪后的信号x'[n]。M路信号都经过降噪处理后得到消噪后的肌电信号X'=[x1',x2',…,xM']T。
步骤四,对步骤三得到的X'进行去均值处理:
X”=X'-E(X')
其中,去均值后的观测数据X”的协方差矩阵可以分解为Cx=E[X”X”T]=UλUT,其中U为Cx的特征向量矩阵,λ为Cx的特征值矩阵,λ=Diag(λ1,λ2,...,λN)。
对X”进行白化处理:
Z=BX”
其中,白化矩阵B=λ-1/2UT。去均值和白化过程使信号成为零均值和具有单位方差且各分量互不相关的矢量Z=[z1,z2,…,zM]T。
步骤五,解混矩阵W的第i列记为解混向量wi。将步骤四得到的Z带入FastICA算法的迭代式得:
wi(k+1)=E{Z(wi(k)T)3Z}-3γE{ZZT}wi(k)
k表示迭代次数。随机产生wi(0),且||wi(0)||2=1。γ为一大于零的可变参数。
步骤六,为保证每次提取出来的都是之前从未提取过的独立分量,在计算wi且i≥2时,需要添加正交化步骤,并把步骤五中得到的wi(k+1)减去已提取的独立分量,即按下式处理:
步骤七,为确保分离出的独立分量具有单位能量,需要对步骤六得到的wi(k+1)进行归一化处理:
步骤八,如果|wi(k+1)Twi(k)|收敛于1,则算法收敛,得到一个解混向量wi,i=i+1,否则返回步骤五,直至算法收敛。
步骤九,重复步骤五,六,七,八得到解混矩阵W=[w1,w2,...,wk]。可求出解混信号:
Y=WZ
即Y=[y1,y2,...,yk]为最终降噪去混迭的信号。
本发明的有益效果:在基于肌电信号的人体动作模式识别的研究中,往往需要用多个电极同时进行多路信号的采集。由于人体本身存在天然电流,各神经元之间又相互连接,以及采样电极相邻通道之间存在耦合串扰等因素,使得多通道肌电信号采集到的信号存在相邻肌电信号的混迭成分。利用本发明对含混信号进行ICA分离,可快速有效的去除信号中的混迭成分。同时,肌电信号本身具有非线性、非平稳,非常微弱极易受到噪声干扰的特点,噪声的加入会影响模式识别的结果。本发明在对sEMG作去混迭处理前采用二代小波降噪算法先滤除信号中的噪声。二代小波降噪算法运用于肌电信号,使降噪后的信号变得平滑,波形当中不必要的震荡得到抑制,肌电信号特征更加明显。利用本发明对肌电信号做以上预处理,可大大去除信号中的干扰,以提高模式识别的准确度。
附图说明
图1为本发明的实施流程图;
图2为本发明实施例的3路原始含有噪声和混迭的肌电信号图;
图3为本发明实施例经过二代小波降噪的结果图;
图4为本发明实施例的降噪去混迭结果图。
具体实施方式
下面结合附图对本发明的实施例作详细说明:本实施例在以本发明技术方案为前提下进行实施,给出了详细的实施方式和具体的操作过程,但本发明的保护范围不限于下述的实施例。
如图1所示,本实施例包括如下步骤:
步骤一,对仪器采集到的含噪线性混迭的M路观测信号X=[x1,x2,…,xM]T进行二代小波分解,在本实施例中M=3,肌电信号如图2所示。
对每一路信号x[n]进行二代小波分解的具体步骤如下:
①分裂:将信号序列x[n]分成两个互不相交的子集,通常按奇样本xodd[n]、偶样本xeven[n]分成两个等长的序列。
②预测:采用一个与数据结构无关的预测算子P,保持偶样本不变,由偶样本估计奇样本xodd[n]=P(xeven[n])。如果信号具有局部相关结构,偶样本集合和奇样本集合一定是高度相关的,可以以预测算子P用一个子集估计另一个子集。原值与预测值之差作为小波系数d[n],表示为d[n]=xodd[n]-P(xeven[n])。
③修正:为了保持存在于原始信号中的低频尺度特性,引入修正算子U,以得到下一尺度上的尺度系数:c[n]=xeven[n]-U(d[n])。
上面这三步构成一个提升过程,对输出c[n]重复以上提升过程,可构成一个完整的离散小波变换。得到表示原始数据的低频成分的尺度系数cj[n]和表示原始数据的高频成分的小波系数dj[n]。j为二代小波分解的层级,最高层数记为L,在本实施例中L=4。
步骤二,对步骤一所得到的小波系数dj[n]进行阈值处理得到消噪后的高频系数d'j[n],阈值函数如下:
式中α∈[0,1]为一可变参数。当α=0时,阈值函数等效于原始信号。阈值σ=1/0.6745med(|x|),只要选取合适的参数α就可以得到平滑的降噪信号,抑制波形当中不必要的振荡。在本实施例中α=0.5,N=2020。
步骤三,将步骤二所得的d'j[n]与步骤一所得的cj[n]进行小波重构得到消噪后的信号x'[n]。M路信号都经过降噪处理后得到X'=[x1',x2',…,xM']T,如图3所示。
步骤四,对步骤三得到的X'进行去均值处理:
X”=X'-E(X')
其中,去均值后的观测数据X”的协方差矩阵可以分解为Cx=E[X”X”T]=UλUT,其中U为Cx的特征向量矩阵,λ为Cx的特征值矩阵,λ=Diag(λ1,λ2,...,λN)。
对X”进行白化处理:
z=BX”
其中,白化矩阵B=λ-1/2UT。去均值和白化过程使信号成为零均值和具有单位方差且各分量互不相关的矢量Z=[z1,z2,…,zM]T。
步骤五,解混矩阵W的第i列记为解混向量wi。将步骤四得到的Z带入FastICA算法的迭代式得:
wi(k+1)=E{Z(wi(k)T)3Z}-3γE{ZZT}wi(k)
k表示迭代次数。随机产生wi(0),且||wi(0)||2=1。γ为一大于零的可变参数。在本实施例中γ=0.8。
步骤六,为保证每次提取出来的都是之前从未提取过的独立分量,在计算wi且i≥2时,需要添加正交化步骤,并把步骤五中得到的wi(k+1)减去已提取的独立分量,即按下式处理:
步骤七,为确保分离出的独立分量具有单位能量,需要对步骤六得到的wi(k+1)进行归一化处理:
步骤八,如果|wi(k+1)Twi(k)|收敛于1,则算法收敛,得到一个解混向量wi,i=i+1,否则返回步骤五,直至算法收敛。
步骤九,重复步骤五,六,七,八得到解混矩阵W=[w1,w2,...,wk]。可求出解混信号:
Y=WZ
即Y=[y1,y2,...,yk]为最终降噪去混迭的信号,如图4所示。
Claims (1)
1.基于二代小波和独立分量分析肌电信号降噪与去混迭方法,其特征在于该方法包括如下步骤:
步骤一,对仪器采集到的含噪声和线性混迭的M路观测信号,X=[x1,x2,…,xM]T进行二代小波分解;
对每一路信号x[n]进行二代小波分解的具体步骤如下:
①分裂:将信号序列x[n]分成两个互不相交的子集,按奇样本xodd[n]、偶样本xeven[n]分成两个等长的序列;
②预测:采用一个与数据结构无关的预测算子P,保持偶样本不变,由偶样本估计奇样本xodd[n]=P(xeven[n]);如果信号具有局部相关结构,偶样本集合和奇样本集合一定是高度相关的,以预测算子P用一个子集估计另一个子集;原值与预测值之差作为小波系数d[n],表示为d[n]=xodd[n]-P(xeven[n]);
③修正:为了保持存在于原始信号中的低频尺度特性,引入修正算子U,以得到下一尺度上的尺度系数:c[n]=xeven[n]-U(d[n]);
上面这三步构成一个提升过程,对输出c[n]重复以上提升过程,构成一个完整的离散小波变换;得到表示原始数据的低频成分的尺度系数cj[n]和表示原始数据的高频成分的小波系数dj[n];j为二代小波分解的层级,最高层数记为L;
步骤二,对步骤一所得到的小波系数dj[n]进行阈值处理得到消噪后的高频系数d'j[n],阈值函数如下:
式中α∈[0,1]为一可变参数;当α=0时,阈值函数等效于原始信号;阈值
步骤三,将步骤二所得的d'j[n]与步骤一所得的cj[n]进行小波重构得到消噪后的信号x'[n];M路信号都经过降噪处理后得到消噪后的肌电信号X'=[x1',x2',…,xM']T;
步骤四,对步骤三得到的X'进行去均值处理:
X”=X'-E(X')
其中,去均值后的观测数据X”的协方差矩阵分解为Cx=E[X”X”T]=UλUT,其中U为Cx的特征向量矩阵,λ为Cx的特征值矩阵,λ=Diag(λ1,λ2,...,λN);
对X”进行白化处理:
Z=BX”
其中,白化矩阵B=λ-1/2UT;去均值和白化过程使信号成为零均值和具有单位方差且各分量互不相关的矢量Z=[z1,z2,…,zM]T;
步骤五,解混矩阵W的第i列记为解混向量wi;将步骤四得到的Z带入FastICA算法的迭代式得:
wi(k+1)=E{Z(wi(k)T)3Z}-3γE{ZZT}wi(k)
k表示迭代次数;随机产生wi(0),且||wi(0)||2=1;γ为一大于零的可变参数;
步骤六,为保证每次提取出来的都是之前从未提取过的独立分量,在计算wi且i≥2时,添加正交化步骤,并把步骤五中得到的wi(k+1)减去已提取的独立分量,即按下式处理:
步骤七,对步骤六得到的wi(k+1)进行归一化处理:
步骤八,如果|wi(k+1)Twi(k)|收敛于1,则算法收敛,得到一个解混向量wi,i=i+1,否则返回步骤五,直至算法收敛;
步骤九,重复步骤五,六,七,八得到解混矩阵W=[w1,w2,...,wk];可求出解混信号:
Y=WZ
即Y=[y1,y2,...,yk]为最终降噪去混迭的信号。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201310755427.8A CN103761424B (zh) | 2013-12-31 | 2013-12-31 | 基于二代小波和独立分量分析肌电信号降噪与去混迭方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201310755427.8A CN103761424B (zh) | 2013-12-31 | 2013-12-31 | 基于二代小波和独立分量分析肌电信号降噪与去混迭方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN103761424A CN103761424A (zh) | 2014-04-30 |
CN103761424B true CN103761424B (zh) | 2016-09-14 |
Family
ID=50528660
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201310755427.8A Active CN103761424B (zh) | 2013-12-31 | 2013-12-31 | 基于二代小波和独立分量分析肌电信号降噪与去混迭方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN103761424B (zh) |
Families Citing this family (10)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104367317B (zh) * | 2014-10-15 | 2017-04-12 | 北京理工大学 | 多尺度特征融合的心电图心电信号分类方法 |
CN105741305A (zh) * | 2016-03-02 | 2016-07-06 | 深圳竹信科技有限公司 | 一种基于平稳小波变换滤除肌电干扰的方法和系统 |
CN106096579A (zh) * | 2016-06-22 | 2016-11-09 | 天津理工大学 | 一种心电信号预处理的方法 |
CN106264521A (zh) * | 2016-09-22 | 2017-01-04 | 小菜儿成都信息科技有限公司 | 多通道脑电信号中下颌干扰的自动去除方法 |
CN108319889A (zh) * | 2017-01-18 | 2018-07-24 | 上海诺师信息科技有限公司 | 利用独立分量分析进行emg信号特征提取的方法 |
CN107468260A (zh) * | 2017-10-12 | 2017-12-15 | 公安部南昌警犬基地 | 一种判断动物心理状态的脑电分析装置及分析方法 |
CN108836321A (zh) * | 2018-05-03 | 2018-11-20 | 江苏师范大学 | 一种基于自适应噪声消除系统的脑电信号预处理方法 |
CN108918927B (zh) * | 2018-08-24 | 2020-11-10 | 广东石油化工学院 | 一种能量分解中功率信号滤波方法及系统 |
CN109567798A (zh) * | 2018-12-26 | 2019-04-05 | 杭州电子科技大学 | 基于肌电小波相干性和支持向量机的日常行为识别方法 |
CN110132866A (zh) * | 2019-05-27 | 2019-08-16 | 东北大学 | 一种土壤高光谱解混方法及系统 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102697495A (zh) * | 2012-05-23 | 2012-10-03 | 杭州电子科技大学 | 基于总体平均经验模式分解的二代小波肌电信号消噪方法 |
CN102832907A (zh) * | 2012-08-29 | 2012-12-19 | 杭州电子科技大学 | 一种消除触觉传感器工频噪声的方法 |
-
2013
- 2013-12-31 CN CN201310755427.8A patent/CN103761424B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102697495A (zh) * | 2012-05-23 | 2012-10-03 | 杭州电子科技大学 | 基于总体平均经验模式分解的二代小波肌电信号消噪方法 |
CN102832907A (zh) * | 2012-08-29 | 2012-12-19 | 杭州电子科技大学 | 一种消除触觉传感器工频噪声的方法 |
Also Published As
Publication number | Publication date |
---|---|
CN103761424A (zh) | 2014-04-30 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN103761424B (zh) | 基于二代小波和独立分量分析肌电信号降噪与去混迭方法 | |
Chen et al. | Removal of muscle artifacts from the EEG: A review and recommendations | |
CN102697493B (zh) | 一种快速的脑电信号中眼电伪迹自动识别和去除的方法 | |
CN103961092B (zh) | 基于自适应阈值处理的脑电信号去噪方法 | |
CN102697495B (zh) | 基于总体平均经验模式分解的二代小波肌电信号消噪方法 | |
CN101221554A (zh) | 基于小波变换和bp神经网络的脑电特征提取方法 | |
CN109498370B (zh) | 基于肌电小波关联维的下肢关节角度预测方法 | |
CN105342605A (zh) | 一种去除脑电信号中肌电伪迹的方法 | |
CN104688220A (zh) | 一种去除脑电信号中眼电伪迹的方法 | |
CN109359619A (zh) | 一种基于卷积盲源分离的高密度表面肌电信号分解方法 | |
Hettiarachchi et al. | Motor imagery data classification for BCI application using wavelet packet feature extraction | |
Abdullah et al. | Surface EMG signal classification by using WPD and ensemble tree classifiers | |
CN104905786A (zh) | 一种心电伪迹在线去除算法 | |
Sheoran et al. | Methods of denoising of electroencephalogram signal: A review | |
CN108836301A (zh) | 一种基于奇异谱分析和稀疏表示的单次诱发电位提取方法 | |
Ferdowsi et al. | Removing ballistocardiogram artifact from EEG using short-and long-term linear predictor | |
CN113842115A (zh) | 一种改进的eeg信号特征提取方法 | |
CN116584960A (zh) | 一种膈肌肌电信号降噪方法 | |
CN114757236A (zh) | 基于tqwt与svmd的脑电信号去噪优化方法及系统 | |
Ma et al. | Classification of motor imagery EEG signals based on wavelet transform and sample entropy | |
CN1744073A (zh) | 利用小波神经网络提取想象动作电位的方法 | |
CN104000587A (zh) | 一种基于边缘小波特征的脑电波(eeg)信号识别系统 | |
Bustomi et al. | Analyzing power spectral of electroencephalogram (EEG) signal to identify motoric arm movement using EMOTIV EPOC+ | |
CN114358090A (zh) | 一种基于psd和csp的运动想象脑电信号分类方法 | |
Kaur et al. | Wavelet based machine learning technique to classify the different shoulder movement of upper limb amputee |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
C14 | Grant of patent or utility model | ||
GR01 | Patent grant | ||
EE01 | Entry into force of recordation of patent licensing contract | ||
EE01 | Entry into force of recordation of patent licensing contract |
Application publication date: 20140430 Assignee: HANGZHOU DUKANG TECHNOLOGY CO.,LTD. Assignor: HANGZHOU DIANZI University Contract record no.: X2022330000025 Denomination of invention: EMG signal denoising and de aliasing method based on second generation wavelet and independent component analysis Granted publication date: 20160914 License type: Common License Record date: 20220128 |