本发明涉及肌肉电信号的预处理算法,尤其涉及一种基于经验模式分解的肌电信号滤波方法,具体涉及肌肉电信号中工频干扰、运动伪迹以及基线漂移的滤除算法,应用领域包括肌电信号处理、手势识别、肌电义肢、康复工程、生物机械等领域。
背景技术:
可穿戴式设备,尤其是基于肌肉电信号原理的可穿戴式设备,已逐渐被大众关注和接受。肌肉电信号在被应用于控制领域以来,经历了近百年的发展,一直被研究应用于医学诊断和生物机械领域。随着生物医学技术、人工智能技术的发展,使用肌电信号进行手势识别的方法被提出并不断探索。肌肉电信号比较微弱,幅值在μV~mV级。肌肉电信号的采集过程中,主要受到如下因素的干扰:周围电子设备的电磁辐射产生的工频干扰;皮肤与电极接触位置的移动等原因造成的运动伪迹;电子设备长时间工作造成的基线漂移。上述干扰会严重影响肌电信号的后续处理,因此对肌电信号滤波就存在着必然性。
现有的滤波算法多采用数字滤波器,采用高通数字滤波器去除基线漂移以及低频的运动伪迹,采用多个带通数字滤波器去除50Hz工频干扰及其谐波成分。数字滤波器方法成熟,但是缺点也是显而易见的:数字滤波器滤波频带较宽,滤除干扰的同时,势必会滤除有用的肌电信号;此外,人为设置参数,受主观影响较大,可能出现滤波过深或不足。这些缺陷都将影响信号后续的分析识别。
因此,需要设计一种数据驱动的自适应滤波算法,最大限度去除干扰的同时,保留原始有用信号。
技术实现要素:
因此,为解决现有技术存在的技术缺陷和不足,本发明提供一种基于经验模式分解进行滤波的肌电信号预处理方法。
本发明提供一种基于经验模式分解滤波的肌电信号预处理方法,包括如下步骤:
采集多个肌电信号,按照一定时间长度的采样窗口和采样步长对处于多个检测通道的肌电信号进行采样形成多个肌电信号样本;
对肌电信号样本根据经验模式分解方法分解成多个固有模态函数IMF和分解余项;
从多个固有模态函数IMF中选取至少一个进行相加重构,得到去除干扰后的肌电信号。
进一步地,在本发明提供的一种基于经验模式分解滤波的肌电信号预处理方法中,在采样窗口内包括64、128、256、512或1024个采样点。
进一步地,在本发明提供的一种基于经验模式分解滤波的肌电信号预处理方法中,在底层嵌入式实现时,采样窗口内包括的采样点个数与采样频率数值的比值不超过0.3。
进一步地,在本发明提供的一种基于经验模式分解滤波的肌电信号预处理方法中,从多个固有模态函数IMF中选取至少一个进行相加重构,包括:
将肌电信号样本分解后的固有模态函数IMF按照频率的高低依次排列;
去除含有工频干扰频率的所述固有模态函数IMF及所有小于所述工频干扰频率的所述固有模态函数IMF,对剩余的所述固有模态函数IMF相加进行重构;或者,
保留所有高于工频干扰频率的所述固有模态函数IMF,并且在所有小于所述工频干扰频率的所述固有模态函数IMF中,保留频率最高的所述固有模态函数IMF,对保留的所述固有模态函数IMF相加进行重构。
进一步地,在本发明提供的一种基于经验模式分解滤波的肌电信号预处理方法中,对剩余的固有模态函数IMF进行相加重构,包括:
使用固有模态函数IMF中的第1个固有模态函数IMF对信号进行相加重构;或者,
使用固有模态函数IMF中的第1个和第2个固有模态函数IMF对信号进行相加重构;或者,
使用固有模态函数IMF中的第1个、第2个和第4个固有模态函数IMF对信号进行相加重构。
进一步地,在本发明提供的一种基于经验模式分解滤波的肌电信号预处理方法中,在从多个固有模态函数IMF中选取至少一个进行相加重构之后,还包括:
使用数字滤波器的方案对肌电信号样本进行滤波,得到数字滤波信号;
对采集到的肌电信号不经过任何处理,得到原始肌电信号;
分别统计数字滤波信号、原始肌电信号及经过经验模式分解相加重构操作之后的肌电信号的识别率。
进一步地,在本发明提供的一种基于经验模式分解滤波的肌电信号预处理方法中,分别统计数字滤波信号、原始肌电信号及经过经验模式分解相加重构操作之后的肌电信号的识别率,包括:
对数字滤波信号进行特征提取、训练和识别处理,得到手势的识别结论,并将识别结论与产生肌电信号的实际手势进行比较;如果识别结论与实际手势相一致,则识别成功;如果识别结论与实际手势不一致,则识别失败;
对原始肌电信号进行特征提取、训练和识别处理,得到手势的识别结论,并将识别结论与产生肌电信号的实际手势进行比较;如果识别结论与实际手势相一致,则识别成功;如果识别结论与实际手势不一致,则识别失败;
对经过经验模式分解相加重构操作之后的肌电信号进行特征提取、训练和识别处理,得到手势的识别结论,并将识别结论与产生肌电信号的实际手势进行比较;如果识别结论与实际手势相一致,则识别成功;如果识别结论与实际手势不一致,则识别失败;
进一步地,在本发明提供的一种基于经验模式分解滤波的肌电信号预处理方法中,按照一定时间长度的采样窗口和采样步长对处于多个检测通道的肌电信号进行采样形成多个肌电信号样本,包括:
步骤a、设置采样窗口长为L,时间序列为x(t),步长为L;
步骤b、初始化处理,令ri-1(t)=x(t),i=1;其中,ri(t)为第i个分解余项;
步骤c、初始化,令hj-1(t)=ri-1(t),j=1;
步骤d、确定hj-1(t)的所有极大值点和极小值点;
步骤e、用三次样条曲线分别对极大值点和极小值点进行插值得到信号的上下包络线,计算上下包络的平均值mj-1(t);
步骤f、计算hj(t)=hj-1(t)-mj-1(t);
步骤g、若则ci(t)=hj(t);否则,j=j+1,并重复执行步骤d。
进一步地,在本发明提供的一种基于经验模式分解滤波的肌电信号预处理方法中,对肌电信号样本根据经验模式分解方法分解成多个固有模态函数IMF和分解余项,包括如下步骤:
步骤h、ri(t)=ri-1(t)-ci(t),若ri(t)极值点数不小于2个,则i=i+1,转到步骤c重复执行;否则,分解结束;其中,ri(t)是第i个余项,hj(t)为第j个过渡量,mj(t)为hj(t)上下包络的平均值,σ为判断固有模态函数IMF的停止条件的阈值,且0.2≤σ≤0.3,ci(t)为第i个IMF。
进一步地,在本发明提供的一种基于经验模式分解滤波的肌电信号预处理方法中,判断IMF的停止条件的阈值σ取值为0.2。
本发明克服了数字滤波器滤波势必会滤除有用的肌电信号、人为主观设置滤波参数导致的滤波过深或滤波不足的缺陷。
附图说明
为了更清晰地说明本发明或现有技术的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作一简单地介绍。显而易见地,下面描述中的附图是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他附图。下面将结合附图,对本发明的具体实施方式进行详细的说明。
图1为本发明实施例提供的一种基于经验模式分解滤波的肌电信号预处理方法的流程示意图;
图2为本发明实施例提供的一种基于经验模式分解滤波的肌电信号预处理方法的时序示意图;
图3为本发明实施例提供的经验模式分解的算法流程示意图。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚,下面将结合本发明的附图,对本发明的技术方案进行清楚、完整的描述。显然,所描述的实施例是本发明一部分实施例,而不是全部的实施例。基于本发明的实施例,本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例,都属于本发明的保护范围。
请参见图1,图1为本发明实施例提供的一种基于经验模式分解滤波的肌电信号预处理方法的流程示意图,该方法可以包括如下步骤:
S101、采集多个肌电信号,按照一定时间长度的采样窗口和采样步长对处于多个检测通道的肌电信号进行采样形成多个肌电信号样本;
S102、对肌电信号样本根据经验模式分解方法分解成多个固有模态函数IMF和分解余项;
S103、从多个固有模态函数IMF中选取至少一个进行相加重构,得到去除干扰后的肌电信号。
本发明提出了一种数据驱动的基于经验模式分解(Empirical Mode Decomposition,简称EMD)的自适应肌电滤波算法。该算法基于数据的局部时间特征尺度,把采集到的样本信号分解为一系列的固有模态函数(Intrinsic Mode Function,简称IMF)和分解余项,通过重构固有模态函数能够有效的去除工频干扰、运动伪迹干扰以及基线漂移。
进一步地,因为可穿戴设备对数据处理的实时性有较高要求,对数据的采样窗口的长度有一定的要求。在本发明的实施例中,在采样窗口内包括64、128、256、512或1024个采样点,也就是说在底层嵌入式实现时,窗口长度L可以是64、128、256、512或1024。
窗口选择越小,响应时间越快,但是相同数据,处理的段数就越多,计算量和计算的复杂度就越大。同时滤波完还有后续处理,考虑到实时性,本发明的实施者可以灵活设置窗口长度L的具体数值,在采样数据的实时性和计算量以及计算复杂度之间达到平衡。
进一步地,在本发明的实施例中,采样窗口内包括的采样点个数与采样频率数值的比值不超过0.3。这是因为人机交互的实时性:人在延迟时间大于300ms时,会有一定的延迟感。因此在采样频率为F的情况下,L/F应小于0.3,即小于0.3s。例如,窗口内有128个点,采用频率为1000Hz,此时,L/F=128/1000=0.128<0.3,就能满足人机交互的实时性要求。
进一步地,在本发明的实施例中,从多个固有模态函数IMF中选取至少一个进行相加重构,具体方法可以为:
将肌电信号样本分解后的固有模态函数IMF按照频率的高低依次排列;去除含有工频干扰频率的所述固有模态函数IMF及所有小于所述工频干扰频率的所述固有模态函数IMF,对剩余的所述固有模态函数IMF相加进行重构;或者,
保留所有高于工频干扰频率的所述固有模态函数IMF,并且在所有小于所述工频干扰频率的所述固有模态函数IMF中,保留频率最高的所述固有模态函数IMF,对保留的所述固有模态函数IMF相加进行重构。
在本发明的实施例中,50Hz的工频干扰是规律的振荡,必然存在于一个或多个固有模态函数中;运动伪迹属于不规律的低频信号,存在于多个低频固有模态函数中;基线漂移代表了信号的趋势走向,存在于分解后的余项。
通过对采集的肌电信号进行经验模式分解,把信号分解为一系列的由高频到低频的固有模态函数,去除工频干扰和运动伪迹干扰所在的固有模态函数,同时去除余项,把剩余固有模态函数重构,就可以实现对肌电信号的滤波。
进一步地,在本发明的实施例中,在具体选择需要重构的固有模态函数IMF的方式如下:
50Hz的工频干扰信号存在于分解后的第2个或者第3个或者两者都有的固有模态函数IMF中。
使用分解后的第1个固有模态函数IMF对信号进行相加重构。
或者,使用分解后的第1个和第2个固有模态函数IMF对信号进行相加重构。
或者,使用分解后的第1个、第2个和第4个固有模态函数IMF对信号进行相加重构。
进一步地,本发明提供的实施例中,在从多个固有模态函数IMF中选取至少一个进行相加重构之后,还可以使用数字滤波器的方案对肌电信号样本进行滤波,得到数字滤波信号;
对采集到的肌电信号不经过任何处理,得到原始肌电信号;
分别统计数字滤波信号、原始肌电信号及经过经验模式分解相加重构操作之后的肌电信号的识别率。具体的,可以采用如下方式:
对数字滤波信号进行特征提取、训练和识别后处理,得到手势的识别结论,并将识别结论与产生肌电信号的实际手势进行比较;如果识别结论与实际手势相一致,则识别成功;如果识别结论与实际手势不一致,则识别失败;
对原始肌电信号进行特征提取、训练和识别处理,得到手势的识别结论,并将识别结论与产生肌电信号的实际手势进行比较;如果识别结论与实际手势相一致,则识别成功;如果识别结论与实际手势不一致,则识别失败;
对经过经验模式分解相加重构操作之后的肌电信号进行特征提取、训练和识别处理,得到手势的识别结论,并将识别结论与产生肌电信号的实际手势进行比较;如果识别结论与实际手势相一致,则识别成功;如果识别结论与实际手势不一致,则识别失败。
例如,肌肉执行握拳动作,对握拳动作产生的肌肉电信号分别采用数字滤波器的去干扰方式、不去任何干扰的方式及经过经验模式分解滤波得到并重构固有模态函数IMF的方式,分别进行特征提取、训练和识别后续处理,分别得到识别结论。如果识别结论就是握拳,则识别成功;否则识别失败。
图2为本发明实施例提供的一种基于经验模式分解滤波的肌电信号预处理方法的时序示意图。在本发明提供的实施例中,按1kHz的采样频率对信号进行采样,信号的采样窗口大小设置为L=128。在1kHz的采样频率下,L耗时128ms,较好地满足了对实时性的要求。
本实施例中选择了5个实验者(LHC、YQ、LN、LZH和YG)分别采集了休息、伸掌、握拳、内摆及外摆5个手势动作8个通道的肌电信号,并分别进行经验模式分解,得到各个样本的固有模态函数IMF。
随后统计这5个实验者休息和外摆两个手势状态肌电数据的50Hz工频干扰以及小于50Hz的信号存在的固有模态函数IMF所在的位置,如表1所示:
表1
如表1所示,趋势项是经验模式分解的余项,表示信号的基线漂移;运动伪迹存在于小于50Hz的固有模态函数IMF;标记“宽”的位置表示所在的固有模态函数IMF的频带比较宽。
另外取窗口长度L分别为256、512或1024个点,重复以上步骤,50Hz的工频信号干扰同样存在于第2个,或第3个,或第2和第3固有模态函数IMF中,余项为基线漂移干扰,只不过小于50Hz的运动伪迹存在的固有模态函数IMF增加了。由于50Hz的工频信号所在的固有模态函数IMF位置不变,因此,在本实施例中,窗口长度L的选择不影响最终的干扰去除。
在本实施例中,对5个手势动作的肌电信号进行经验模式分解,用第1个固有模态函数IMF对信号进行重构,并进行特征提取、训练和识别后续处理;用第1个和第2个固有模态函数IMF对信号进行重构,并进行特征提取、训练和识别后续处理;用第1个、第2个和第4个固有模态函数IMF对信号进行重构,并进行特征提取、训练和识别后续处理;用通过数字滤波器滤波后的信号进行特征提取、训练和识别后续处理;用没有经过去除干扰处理的信号进行特征提取、训练和识别后续处理。分别统计上述不同来源的信号进行特征提取、训练和识别后续处理之后5个手势的最高识别率,如下表2。
表2
如表2所示,用第1个固有模态函数IMF重构的滤波信号,识别率高于第1,2个固有模态函数IMF、第1,2,4个固有模态函数IMF重构的滤波信号,并明显高于数字滤波器滤波后信号和没有滤波的信号。本实施例进一步说明了经验模式分解滤波的有效性。
进一步地,图3为本发明实施例提供的经验模式分解的算法流程示意图。将采样后的信号按照经验模式分解的方法,分解成多个固有模态函数IMF和分解余项,分解步骤如下:
步骤a、设置采样窗口长为L,时间序列为x(t),步长为L;
步骤b、初始化处理,令ri-1(t)=x(t),i=1;其中,ri(t)为第i个分解余项;
步骤c、初始化,令hj-1(t)=ri-1(t),j=1;
步骤d、确定hj-1(t)的所有极大值点和极小值点;
步骤e、用三次样条曲线分别对极大值点和极小值点进行插值得到信号的上下包络线,计算上下包络的平均值mj-1(t);
步骤f、计算hj(t)=hj-1(t)-mj-1(t);
步骤g、若则ci(t)=hj(t);否则,j=j+1,并重复执行步骤d。
进一步地,在本发明提供的一种基于经验模式分解滤波的肌电信号预处理方法中,对肌电信号样本根据经验模式分解方法分解成多个固有模态函数IMF和分解余项,包括如下步骤:
步骤h、ri(t)=ri-1(t)-ci(t),若ri(t)极值点数不小于2个,则i=i+1,转到步骤c重复执行;否则,分解结束;其中,ri(t)是第i个余项,hj(t)为第j个过渡量,mj(t)为hj(t)上下包络的平均值,σ为判断固有模态函数IMF的停止条件的阈值,且0.2≤σ≤0.3,ci(t)为第i个IMF。
进一步地,在本发明提供的一种基于经验模式分解滤波的肌电信号预处理方法中,判断IMF的停止条件的阈值σ取值为0.2,可以得到很好的分解效果。
综上所述,本文中应用了具体个例对本发明的原理及实施方式进行了阐述,以上实施例的说明只是用于帮助理解本发明的方法及其核心思想;同时,对于本领域的一般技术人员,依据本发明的思想,在具体实施方式及应用范围上均会有改变之处,综上所述,本说明书内容不应理解为对本发明的限制,本发明的保护范围应以所附的权利要求为准。