CN102930284B

CN102930284B - 基于经验模态分解和分形的表面肌电信号模式识别方法

Info

Publication number: CN102930284B
Application number: CN201210337353.1A
Authority: CN
Inventors: 张启忠; 席旭刚; 罗志增; 佘青山; 高云园
Original assignee: Hangzhou Dianzi University
Current assignee: Hangzhou Dianzi University
Priority date: 2012-09-13
Filing date: 2012-09-13
Publication date: 2015-06-17
Anticipated expiration: 2032-09-13
Also published as: CN102930284A

Abstract

本发明涉及一种基于经验模态分解和分形的表面肌电信号模式识别方法。现有的方法大都采用单分形理论，对肌电信号仅做整体的奇异性评价，而没有研究信号的局部奇异性特征。本发明首先从相关肌肉组上采集相应的表面肌电信号，然后用经验模态分解的方法提取肌电信号的多层内在模态函数，利用多重分形分析的方法提取各层内在模态函数上的广义维数谱。最后，以各层模态函数上的广义维数谱作为模式识别的特征向量，以支持向量机为分类器实现多运动模式的分类与识别。本发明提出利用多重分形分析的方法提取各层内在模态函数上的广义维数谱作为表面肌电信号的特征，具有较好的鲁棒性，能从信噪比相对较低的肌电信号中计算出稳定的特征数据。

Description

基于经验模态分解和分形的表面肌电信号模式识别方法

技术领域

本发明属于模式识别领域，涉及一种肌电信号模式识别方法，特别涉及一种基于肌电信号的手部多运动模式识别方法。

背景技术

肌电信号（Electromyography，EMG）是一种伴随肌肉活动的生物电信号，是众多肌纤维中运动单元动作电位的叠加，蕴涵了肌肉活动的各种信息。而表面肌电信号（Surface electromyography ，sEMG）则是浅层肌肉EMG和神经干上电活动在皮肤表面的综合效应。表面肌电信号在生理基础上具有高度非线性特征，是一个非线性的动力学系统。因而，对表面肌电信号的研究采用非线性分析方法，如混沌、分形理论是一个值得重视的方向。

目前，已有一些科研机构的学者以混沌、分形理论作为数学工具，对肌电信号进行各种目的的研究。国外，克罗地亚萨格勒布大学的Cifrek, Mario用表面肌电信号的分维数对肌肉疲劳程度进行参数估计。澳大利亚墨尔本大学的Naik,.Ganesh R、泰国宋卡大学的phinyomark Angkoon等用分形的方法研究了上臂弱肌电信号的多模式分类问题，并把研究成果应用到人机接口技术中。国内，清华大学的王人成，上海交通大学的王立中等、浙江大学的刘加海，王健等在混沌、分形理论应用于肌电信号信息处理方面做了较多的工作。然而，综合这些学者的研究成果发现，现有的对肌电信号的研究，大都采用单分形理论，对所研究的肌电信号仅做整体的奇异性评价，如在模式识别时，提取整体信号的Lyapunov指数、分形维等，而没有研究信号的局部奇异性特征。

发明内容

表面肌电信号具有非线性非平稳特性，为了有效地提取表面肌电信号的特征和提高手部多运动模式识别率，本发明提出了一种经验模态分解（Empirical Mode Decomposition,EMD）与多重分形分析相结合的特征提取方法。首先从相关肌肉组上采集相应的表面肌电信号，然后用经验模态分解的方法提取肌电信号的多层内在模态函数(intrinsic mode function, IMF)，利用多重分形分析的方法提取各层内在模态函数上的广义维数谱，从而实现了对表面肌电信号在不同层次上的奇异性测度，提高了对信号几何特征和局部尺度行为刻画的精细程度。最后，以各层模态函数上的广义维数谱作为模式识别的特征向量，以支持向量机为分类器实现多运动模式的分类与识别。

为了实现以上目的，本发明方法主要包括以下步骤：

步骤(1).获取人体上肢肌电信号样本数据，具体是：首先通过肌电信号采集仪拾取人体上肢肌电信号，再采用空域相关滤波方法对含有干扰噪声的肌电信号进行消噪。

所述的空域相关滤波方法，具体如下：

用空域相关滤波对采集的SEMG信号进行消噪预处理，由于空域相关滤波中涉及噪声能量阈值的设定，且噪声能量阈值的设定没有通用的算法，因此给出一种针对SEMG信号的噪声能量阈值设定算法。算法的具体实现如下：

由于手部未动作时采集的SEMG信号对应为信号的噪声，若对该含有噪声的SEMG信号进行多尺度小波分解，各尺度的高频系数的首尾部分对应手部无动作时的SEMG信号，因此可以用手部无动作时的SEMG信号的高频系数来估计各层的噪声能量门限。以伸腕动作SEMG信号为例，取伸腕动作SEMG信号小波变换高频系数的前个点，取到的这些点对应于手部未动作时SEMG信号的高频系数，用这些点的方差来估计SEMG信号的噪声能量阈值，此时滤波算法中用下式控制迭代过程：

; （1）

其中为比例系数，是根据经验和实验确定的；为迭代后的高频系数的方差；为估计的噪声能量阈值。当上式成立时迭代继续，否则停止。

空域相关法利用真实信号的小波系数在各尺度上的相关性，使得小尺度上的真实信号的小波系数得以突显，从而实现了真实信号的小波系数与噪声小波系数的分离，取得了比较好的滤波效果。

步骤(2).将步骤（1）获取的肌电信号进行经验模态分解，得到表示肌电信号细节构成的多层窄带信号即内在模态函数。

经验模态分解（EMD）方法是Huang等人于1998年提出的一种时间序列处理方法。该方法根据信号本身的局部特征将信号分解成一组内在模态函数(IMF)。其实质是将非线性非平稳信号分解成一组单分量近似窄带信号。每个IMF满足下面条件：

1）整个IMF中极点数（包括极大值和极小值）和零点数相等或至多相差1。

2）信号关于时间轴对称，即信号由局部极大值确定的包络线和局部极小值确定的包络线均值为零。

EMD的分解步骤如下：

1）求信号的局部极大、极小值。通过插值函数获得信号的上下包络线。并对上下包络线求平均，记为。

2）设原信号为，用原信号减去，记为：。

3）判断是不是符合基本模态分量的两个要求，是则将赋值给，不是则把看作，重复前二步先得到，再得到。判断是否满足IMF的要求，不满足则不断循环直至次循环后得到的满足IMF的要求。记。可以看成是信号的一个基本模态。

4）记。

将作为信号，重复1)到4)步求出，由此循环可继而求出第二阶IMF、……、第N-1节IMF，第阶剩余信号可表示为，当满足某个控制条件时，筛选结束。原始信号表示为：

（2）

步骤(3). 对对步骤(2)所获得的内在模态函数进行多重分形分析，提取其广义维数谱，作为肌电信号多模式识别的特征向量。

所述的多重分形分析，提取其广义维数谱，具体如下：

设肌电信号时间序列为，嵌入维数，时间延迟为，则重构相空间为：

（3）

式中，为相空间中的第个相点，相空间中的点数为。构造其相空间后，以维相空间集中的一点作为参考点，计算另外点与它的距离，则可统计出落于以点为中心，以小标量为半径的体积元中的点的个数，从而得到q阶关联积分

(4)

其中为Heaviside阶跃函数。

计算广义维数谱算法步骤如下：

1）把一维时间序列信号转换到相空间，求出相点间的最大伸展距离。

2）依据的大小确定计算关联积分时增幅，确定需计算的小标量初始个数,,为单调上升的。根据公式(4)求得,计算数据点对，为表达方便，令为。

3）根据下式计算数据组的斜率：

(5)

4) 增加的值，根据公式(5)重新计算值，若连续多次，k的值变化极小，则停止计算，以最后几次所得值的平均，作为肌电信号的的广义维数谱。

步骤(4). 以步骤（3）获取的特征向量输入支持向量机进行分类识别，获得识别结果。

支持向量机是基于统计学习理论和结构风险最小化原则，其基本思想是把输入空间的样本通过非线性变换映射到高维特征空间，然后在特征空间中求取把样本线性分开的最优分类面。算法使用分类间隔控制线性学习机器的容量,从而使结构风险最小，也使其在有限样本下具有了较强的泛化能力。根据结构风险最小化原则，为了最小化期望风险的上界，SVM通过最优超平面的构造，在固定学习机经验风险的条件下最小化VC (Vapnik和Chervonnenkis)置信范围。这里，最优超平面的构造问题实质上是约束条件下求解一个二次规划问题，以得到一个最优分类函数：

（6）

这里是阈值，是第个训练样本，是测试样本，为总的训练样本数，是拉格朗日系数，为满足Mercer 条件的核函数。不同的核函数变换到不同的特征空间,使用核函数避免了在高维的特征空间中直接计算，一般应用较多的核函数有下面3种：多项式核函数、径向基核函数和神经网络核函数。

SVM的思想是针对两类问题的分类，对于多类问题，必须重新构造SVM分类器来求解，目前主要有两种方法。一种是以Weston在1998年提出的多类算法为代表，这个算法是在经典SVM理论的基础上，重新构造多值分类模型实现多值分类。该算法选择的目标函数十分复杂，实现困难，计算复杂度也非常高，因而较少使用。另一种构造方法是通过组合多个二值子分类器实现对多值分类器的构造，该方法有两种分支算法，即“一对多”和“一对一”算法。“一对多”算法由Vapnik提出，它的基本思想是对于类问题构造个两类分类器，第个SVM用第类中的训练样本作为正的训练样本，而将其它的样本作为负的训练样本，最后的输出是两类分类器输出为最大的那一类。“一对一”算法是在类样本中构造所有可能的两类分类器，每次仅仅在类中的两类样本上训练，结果共构造个分类器，组合这些两类分类器并使用投票法，得票最多的类为新点所属的类。本发明选用“一对一”算法的基础上，运用聚类分析中的类间距离和二叉树的方法来构造多类分类器。

本发明与已有的诸多手部肌电信号动作识别方法相比，具有如下特点：

本发明提出利用多重分形分析的方法提取各层内在模态函数上的广义维数谱作为表面肌电信号的特征，具有较好的鲁棒性，能从信噪比相对较低的肌电信号中计算出稳定的特征数据。

鉴于单重分形分析仅能从整体上标志信号的构造规律，而不能揭示信号局部或细节上的构造规律及生成过程，本发明用经验模态分解的方法得到代表肌电信号细节构成的近似单分量窄带信号-多层内在模态函数，然后在模态函数上进行多重分形分析，提取其广义维数谱，作为肌电信号动作多模式识别的特征向量，该方法克服了单重分形分析只能反映信号的整体特性缺乏对信号局部奇异性描述的缺点。

附图说明

图1为本发明的实施流程图；

图2为本发明实施例展拳动作伸肌的表面肌电信号EMD分解；

图3为本发明实施例展拳动作屈肌的肌电信号EMD分解。

具体实施方式

下面结合附图对本发明的实施例作详细说明：本实施例在以本发明技术方案为前提下进行实施，给出了详细的实施方式和具体的操作过程。

如图1所示，本实施例包括如下步骤：

步骤一，获取人体上肢肌电信号样本数据，具体是：首先通过肌电信号采集仪拾取人体上肢肌电信号，再采用空域相关滤波方法对含有干扰噪声的肌电信号进行消噪。

（1）采集人体上肢的肌电信号。受试者分别进行握拳、展拳、屈腕、伸腕4个动作各60组，选用上肢尺侧腕伸肌和尺侧腕屈肌作为表面肌电信号来源。实验时先用酒精分别在受试者的尺侧腕伸肌和尺侧腕屈肌上擦涂去污，以增强信号拾取能力，采用MyoTrace 400肌电信号采集仪来拾取尺侧腕伸肌和尺侧腕屈肌对应的表面肌电信号。

（2）采用空域相关滤波方法进行消噪处理。对采集的原始肌电信号进行5层小波分解，基小波选用双正交样条小波bior1.5。在空域相关滤波消噪中需要设定某一噪声阈值，取伸腕肌电信号的小波变换高频系数的前（本实验中取200）个点，取到的这些点对应于手部未动作时肌电信号的高频系数，用这些点的方差来估计肌电信号的噪声能量阈值。

步骤二，将步骤一获取的表面肌电信号进行经验模态分解，得到表示肌电信号细节构成的多层窄带信号即内在模态函数。本实施例分解到第6层，图2和图3分别是展拳动作尺侧腕伸肌及尺侧腕屈肌表面肌电信号分解得到6层IMFs。

由于各种内、外源干扰的影响，并不能确定EMD分解所得各模态函数都是有效的肌电信号，因此要对各IMFs做有效性判别。本实施例采用简单有效的功率谱估计的方法，先求得各IMF的功率谱密度，然后统计频率在15～300Hz区间信号能量的百分比，若其大于阈值80%，则认为该层IMF为有效的。由此，可得到某一动作屈肌及伸肌上IMFs有效层数。以上述展拳动作为例，伸肌肌电信号IMFs前4层，屈肌的前3层都是有效层。在不影响模式识别效果的情况下，本实施例对所有动作取IMFs的前3层作为特征向量提取的信号源。

步骤三，对步骤二所获得的内在模态函数进行多重分形分析，提取其广义维数谱，作为肌电信号多模式识别的特征向量。具体步骤如下：

1)相空间的重构：根据公式（3）重构相空间，重构过程中需确定2个非常重要的参数，延迟时间和嵌入维数。它们的选择直接关系到相空间重构质量。本实施例采用Kim H S（Kim H S, Eykholt R, Salas J D. Nonlinear dynamics, delay times, and embedding windows [J]. Physica D, 1999,127(1):48-60.）提出的基于关联积分的C-C算法来同时估计延迟时间和嵌入维数。最后，所得结果为,。

2)把一维时间序列信号转换到相空间后，先求得相点的中心点，以及相点间的最大伸展距离，按等差数列方式设定,。令,求得前80个小标量点的关联积分，得到数据组，用公式（5）求取斜率，增加的值，重新计算，若连续5次计算其数值差不超过0.01，5次运算数值差的和不超过0.1，则停止计算，以最后5次运算的平均作为信号的阶广义维数谱。实验发现，当达到87~95之间，数值最为稳定。本实施选择计算的广义维数谱。每类动作做60次，测定其肌电信号，分别计算其EMD分解后3层内在模态函数上的广义维数谱，从而实现对肌电信号在不同层次上的多重分形分析，表1为动作样本点的统计数据。

表1 四种动作的广义维数谱

从表1数据中可以看到四种动作的广义维数谱有明显的差别。因而广义维数谱，可以作为肌电信号的特征用于不同动作的分类。

步骤四，以步骤三获取的特征向量输入支持向量机进行分类识别，获得识别结果。

每个动作对应60组六维特征向量，从中任取20组特征向量作为训练集输入SVM进行训练，伸腕、屈腕、展拳、握拳4个动作对应的训练目标分别为1，2，3和4；其余40组作为测试集输入支持向量机，这4个动作对应的测试目标也分别为1，2，3和4。若作为测试的特征向量输入SVM后的输出结果与测试目标一致，说明对测试的动作进行了正确归类，否则是错误的。

对于支持向量机分类器，本实施例采用了3种不同的核函数，即linear核、rbf核和poly核。分别用三种不同核函数对同样的特征值进行了实验，结果表明采用 rbf核（径向基核）的分类效果最好。因此本发明采用rbf核函数设计SVM分类器，这里还涉及到参数的选择，包括控制对错分样本惩罚的程度的可调参数，以及rbf核函数中的待定参数，这两个参数的选择直接影响分类效果。为了有效确定这两个参数的值，、以固定间隔逐步增加参数的值，来搜索使模式识别率最高的参数值。

Claims

1.基于经验模态分解和分形的表面肌电信号模式识别方法，其特征在于该方法包括如下步骤：

步骤(1).获取人体上肢肌电信号样本数据，具体是：首先通过肌电信号采集仪拾取人体上肢肌电信号，再采用空域相关滤波方法对含有干扰噪声的肌电信号进行消噪；

步骤(2).将步骤(1)获取的肌电信号进行经验模态分解，得到代表肌电信号细节构成的多层窄带信号，即内在模态函数；

步骤(3).对步骤(2)所获得的内在模态函数进行多重分形分析，提取其广义维数谱，作为肌电信号多模式识别的特征向量；

步骤(4).以步骤(3)获取的特征向量输入对支持向量机进行分类识别，获得识别结果；

所述的多重分形分析，提取其广义维数谱，具体如下：

设肌电信号时间序列为{x₁,x₂,...x_n}，嵌入维数m，时间延迟为τ，则重构相空间为：

X_i＝(x_i,x_i+τ,…,x_i+(m-1)τ) i＝1,2,…,n-(m-1)τ (1)

式中，X_i为相空间中的第i个相点，相空间中的点数为N＝n-(m-1)τ；构造其相空间后，以m维相空间集{X_i}中的一点X_I＝{x_I1,x_I2,...,x_Im}作为参考点，计算另外N-1点与它的距离，则可统计出落于以点X_I为中心，以小标量r为半径的体积元中的点的个数，从而得到q阶关联积分C^q(m,r)

其中f(x)为Heaviside阶跃函数；

计算广义维数谱D_q算法步骤如下：

1)把一维时间序列信号转换到相空间，求出相点间的最大伸展距离d_max；

2)依据d_max的大小确定计算关联积分C^q(m,r)时r增幅，确定需计算的小标量r初始个数l,r_h为单调上升，h＝1,2,...l；根据公式(2)求得C^q(m,r_h),计算数据点对(ln(r_h),ln(C^q(m,r_h)))，为表达方便，令为(u_h,v_h)；

3)根据下式计算数据组的斜率：

4)增加l的值，根据公式(3)重新计算k值，若连续多次，k的值变化极小，则停止计算，以最后几次所得k值的平均，作为肌电信号的广义维数谱D_q。