本发明涉及一种肌电信号特征提取领域，特别涉及一种融合小波包和双谱分析的肌电信号特征提取方法。
背景技术：
表面肌电(surface electromyography,sEMG)信号是肌肉收缩时所产生的电波动，既与肌肉本身的组织生理特性有关，也与神经控制系统有关，它反映了神经肌肉的活动和功能状态。因此，肌电信号已广泛应用于生理医学、康复医学及运动医学等领域的研究，并且成为驱动机器人、控制假肢运动以及功能性电刺激的理想控制信号。
sEMG具有非平稳特性，同时，在采集sEMG信号时易受环境影响而引入高斯噪声。单采用小波变换、AR模型或小波包变换等传统方法进行特征提取时，其正确识别率依然不高，导致指令控制智能硬件时出现误动作。鉴于此，本文结合小波包的时频分析能力和双谱对高斯噪声不敏感的优势，利用小波包对sEMG的各频带进行分解，得到sEMG在不同频带内的能量分布，然后通过对小波包重构的时域信号分析，提取出双谱特征，最后结合这两种方法提取出特征向量。不仅能得到sEMG信号的特征信息，还能抑制噪声，对sEMG信号识别具有很重要的意义。
技术实现要素：
本发明旨在解决以上现有技术的问题。提出了一种提高sEMG信号的识别率的融合小波包和双谱分析的肌电信号特征提取方法。本发明的技术方案如下：
一种融合小波包和双谱分析的肌电信号特征提取方法，其特征在于，包括以下步骤：
S1、输入表面肌电sEMG信号，采用小波包对sEMG信号进行j层分解；S2、选取n‘个有效频带，求出每个频带的能量；S3、将步骤S2求出的频带能量进行归一化处理，作为表面肌电sEMG信号的特征向量；S4、小波包重构sEMG信号的时域信号；S5、对表面肌电sEMG的时域信号进行双谱分析，提取出双谱特征；S6、融合小波包和双谱分析的肌电信号特征向量。
步骤S1中采用小波包对表面肌电sEMG信号进行j层分解，小波包对sEMG进行分解的算法为：
式中为上层小波包分解结果，和为下一级分解结果，j是尺度指标，l为位置指标，n为频率指标，k为变量，h0和h1为分解使用的多分辨率滤波器系数。
进一步的，所述步骤S2：选取n个有效频带，求出每个频带的能量具体包括步骤：假设信号的采样频率为2f，如果对信号进行j层小波包分解，就可以形成2j个等宽频带，每个区间频宽为f/(2j)，小波包分解之后，得到j层小波包系数Cj,m,k，k＝0,1,…,2j-1,m为小波包空间位置标识；
由Parseval能量积分等式，sEMG信号x(t)在时域上的能量为：
根据该能量积分等式关联起来可知，式子与x(t)的小波包变换系数Cj,k具有能量量纲，因此原始信号的能量分布用小波包能量谱来表示是可行的。
进一步的，将步骤S2求出的频带能量进行归一化处理的公式为：
Tj,k表示频带能量；T′j,k表示归一化处理的频带能量。
并将其当作sEMG信号的特征向量，即：
T＝[T′1,T′2,…,T′n]。
T＝[T′1,T′2,…,T′n]。
进一步的，步骤S4采用小波包重构表面肌电sEMG的时域信号，其中小波包重构算法为：
式中为上层小波包分解结果，和为下一级分解结果，j是尺度指标，l为位置指标，n为频率指标，k为变量，h0和h1为分解使用的多分辨率滤波器系数。
进一步的，S5、对表面肌电sEMG的时域信号进行双谱分析，提取出双谱特征，具体包括步骤：
设一离散、平稳时间序列{x(i)}，t＝1,2,…,N，根据信号的离散傅里叶变换，双谱定义为：
Bx(f1,f2)＝E[X(f1)X(f2)X*(f1+f2)]，
式中fi为频率变量；X(fi)为信号的离散傅里叶变换；i＝1,2；E[]为数学期望，双谱估计采用双谱估计间接法，首先估计实验数据{x(i)}的三阶累积量，然后对该累积量进行二维DFT变换，这样就得到随机序列的双谱估计。
进一步的，所述双谱估计间接法算法的具体描述如下：
a、将长度为N的实验数据{x(i)}分成K段，每段有M个数据，即N＝KM；
b、进行去均值操作，去除每段数据的均值，使将分析数据的均值为0；
c、假设{xj(i)}(i＝1,2,…,M；j＝1,2,…,K)为第j段的数据，估计每段数据的三阶累积量：
表示第j段数据的三阶积累量，
k1＝max{1,-m,-n}，k1表示1、-m、-n中最大值；
k2＝min{M,M-m,M-n}；k2表示M、M-m、M-n中最小值；
d、对进行统计平均，得到K组数据的累积量估计，即：
e、对三阶累积量估计进行二维DFT变换，即得双谱估计。
进一步的，所述融合小波包和双谱分析的肌电信号特征向量，其具体步骤为：
每个样本的双谱值为Bj(f1,f2)，j＝1,2,…,n，n为样本数，对每个样本分别进行双谱计算，取平均可得样本双谱为：
接下来再对每个样本求极大值Max(Bj(f1,f2))所对应的即为特征向量，
然后与所述能量特征组合作为新的特征向量：
本发明的优点及有益效果如下：
本发明针对表面肌电信号(sEMG)的非平稳特性导致信号的正确识别率低的问题，提出一种融合小波包和双谱分析相结合的特征提取方法。首先，采用小波包对sEMG的各频带进行分解，得到sEMG在不同频带内的能量分布；其次，通过对小波包重构的时域信号进行双谱分析，提取出双谱特征；最后，融合这两种方法提取出特征向量。这种融合小波包和双谱分析的肌电信号特征提取方法不仅能提取到高时频分辨率的sEMG信号特征信息，还能抑制噪声，从而提高sEMG信号的识别率。
附图说明
图1是本发明提供优选实施例融合小波包和双谱分析的肌电信号特征提取方法流程图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、详细地描述。所描述的实施例仅仅是本发明的一部分实施例。
本发明解决上述技术问题的技术方案是，
图1所示，本发明提供了一种融合小波包和双谱分析的肌电信号特征提取方法，其包括以下步骤：
第一步，采用小波包对sEMG进行j层分解。小波包对sEMG进行分解算法为：
式中为上层小波包分解结果，和为下一级分解结果，j是尺度指标，l为位置指标，n为频率指标，k为变量，h0和h1为分解使用的多分辨率滤波器系数。
第二步，选取n个有效频带，求出每个频带的能量。假设信号的采样频率为2f，如果对信号进行j层小波包分解，就可以形成2j个等宽频带，每个区间频宽为f/(2j)。小波包分解之后，得到j层小波包系数Cj,m,k，k＝0,1,…,2j-1,m为小波包空间位置标识。
由Parseval能量积分等式，sEMG信号x(t)在时域上的能量为：
根据该能量积分等式关联起来可知，式子与x(t)的小波包变换系数Cj,k具有能量量纲，因此原始信号的能量分布用小波包能量谱来表示是可行的。
第三步，将上述的频带能量进行归一化处理：
并将其当作sEMG信号的特征向量，即：
T＝[T′1,T′2,…,T′n]。
T＝[T′1,T′2,…,T′n]。
第四步，小波包重构sEMG的时域信号。小波包重构算法为：
式中为上层小波包分解结果，和为下一级分解结果，j是尺度指标，l为位置指标，n为频率指标，k为变量，h0和h1为分解使用的多分辨率滤波器系数。
第五步，对重构的sEMG的时域信号进行双谱分析，提取出双谱特征。具体算法和步骤为：
设一离散、平稳时间序列{x(i)}，t＝1,2,…,N，根据信号的离散傅里叶变换，双谱定义为：
Bx(f1,f2)＝E[X(f1)X(f2)X*(f1+f2)]，
式中fi为频率变量；X(fi)为信号的离散傅里叶变换；i＝1,2；E[]为数学期望。双谱估计有直接法和间接法，本文采用双谱估计间接法，其核心思想是首先估计实验数据{x(i)}的三阶累积量，然后对该累积量进行二维DFT变换，这样就得到随机序列的双谱估计。该算法的具体描述如下：
a、将长度为N的实验数据{x(i)}分成K段，每段有M个数据，即N＝KM；
b、进行去均值操作，去除每段数据的均值，使将分析数据的均值为0；
c、假设{xj(i)}(i＝1,2,…,M；j＝1,2,…,K)为第j段的数据，估计每段数据的三阶累积量：
k1＝max{1,-m,-n}，
k2＝min{M,M-m,M-n}；
d、对进行统计平均，得到K组数据的累积量估计，即：
e、对三阶累积量估计进行二维DFT变换，即得双谱估计。
第六步，融合小波包和双谱分析的肌电信号特征向量。其具体步骤为：
每个样本的双谱值为Bj(f1,f2)，j＝1,2,…,n，n为样本数，对每个样本分别进行双谱计算，取平均可得样本双谱为：
接下来再对每个样本求极大值Max(Bj(f1,f2))所对应的即为特征向量。
然后与步骤4中的能量特征组合作为新的特征向量：
以上这些实施例应理解为仅用于说明本发明而不用于限制本发明的保护范围。在阅读了本发明的记载的内容之后，技术人员可以对本发明作各种改动或修改，这些等效变化和修饰同样落入本发明权利要求所限定的范围。