CN106108897A - 一种基于经验模式分解的肌电信号滤波方法 - Google Patents
一种基于经验模式分解的肌电信号滤波方法 Download PDFInfo
- Publication number
- CN106108897A CN106108897A CN201610573087.0A CN201610573087A CN106108897A CN 106108897 A CN106108897 A CN 106108897A CN 201610573087 A CN201610573087 A CN 201610573087A CN 106108897 A CN106108897 A CN 106108897A
- Authority
- CN
- China
- Prior art keywords
- electromyographic signal
- mode function
- signal
- intrinsic mode
- function imf
- 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.)
- Pending
Links
- 238000000354 decomposition reaction Methods 0.000 title claims abstract description 58
- 238000000034 method Methods 0.000 title claims abstract description 57
- 238000001914 filtration Methods 0.000 title claims abstract description 41
- 238000005070 sampling Methods 0.000 claims abstract description 29
- 238000000605 extraction Methods 0.000 claims description 16
- 230000008569 process Effects 0.000 claims description 9
- 238000012549 training Methods 0.000 claims description 7
- 238000012545 processing Methods 0.000 claims description 6
- FGUUSXIOTUKUDN-IBGZPJMESA-N C1(=CC=CC=C1)N1C2=C(NC([C@H](C1)NC=1OC(=NN=1)C1=CC=CC=C1)=O)C=CC=C2 Chemical compound C1(=CC=CC=C1)N1C2=C(NC([C@H](C1)NC=1OC(=NN=1)C1=CC=CC=C1)=O)C=CC=C2 FGUUSXIOTUKUDN-IBGZPJMESA-N 0.000 claims description 3
- 230000007704 transition Effects 0.000 claims description 3
- 239000007787 solid Substances 0.000 claims description 2
- 230000015572 biosynthetic process Effects 0.000 claims 1
- 230000007547 defect Effects 0.000 abstract description 2
- 210000003205 muscle Anatomy 0.000 description 8
- 238000010586 diagram Methods 0.000 description 4
- 238000005516 engineering process Methods 0.000 description 3
- 230000002452 interceptive effect Effects 0.000 description 3
- 230000009471 action Effects 0.000 description 2
- 238000004364 calculation method Methods 0.000 description 2
- 230000007812 deficiency Effects 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 230000003993 interaction Effects 0.000 description 2
- 230000000717 retained effect Effects 0.000 description 2
- 238000004804 winding Methods 0.000 description 2
- 241000208340 Araliaceae Species 0.000 description 1
- 235000005035 Panax pseudoginseng ssp. pseudoginseng Nutrition 0.000 description 1
- 235000003140 Panax quinquefolius Nutrition 0.000 description 1
- 230000006978 adaptation Effects 0.000 description 1
- 230000003044 adaptive effect Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 238000013473 artificial intelligence Methods 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000003745 diagnosis Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000005670 electromagnetic radiation Effects 0.000 description 1
- 235000008434 ginseng Nutrition 0.000 description 1
- 238000002372 labelling Methods 0.000 description 1
- 230000035800 maturation Effects 0.000 description 1
- 238000012805 post-processing Methods 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 238000003672 processing method Methods 0.000 description 1
- 230000004044 response Effects 0.000 description 1
Classifications
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/24—Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
- A61B5/316—Modalities, i.e. specific diagnostic methods
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/24—Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
- A61B5/316—Modalities, i.e. specific diagnostic methods
- A61B5/389—Electromyography [EMG]
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/68—Arrangements of detecting, measuring or recording means, e.g. sensors, in relation to patient
- A61B5/6801—Arrangements of detecting, measuring or recording means, e.g. sensors, in relation to patient specially adapted to be attached to or worn on the body surface
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
- A61B5/7203—Signal processing specially adapted for physiological signals or for diagnostic purposes for noise prevention, reduction or removal
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
- A61B5/7225—Details of analog processing, e.g. isolation amplifier, gain or sensitivity adjustment, filtering, baseline or drift compensation
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
- A61B5/7235—Details of waveform analysis
- A61B5/725—Details of waveform analysis using specific filters therefor, e.g. Kalman or adaptive filters
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61F—FILTERS IMPLANTABLE INTO BLOOD VESSELS; PROSTHESES; DEVICES PROVIDING PATENCY TO, OR PREVENTING COLLAPSING OF, TUBULAR STRUCTURES OF THE BODY, e.g. STENTS; ORTHOPAEDIC, NURSING OR CONTRACEPTIVE DEVICES; FOMENTATION; TREATMENT OR PROTECTION OF EYES OR EARS; BANDAGES, DRESSINGS OR ABSORBENT PADS; FIRST-AID KITS
- A61F2/00—Filters implantable into blood vessels; Prostheses, i.e. artificial substitutes or replacements for parts of the body; Appliances for connecting them with the body; Devices providing patency to, or preventing collapsing of, tubular structures of the body, e.g. stents
- A61F2/50—Prostheses not implantable in the body
- A61F2/68—Operating or control means
- A61F2/70—Operating or control means electrical
- A61F2/72—Bioelectric control, e.g. myoelectric
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- General Health & Medical Sciences (AREA)
- Biomedical Technology (AREA)
- Veterinary Medicine (AREA)
- Public Health (AREA)
- Animal Behavior & Ethology (AREA)
- Heart & Thoracic Surgery (AREA)
- Medical Informatics (AREA)
- Pathology (AREA)
- Physics & Mathematics (AREA)
- Molecular Biology (AREA)
- Surgery (AREA)
- Biophysics (AREA)
- Signal Processing (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Artificial Intelligence (AREA)
- Psychiatry (AREA)
- Physiology (AREA)
- Power Engineering (AREA)
- Cardiology (AREA)
- Oral & Maxillofacial Surgery (AREA)
- Transplantation (AREA)
- Vascular Medicine (AREA)
- Measurement And Recording Of Electrical Phenomena And Electrical Characteristics Of The Living Body (AREA)
Abstract
本发明提供一种基于经验模式分解的肌电信号滤波方法，包括如下步骤：采集多个肌电信号，按照一定时间长度的采样窗口和采样步长对处于多个检测通道的所述肌电信号进行采样形成多个肌电信号样本。对所述肌电信号样本根据经验模式分解方法分解成多个固有模态函数IMF和分解余项。从多个所述固有模态函数IMF中选取至少一个进行相加重构，得到去除干扰后的所述肌电信号。本发明克服了数字滤波器滤波势必会滤除有用的肌电信号、人为主观设置滤波参数导致的滤波过深或滤波不足的缺陷。
Description
技术领域
本发明涉及肌肉电信号的预处理算法，尤其涉及一种基于经验模式分解的肌电信号滤波方法，具体涉及肌肉电信号中工频干扰、运动伪迹以及基线漂移的滤除算法，应用领域包括肌电信号处理、手势识别、肌电义肢、康复工程、生物机械等领域。
背景技术
可穿戴式设备，尤其是基于肌肉电信号原理的可穿戴式设备，已逐渐被大众关注和接受。肌肉电信号在被应用于控制领域以来，经历了近百年的发展，一直被研究应用于医学诊断和生物机械领域。随着生物医学技术、人工智能技术的发展，使用肌电信号进行手势识别的方法被提出并不断探索。肌肉电信号比较微弱，幅值在μ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 ModeDecomposition，简称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，可以得到很好的分解效果。
综上所述，本文中应用了具体个例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处，综上所述，本说明书内容不应理解为对本发明的限制，本发明的保护范围应以所附的权利要求为准。
Claims (10)
1.一种基于经验模式分解的肌电信号滤波方法，其特征在于，包括如下步骤：
采集多个肌电信号，按照一定时间长度的采样窗口和采样步长对处于多个检测通道的所述肌电信号进行采样形成多个肌电信号样本；
对所述肌电信号样本根据经验模式分解方法分解成多个固有模态函数IMF和分解余项；
从多个所述固有模态函数IMF中选取至少一个进行相加重构，得到去除干扰后的所述肌电信号。
2.如权利要求1所述的方法，其特征在于，在所述采样窗口内包括64、128、256、512或1024个采样点。
3.如权利要求1或2所述的方法，其特征在于，在底层嵌入式实现时，所述采样窗口内包括的采样点个数与所述采样频率数值的比值不超过0.3。
4.如权利要求1所述的方法，其特征在于，从多个所述固有模态函数IMF中选取至少一个进行相加重构，包括：
将所述肌电信号样本分解后的所述固有模态函数IMF按照频率的高低依次排列；
去除含有工频干扰频率的所述固有模态函数IMF及所有小于所述工频干扰频率的所述固有模态函数IMF，对剩余的所述固有模态函数IMF相加进行重构；或者，
保留所有高于工频干扰频率的所述固有模态函数IMF，并且在所有小于所述工频干扰频率的所述固有模态函数IMF中，保留频率最高的所述固有模态函数IMF，对保留的所述固有模态函数IMF相加进行重构。
5.如权利要求4所述的方法，其特征在于，对所述固有模态函数IMF相加进行重构，包括：
使用所述分解后的第1个所述固有模态函数IMF对信号进行相加重构；或者，
使用所述分解后的第1个和第2个所述固有模态函数IMF对信号进行相加重构；或者，
使用所述分解后的第1个、第2个和第4个所述固有模态函数IMF对信号进行相加重构。
6.如权利要求1所述的方法，其特征在于，在从多个所述固有模态函数IMF中选取至少一个进行相加重构之后，还包括：
使用数字滤波器的方案对所述肌电信号样本进行滤波，得到数字滤波信号；
对采集到的所述肌电信号不经过任何处理，得到原始肌电信号；
分别统计所述数字滤波信号、所述原始肌电信号及经过经验模式分解相加重构操作之后的所述肌电信号的识别率。
7.如权利要求6所述的方法，其特征在于，分别统计所述数字滤波信号、所述原始肌电信号及经过经验模式分解相加重构操作之后的所述肌电信号的识别率，包括：
对所述数字滤波信号进行特征提取、训练和识别处理，得到手势的识别结论，并将所述识别结论与产生所述肌电信号的实际手势进行比较；如果所述识别结论与所述实际手势相一致，则识别成功；如果所述识别结论与所述实际手势不一致，则识别失败；
对所述原始所述肌电信号进行特征提取、训练和识别处理，得到手势的识别结论，并将所述识别结论与产生所述肌电信号的实际手势进行比较；如果所述识别结论与所述实际手势相一致，则识别成功；如果所述识别结论与所述实际手势不一致，则识别失败；
对所述经过所述经验模式分解相加重构操作之后的所述肌电信号进行特征提取、训练和识别处理，得到手势的识别结论，并将所述识别结论与产生所述肌电信号的实际手势进行比较；如果所述识别结论与所述实际手势相一致，则识别成功；如果所述识别结论与所述实际手势不一致，则识别失败。
8.如权利要求1所述的方法，其特征在于，按照一定时间长度的采样窗口和采样步长对处于多个检测通道的所述肌电信号进行采样形成多个肌电信号样本，包括：
步骤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。
9.如权利要求8所述的方法，其特征在于，对所述肌电信号样本根据经验模式分解方法分解成多个固有模态函数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。
10.如权利要求9所述的方法，其特征在于，所述判断IMF的停止条件的阈值σ取值为0.2。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201610573087.0A CN106108897A (zh) | 2016-07-20 | 2016-07-20 | 一种基于经验模式分解的肌电信号滤波方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201610573087.0A CN106108897A (zh) | 2016-07-20 | 2016-07-20 | 一种基于经验模式分解的肌电信号滤波方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN106108897A true CN106108897A (zh) | 2016-11-16 |
Family
ID=57290816
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201610573087.0A Pending CN106108897A (zh) | 2016-07-20 | 2016-07-20 | 一种基于经验模式分解的肌电信号滤波方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN106108897A (zh) |
Cited By (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107049310A (zh) * | 2017-02-22 | 2017-08-18 | 广东工业大学 | 一种基于经验模态分解的对肌电图的预处理方法 |
| CN108742611A (zh) * | 2018-04-08 | 2018-11-06 | 南京医科大学 | 一种随机动态参数电刺激条件下的自主肌电提取方法 |
| CN109558006A (zh) * | 2018-11-23 | 2019-04-02 | 武汉灏存科技有限公司 | 无线分布式肢体动作捕捉设备 |
| CN109998541A (zh) * | 2019-04-19 | 2019-07-12 | 中国地质大学（武汉） | 一种基于经验模态全变分的表面肌电信号去噪方法及系统 |
| CN110974180A (zh) * | 2019-12-25 | 2020-04-10 | 索思（苏州）医疗科技有限公司 | 基于母体生理电信号的宫缩检测装置及方法 |
| CN115251957A (zh) * | 2022-09-26 | 2022-11-01 | 深圳市心流科技有限公司 | 一种逐级调整采样频率的方法、装置及存储介质 |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102961203A (zh) * | 2012-12-10 | 2013-03-13 | 杭州电子科技大学 | 基于emd样本熵的表面肌电信号识别方法 |
| CN103393484A (zh) * | 2013-07-31 | 2013-11-20 | 刘洪运 | 用于电子耳蜗的语音处理方法 |
| CN103610466A (zh) * | 2013-10-17 | 2014-03-05 | 杭州电子科技大学 | 基于emd排列组合熵的肌电信号跌倒识别方法 |
| CN104702244A (zh) * | 2013-12-05 | 2015-06-10 | 中国科学院深圳先进技术研究院 | 基于eemd算法滤除肌电信号检测中工频干扰的自适应滤波器 |
| CN105326500A (zh) * | 2014-08-13 | 2016-02-17 | 华为技术有限公司 | 基于表面肌电信号的动作识别方法和设备 |
-
2016
- 2016-07-20 CN CN201610573087.0A patent/CN106108897A/zh active Pending
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102961203A (zh) * | 2012-12-10 | 2013-03-13 | 杭州电子科技大学 | 基于emd样本熵的表面肌电信号识别方法 |
| CN103393484A (zh) * | 2013-07-31 | 2013-11-20 | 刘洪运 | 用于电子耳蜗的语音处理方法 |
| CN103610466A (zh) * | 2013-10-17 | 2014-03-05 | 杭州电子科技大学 | 基于emd排列组合熵的肌电信号跌倒识别方法 |
| CN104702244A (zh) * | 2013-12-05 | 2015-06-10 | 中国科学院深圳先进技术研究院 | 基于eemd算法滤除肌电信号检测中工频干扰的自适应滤波器 |
| CN105326500A (zh) * | 2014-08-13 | 2016-02-17 | 华为技术有限公司 | 基于表面肌电信号的动作识别方法和设备 |
Non-Patent Citations (2)
| Title |
|---|
| 席旭刚等: "《基于EMD自相关的表面肌电信号消噪方法》", 《仪器仪表学报》 * |
| 时改杰 等: "《基于EMD分解的表面肌电信号动作模式识别》", 《振动与冲击》 * |
Cited By (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107049310A (zh) * | 2017-02-22 | 2017-08-18 | 广东工业大学 | 一种基于经验模态分解的对肌电图的预处理方法 |
| CN107049310B (zh) * | 2017-02-22 | 2020-01-17 | 广东工业大学 | 一种基于经验模态分解的对肌电图的预处理方法 |
| CN108742611A (zh) * | 2018-04-08 | 2018-11-06 | 南京医科大学 | 一种随机动态参数电刺激条件下的自主肌电提取方法 |
| CN108742611B (zh) * | 2018-04-08 | 2023-10-13 | 南京医科大学 | 一种随机动态参数电刺激条件下的自主肌电提取方法 |
| CN109558006A (zh) * | 2018-11-23 | 2019-04-02 | 武汉灏存科技有限公司 | 无线分布式肢体动作捕捉设备 |
| CN109998541A (zh) * | 2019-04-19 | 2019-07-12 | 中国地质大学（武汉） | 一种基于经验模态全变分的表面肌电信号去噪方法及系统 |
| CN110974180A (zh) * | 2019-12-25 | 2020-04-10 | 索思（苏州）医疗科技有限公司 | 基于母体生理电信号的宫缩检测装置及方法 |
| CN115251957A (zh) * | 2022-09-26 | 2022-11-01 | 深圳市心流科技有限公司 | 一种逐级调整采样频率的方法、装置及存储介质 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN106108897A (zh) | 一种基于经验模式分解的肌电信号滤波方法 | |
| CN106874833B (zh) | 一种振动事件的模式识别方法 | |
| CN108806708A (zh) | 基于计算听觉场景分析和生成对抗网络模型的语音降噪方法 | |
| CN112885372B (zh) | 电力设备故障声音智能诊断方法、系统、终端及介质 | |
| CN106340309A (zh) | 一种基于深度学习的狗叫情感识别方法及装置 | |
| CN101515200B (zh) | 基于瞬态视觉诱发脑电的目标选择方法 | |
| CN105841961A (zh) | 一种基于Morlet小波变换和卷积神经网络的轴承故障诊断方法 | |
| CN110718235B (zh) | 异常声音检测的方法、电子设备及存储介质 | |
| CN108288015A (zh) | 基于时间尺度不变性的视频中人体动作识别方法及系统 | |
| CN106228200A (zh) | 一种不依赖于动作信息采集设备的动作识别方法 | |
| CN108958482B (zh) | 一种基于卷积神经网络的相似性动作识别装置及方法 | |
| CN110276357A (zh) | 一种基于卷积神经网络的验证码识别方法 | |
| CN106328151A (zh) | 一种环噪消除系统及其应用方法 | |
| CN109256118A (zh) | 基于生成式听觉模型的端到端汉语方言识别系统和方法 | |
| CN106095101A (zh) | 基于节能机制的人体行为识别方法及客户端 | |
| CN106889984A (zh) | 一种心电信号自动降噪方法 | |
| CN109308903A (zh) | 语音模仿方法、终端设备及计算机可读存储介质 | |
| CN110322882A (zh) | 一种生成混合语音数据的方法及系统 | |
| CN103761974A (zh) | 一种人工耳蜗 | |
| CN114897025A (zh) | 一种人体姿态识别模型建立方法和人体姿态识别方法 | |
| CN117095694A (zh) | 一种基于标签层级结构属性关系的鸟类鸣声识别方法 | |
| CN104897180A (zh) | 用于桥梁监测信号的预处理方法 | |
| CN104049869B (zh) | 一种数据处理方法及装置 | |
| CN113633296A (zh) | 反应时间预测模型构建方法、装置、设备及可读存储介质 | |
| CN113450816A (zh) | 一种基于深度神经网络的语音主动降噪方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| TA01 | Transfer of patent application right | ||
| TA01 | Transfer of patent application right |
Effective date of registration: 20181207 Address after: 10 000 07 Group B 7, 3rd floor, Zhangzizhong Road, Dongcheng District, Beijing Applicant after: Zhang Wendong Address before: 710000 Room 204, Arc Tower, 60 West Avenue, New Industrial Park, Xi'an High-tech Zone, Shaanxi Province Applicant before: XI'AN ZHONGKE BIQI INNOVATION TECHNOLOGY CO., LTD. |
|
| RJ01 | Rejection of invention patent application after publication | ||
| RJ01 | Rejection of invention patent application after publication |
Application publication date: 20161116 |