一种基于肌音信号的下肢假肢步态识别方法
技术领域
本发明属于模式识别领域，涉及一种下肢假肢步态识别方法，尤其涉及一种基于肌音信号的下肢假肢步态识别方法。
背景技术
佩戴智能假肢是截肢患者恢复基本运动能力的唯一途径，目前较为成熟的智能假肢有德国Otto Bock公司的C-Leg、Genium系列假肢，冰岛Ossur公司的Rheo Knee假肢等；由于下肢运动的特殊性和复杂性，步态识别是智能假肢研究的核心技术之一。
根据信号源的不同，步态识别信号源主要包括两种：
1)采集下肢关节角度、惯导信息、足底压力信息等生物力学信号识别步态，例如加速度传感器、角度传感器、脚底压力传感器、位置传感器等，但是这些信息只有在佩戴者完成一定的动作之后才能采集到，在识别相应运动模式时不可避免的会产生一定的延迟；
2)表面肌电信号(sEMG)作为一种人体生物信号，包含了大量的神经信号，结合模式识别算法可以解析出相应的步态模式；人体运动时肌肉振动产生的肌音信号(MMG)和肌电信号相似，都是可以反映肌肉状态的信号源；虽然对肌音信号的研究相对肌电信号晚，但肌音信号与肌电信号相比具有其自身优势；如：肌音信号的采集利用压电传感器或者加速度传感器即可完成，对传感器要求较低，而肌电信号的采集需要表面导电涂层，成本和要求较高；另外，肌音信号原始生理信号比肌电信号原始生理信号强度大，并且受到外界和体表的电磁干扰小，更易采集。
由于肌音信号优势显著，近些年国内外一些学者对其特点及应用进行了深入研究；Alves等利用上肢两通道肌音信号对3种模式手部动作进行模式识别，准确率达到了89％±2％；目前对肌音信号的研究主要集中于上肢肌力研究和上肢手臂及手部动作的模式识别，而对步态识别的研究不多。
从上世纪90年代，支持向量机、神经网络等模式识别算法开始被应用到步态分类中；支持向量机(SVM)是由Vapnik在1995年提出的一种通过建立决策超平面，实现样本分类的方法；由于其多样性、鲁棒性和免费可行的工具箱等优点被广泛地应用到人机交互的模式识别中；Lucas等人基于多通道sEMG信号的离散小波变换，并利用SVM的方法实现了步态的分类；支持向量机(SVM)是一种识别较好的肌电假肢步态分类器，但SVM也存在一些缺陷：SVM在步态识别过程中，存在惩罚因子难以确定的缺陷；神经网络在肌电模式分类和肌电信号处理方面，进行了相当多的研究，并获得了非常有价值的结果，但是神经网络常常面临权值参数人工设置繁琐、训练速度慢的缺陷。
发明内容
本发明的目的在于提供一种基于肌音信号的下肢假肢步态识别方法，能够克服SVM在步态识别过程中惩罚因子难以确定的缺陷，还能够克服神经网络权值参数人工设置繁琐、训练速度慢的缺陷。
为了实现上述目的，本发明采用以下技术方案：
一种基于肌音信号的下肢假肢步态识别方法，包括以下步骤：
步骤1：采集大腿截肢患者的下肢肌音信号并对所述肌音信号进行预处理；具体的：
所述的大腿截肢患者为从坐骨结节下10cm至膝关节间隙上8cm范围内的大腿截肢患者，所述肌音信号为股直肌、股内侧肌、股外侧肌不同步态下的肌音信号，所述不同步态为平地行走、上楼、下楼、上坡、下坡5种不同步态；
步骤1.1：分别对大腿截肢患者的股直肌、股内侧肌、股外侧肌在平地行走、上楼、下楼、上坡、下坡5种不同步态下的肌音信号进行采集，并记为：X股直肌i、X股内侧肌i、X股外侧肌i，其中，i＝1,2,3,4,5，分别表示平地行走、上楼、下楼、上坡、下坡5种不同步态，X股直肌i、X股内侧肌i、X股外侧肌i分别表示股直肌、股内侧肌、股外侧肌的肌音信号；
步骤1.2：分别对X股直肌i、X股内侧肌i、X股外侧肌i进行预处理，以获得预处理后的X’股直肌i、X’股内侧肌i、X’股外侧肌i，所述的预处理为滤波预处理；
步骤2：利用差分HOC提取X’股直肌i、X’股内侧肌i、X’股外侧肌i的步态识别的特征值样本集，具体的：
步骤2.1：计算构建肌音信号的特征向量W时选择的最高差分阶数θ，具体的：利用10位大腿截肢患者的股直肌肌音信号，计算得到最高差分阶数θ＝13；
步骤2.2：利用θ＝13时的差分HOC提取X’股直肌i、X’股内侧肌i、X’股外侧肌i的5种步态识别特征值样本集，具体的：
步骤3：将步骤2获得的步态识别的特征值样本集随机分为两组，一组为训练集，一组为测试集，利用训练集计算最优的用于相关向量机的核函数，建立优化的相关向量机多分类步态识别模型，具体的，包括以下步骤：
步骤3.1：通过蝙蝠算法优化相关向量机分类参数，得到最优相关向量机分类器，具体的，包括以下步骤：
步骤3.1.1：定义蝙蝠算法,具体的：
定义蝙蝠i在t时刻的位置为
fi＝fmin+(fmax-fmin)β；
其中，fi表示蝙蝠i发出声波的频率，fi∈[fmin,fmax]，fmin、fmax分别表示蝙蝠i发出声波的频率的最小值和最大值；β∈[0,1]，x*表示t时刻全局搜索过程中蝙蝠i的最优位置；
蝙蝠i位置更新公式如下
xnew＝xold+ηAt；
式中：xold表示所选取的最优解；η是一个随机数，η∈[-1,1]，At是所有蝙蝠在时刻t时发出声波的平均响度；xnew表示对最优解施加一个随机扰动之后得到的新解；
蝙蝠i发出声波的响度Ai和发射速率ri的调节计算公式为：
式中：α表示声波响度增加系数，α∈[0,1]；γ3表示一个大于0的常量；
步骤3.1.2：定义相关向量机，具体的：
给定输入训练集
其中，xi是训练集的第i个样本，xi∈R，ti∈{-1,+1}是类别标号.
相关向量机的核函数采用高斯核函数
步骤3.1.3：通过蝙蝠算法(BA)优化相关向量机(RVM)分类参数，计算最优的用于相关向量机(RVM)的核函数，具体的，包含以下步骤：
步骤3.1.3.1：选择训练集，所述的训练集即为从步骤2获得的步态识别的特征值样本集随机分为两组中的训练集，
步骤3.1.3.2：设定蝙蝠算法(BA)参数：最大迭代次数N＝100，d＝5，fmax＝3，fmin＝0，ri0＝0.0001，α＝0.75，γ3＝0.96，η＝0.5；
步骤3.1.3.3：设置核函数参数σ的取值范围，构建初始种群；
步骤3.1.3.4：以核函数参数σ下的相关向量机(RVM)步态识别准确率作为适应度值，所述的适应度值即为种群中每个蝙蝠个体的适应度；
步骤3.1.3.5：定义并计算核函数参数σbest；
进行迭代更新操作，并且计算迭代更新后的相关向量机(RVM)在更新后的核函数参数σ下的步态识别准确率；若更新后的相关向量机(RVM)在迭代更新后的核函数参数σ下的步态识别准确率最高，则将迭代更新后的核函数参数σ记为核函数参数σbest；若迭代更新后的相关向量机(RVM)在迭代更新后的核函数参数σ下的步态识别准确率不是最高，则保留以前的核函数参数σbest；
步骤3.1.3.6：最大迭代次数设为终止条件，若不满足终止条件，则返回步骤3.1.3.5继续进行迭代更新操作；当满足终止条件，即迭代次数达到最大，记录此时的σbest作为最优核函数参数，相关向量机(RVM)采用此时的最优核函数参数σbest；
步骤3.1.3.7：依据最优核函数参数σbest，建立优化的相关向量机多分类步态识别模型；
步骤4：利用步骤3.1.3.7获得的优化的相关向量机多分类步态识别模型对测试集进行识别，获取步态识别结果。
步骤1.2中所述的滤波预处理的方法为采用2-50Hz的FIR带通滤波器对X股直肌i、X股内侧肌i、X股外侧肌i进行滤波预处理。
步骤1中所述的采集大腿截肢患者的下肢肌音信号的方法为利用TD-3压电加速度传感器采集大腿截肢患者的下肢肌音信号，利用USB-6211数据采集卡存储TD-3压电加速度传感器采集的大腿截肢患者的下肢肌音信号。
本发明的有益效果：
本发明所述的一种基于肌音信号的下肢假肢步态识别方法，将表征时间序列的差分HOC方法用于提取不同步态识别的特征值样本集，结合通过蝙蝠算法优化的相关向量机建立了多分类步态识别模型，显著的提高了下肢假肢步态识别的准确率。
附图说明
为了更清楚地说明本发明具体实施方式或现有技术中的技术方案，下面将对具体实施方式或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施方式，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本发明的方法流程图。
具体实施方式
下面将结合附图对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
如图1所示：本发明所述的一种基于肌音信号的下肢假肢步态识别方法，包括以下步骤：
步骤1：采集大腿截肢患者的下肢肌音信号并对所述肌音信号进行预处理；具体的，
由于采集的是大腿截肢患者的残肢的肌音信号，所以需要考虑残肢的存留情况，本发明中选取大部分大腿截肢患者能保存的肌肉；具体的，本发明中所述的大腿截肢患者为从坐骨结节下10cm至膝关节间隙上8cm范围内的大腿截肢患者，所述肌音信号为股直肌、股内侧肌、股外侧肌不同步态下的肌音信号，所述不同步态为平地行走、上楼、下楼、上坡、下坡5种不同步态；
进一步的，选用北京颐松有限公司的TD-3压电加速度传感器采集大腿截肢患者的下肢肌音信号；数据采集卡选用美国国家仪器公司的USB-6211，USB-6211数据采集卡通过USB总线供电，具有5V直流电压输出端口，可直接对传感器模块供电，省去了传感器供电模块；
具体的，包括以下步骤：
步骤1.1：分别对大腿截肢患者的股直肌、股内侧肌、股外侧肌在平地行走、上楼、下楼、上坡、下坡5种不同步态下的肌音信号进行采集，并记为：X股直肌i、X股内侧肌i、X股外侧肌i，其中，i＝1,2,3,4,5，分别表示平地行走、上楼、下楼、上坡、下坡5种不同步态，X股直肌i、X股内侧肌i、X股外侧肌i分别表示股直肌、股内侧肌、股外侧肌的肌音信号；
步骤1.2：分别对X股直肌i、X股内侧肌i、X股外侧肌i进行预处理，以获得预处理后的X’股直肌i、X’股内侧肌i、X’股外侧肌i，所述的预处理为滤波预处理；
由于肌音信号频率主要为低于50Hz的低频部分，在400-450Hz具有高频干扰，因此，为了更好地保留有用信号，采用2-50Hz的FIR带通滤波器对X股直肌i、X股内侧肌i、X股外侧肌i进行滤波预处理，获得预处理后的X’股直肌i、X’股内侧肌i、X’股外侧肌i；
步骤2：利用差分HOC提取X’股直肌i、X’股内侧肌i、X’股外侧肌i的步态识别的特征值样本集，具体的：
为了便于本领域技术人员进一步理解本发明的技术方案，下面将对HOC和差分HOC做进一步说明：
过零计数是对一个原始时间序列振荡特性的一种检测，当这个原始时间序列被一簇滤波器滤波后，再对其进行过零计数，则称这种过零计数为相对于原始时间序列的HOC；因此，HOC即为对原始时间序列和预滤波处理序列的过零计数；本发明用{Lθ(.),θ∈Θ}表示一簇滤波器；参数空间Θ为参数θ的一维或多维可数集；对时间序列滤波的过程由Zl(θ)＝Lθ(Zl)来表示；{Dθ,θ∈Θ}表示相应的θ阶HOC值；HOC值Dθ定义为
其中，N为原始时间序列的长度，Xi(θ)定义为
最后得到肌音信号的特征向量W为
W＝[D1,D2,...,Dθ],1＜θ＜K；
式中：K表示所估计的HOC序列的最高阶数；θ表示构建特征向量W时选择的最高差分阶数；
定义延迟算子B：
BZl＝Zl-1，则差分滤波器可由下式构造
Yl(θ)＝Zl(θ)＝(1-B)θZl＝b0(θ)Zl+b1(θ)Zl-1+bθ-1(θ)Zl-θ+1；
其中，θ＝0,1,...,K，Yl(θ)表示滤波后的时间序列，K表示所估计的HOC序列的最高阶数，θ表示构建肌音信号的特征向量W时选择的最高差分阶数，Zl表示原始时间序列X股直肌i、X股内侧肌i、X股外侧肌i；系数bm(θ)由下式计算
式中：
步骤2.1：计算构建肌音信号的特征向量W时选择的最高差分阶数θ，具体的：
使用差分HOC对肌音信号的特征向量W进行提取时，特征向量W的最高差分阶数θ需要结合送入相关向量机(RVM)后的结果得到；
本发明利用10位大腿截肢患者的股直肌肌音信号，计算得到在最高差分阶数θ＝13时，步态识别准确率最高，因此，选择最高差分阶数θ＝13，具体的：
当θ＜6时，步态识别准确率在65％以下；当θ继续增大时，识别准确率有所提高；10组数据中,8组数据在θ＝13时识别准确率最高；有1位受试者在θ＝14时识别准确率最高；有1位受试者在θ＝15时识别准确率最高；随着θ从θ＝13向上递增，步态识别准确率整体有下降的趋势，因此本发明取θ＝13；
步骤2.2：利用θ＝13时的差分HOC提取X’股直肌i、X’股内侧肌i、X’股外侧肌i的5种步态识别特征值样本集，具体的：
步骤3：将步骤2获得的步态识别的特征值样本集随机分为两组，一组为训练集，一组为测试集，利用训练集计算最优的用于相关向量机(RVM)的核函数，具体的，包括以下步骤：
步骤3.1：通过蝙蝠算法(BA)优化相关向量机(RVM)分类参数，得到最优相关向量机(RVM)分类器，具体的，包括以下步骤：
步骤3.1.1：定义蝙蝠算法(BA),具体的：
蝙蝠算法(BA)是一种新型进化算法，在蝙蝠算法中，每一个蝙蝠表示一个可行解，蝙蝠的猎物表示优化问题的最优解，在本发明中，蝙蝠的猎物表示求解核函数的最优值；定义蝙蝠i在t时刻的位置为xit、速度为vit；在d维的觅食空间里，蝙蝠i在t+1时刻的位置和速度的更新公式如下：
fi＝fmin+(fmax-fmin)β；
其中，fi表示蝙蝠i发出声波的频率，fi∈[fmin,fmax]，fmin、fmax分别表示蝙蝠i发出声波的频率的最小值和最大值；β∈[0,1]，x*表示t时刻全局搜索过程中蝙蝠i的最优位置；
蝙蝠i位置更新公式如下
xnew＝xold+ηAt；
式中：xold表示所选取的最优解；η是一个随机数，η∈[-1,1]，At是所有蝙蝠在时刻t时发出声波的平均响度；xnew表示对最优解施加一个随机扰动之后得到的新解；
觅食过程中，为了平衡全局搜索和局部搜索过程，蝙蝠发射的声波响度和射速率需要自动调节；蝙蝠i发出声波的响度Ai和发射速率ri的调节计算公式为：
式中：α表示声波响度增加系数，α∈[0,1]；γ3表示一个大于0的常量；ri0表示初始发射速率；
步骤3.1.2：定义相关向量机(RVM)，具体的：
给定输入训练集
其中，xi是训练集的第i个样本，xi∈R，ti∈{-1,+1}是类别标号.
本发明采用高斯核函数
步骤3.1.3：通过蝙蝠算法(BA)优化相关向量机(RVM)分类参数，计算最优的用于相关向量机(RVM)的核函数，具体的，包含以下步骤：
步骤3.1.3.1：选择训练集，所述的训练集即为从步骤2获得的步态识别的特征值样本集随机分为两组中的训练集，
步骤3.1.3.2：设定蝙蝠算法(BA)参数：最大迭代次数N＝100，d＝5，fmax＝3，fmin＝0，ri0＝0.0001，α＝0.75，γ3＝0.96，η＝0.5；
步骤3.1.3.3：设置核函数参数σ的取值范围，构建初始种群；
步骤3.1.3.4：以核函数参数σ下的相关向量机(RVM)步态识别准确率作为适应度值，所述的适应度值即为种群中每个蝙蝠个体的适应度；
步骤3.1.3.5：定义并计算核函数参数σbest；
进行迭代更新操作，并且计算迭代更新后的相关向量机(RVM)在更新后的核函数参数σ下的步态识别准确率；若更新后的相关向量机(RVM)在迭代更新后的核函数参数σ下的步态识别准确率最高，则将迭代更新后的核函数参数σ记为核函数参数σbest；若迭代更新后的相关向量机(RVM)在迭代更新后的核函数参数σ下的步态识别准确率不是最高，则保留以前的核函数参数σbest；
步骤3.1.3.6：最大迭代次数设为终止条件，若不满足终止条件，则返回步骤3.1.3.5继续进行迭代更新操作；当满足终止条件，即迭代次数达到最大，记录此时的σbest作为最优核函数参数，相关向量机(RVM)采用此时的最优核函数参数σbest；
步骤3.1.3.7：依据最优核函数参数σbest，建立优化的相关向量机多分类步态识别模型；
步骤4：利用步骤3.1.3.7获得的优化的相关向量机多分类步态识别模型对测试集进行识别，获取步态识别结果。
为了便于本领域技术人员理解本发明的技术方案，下面将以具体实施例对本发明的技术方案作进一步描述：
步骤1：采集大腿截肢患者的肌音信号并对所述肌音信号进行预处理；
步骤1.1：采样对象为10名穿戴智能假肢的大腿截肢患者，年龄(28±5)岁，体重(68.0±5.0)kg，身高(169.0±5.0)cm；采样对象在实验前被告知实验目的，采样对象在实验前24h没有进行过任何形式的剧烈运动；为避免行走速度影响实验结果，采样对象在平地行走、上楼、下楼、上坡、下坡时以正常速度行走；利用TD-3压电加速度传感器采集肌音信号X股直肌i、X股内侧肌i、X股外侧肌i；由于肌音信号的有用信息大多集中在0Hz-100Hz，因此在进行肌音信号采集时，采样率为500Hz；其中，每个采样对象每种步态采集60组数据，则每个采样对象5种步态采集300组数据；
步骤1.2：由于肌音信号频率主要为低于50Hz的低频部分，在400-450Hz具有高频干扰，因此，为了更好地保留有用信号，利用2-50Hz的FIR带通滤波器对X股直肌i、X股内侧肌i、X股外侧肌i进行滤波处理，以获得预处理后的X’股直肌i、X’股内侧肌i、X’股外侧肌i；
步骤2：利用θ＝13时的差分HOC提取X’股直肌i、X’股内侧肌i、X’股外侧肌i的5种步态识别特征值样本集，具体的：
步骤3：随机选择步骤2获得的每种步态识别特征值样本集中的2/3作为训练集，剩余1/3作为测试集；并利用训练集训练得到通过蝙蝠算法(BA)优化相关向量机(RVM)后的相关向量机(RVM)的核函数最优核函数参数σbest；依据最优核函数参数σbest，建立优化的相关向量机多分类步态识别模型；
步骤4：利用步骤3.1.3.7获得的优化的相关向量机多分类步态识别模型对测试集进行识别，获取步态识别结果。
步骤4：将步骤3中的测试集输入步骤3获得的以σbest作为最优核函数参数的相关向量机(RVM)，获得测试集的步态识别结果。
在优化前和优化后的相关向量机(RVM)中的步态识别结果如表1所示：
表1优化前和优化后的相关向量机(RVM)中的步态识别结果对比表
由表1可以看出：通过蝙蝠算法(BA)优化后的相关向量机(RVM)对步态识别的识别率明显优于没有优化的相关向量机(RVM)。
最后应说明的是：以上各实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述各实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的范围。
