基于经验模板的表面肌电运动单元动作电位序列分解方法研究.doc

基于经验模板的表面肌电运动单元动作电位序列分解方法研究重庆大学硕士学位论文（学术学位）学生姓名：罗万国指导教师：侯文生教授专业：生物医学工程学科门类：工学重庆大学生物工程学院二 O一四年四月 Study on Decomposition of Surface EMG into Motor Unit Action Potential Trains Based on Prior Templates A Thesis Submitted to Chongqing Universityin Partial Fulfillment of the Requirement for theMaster s Degree of EngineeringBySupervised by Prof. Hou Wen-ShengSpecialty: Biomedical EngineeringCollege of Bioengineering ofChongqing University, Chongqing, ChinaApril, 2014 重庆大学硕士学位论文中文摘要摘要表面肌电信号（sEMG）是肌肉活动时所有募集运动单元（MU）产生的运动单元动作电位（MUAP）在表面电极处时空综合叠加的结果，包涵了大量的 MU募集和 MUAP发放信息。sEMG信号分解就是从 sEMG信号中提取主体运动单位动作电位序列（MUAPTs）的过程，分解得到的 MUAP发放信息有助于深入研究神经-肌肉控制系统的调控机理，在临床医学、假肢控制、康复医学、运动医学等领域具有良好的应用前景。目前，sEMG信号分解技术可以大致分为两种类型：一类是盲源分离算法或系统辨识法，另一种是 MUAP形态学方法。由于第一类算法应用于 sEMG信号的基本假设条件并不一定满足，且目前的分解效果并不理想，因此本文基于 MUAP形态特征设计了一种的 sEMG信号分解算法。在分析总结相关文献的基础上，本文结合 MUAP波形常见为双相或三相波形特点，利用 Hermite-Rodriguez函数拟合了 4种时间、幅度可伸缩的 MUAP波形模板；为了减小沿时间轴顺序分割 sEMG信号可能对 MUAP叠加波形识别所带来的影响，本文依据 MU募集的“大小原则”，从整段 sEMG信号中，按照从大到小的顺序逐个剥离出相应的 MUAP波形。值得注意的是，本文在 MUAP波形的识别过程未增加 MUAP发放规律的假设，仅限定了 MUAP发放频率的范围。此外，本文所设计的 sEMG分解算法可对单通道 sEMG信号进行独立分解，克服了对其他通道信息的依赖。为了满足 sEMG分解对信号高信噪比的要求，本文分别采用了分 3个步骤对原始 sEMG信号进行了预处理。首先，本文采用了椭圆带通数字滤波器来消除sEMG信号主频带（20500Hz）以外的部分低频和高频噪声。其次，基于快速独立分量分析算法（FastICA）设计算法实现了工频干扰的分离。最后，采用具有双正交、紧支撑性、近似对称性等优点的 coif2母小波对 sEMG信号进行小波包去噪。实测 sEMG信号的分解结果显示，本文所提 sEMG信号预处理算法不仅有效地滤除工频噪声等噪声，而且较好地保留了 MUAP波形的锐度。由于缺乏 sEMG信号中主体 MUAPTs的先验知识，通常需要专门设计相应的算法准确性验证方案。因此，本文构建了简单的 sEMG信号模型来对本文所提算法的准确性进行验证。本文分别对不同信噪比（5dB、10dB、15dB和 20dB）不同叠加程度（0%、10%、20%和 30%）情况下的 5s长的仿真 sEMG信号（采样率为2kHz）进行了分解，每种情况进行 20组。仿真 sEMG信号分解的结果显示，该算法在噪声水平较高（SNR=20dB）、MUAP叠加程度较轻时(10%)，分解的准确性较I 重庆大学硕士学位论文中文摘要高（90.94% 1.27%）；为进一步验证利用本文算法提取出来的主体 MUAPTs与相应的神经肌肉活动是否具有相关性，本文还将该算法应用于 8名受试者（3组/人）不同手指活动模式下的指浅屈肌多通道（12通道）sEMG信号分解；单通道分解结果显示，高力量水平下 sEMG信号中的主体 MUAPt能够被有效检测和分类；统计结果证实，随着力量水平的增加，MUAP的数目增加；不同大小 MUAP的比重的变化与活动手指和力量水平具有显著的相关性。本文的实验结果初步验证了利用先验模板从sEMG中渐进提取 MUAP的可行性，为 sEMG分解和进一步研究 MU发放规律提供一种新的思路。关键词：表面肌电信号（sEMG），分解，先验模板，运动单元动作电位II 重庆大学硕士学位论文英文摘要ABSTRACT Surface electromyography (sEMG) signal is the result of composite superpositionof motor unit action potentials (MUAP), which are generateed by all recruited motorunits when the muscle is excited, it contains lots of information regarding MUrecruitment and MUAP firing. The decomposition of an sEMG signal is the process ofdecompose sEMG signals into their constituent motor unit action potential trains(MUAPTs), the acquired information of MUAP firing contributes to the research of theneuromuscular regulation mechanism, has a good application prospect in the fields ofclinical medicine and prosthesis control, rehabilitation, sports medicine and so on.At present, the technology of the sEMG signal decomposition can be roughlydivided into two types: One is Blind source separation algorithm or SystemIdentification, the other is based on MUAP shape. Due to the basic assumptions of theformer are not always satisfied when applied to sEMG signal decomposition, and Thecurrent decomposition effect is not ideal, so we designed an algorithm of sEMG signaldecomposition based on the shape of MUAP.Based on the analyzing and summarizing of the related literature, in combinationwith the characteristics of MUAPs that they are common for dual phase or three phasewaveform, 4 kinds of MUAP waveforms with variable time during and amplitude werefitted using Hermite - Rodriguez function; In order to reduce the effects of sEMG signalsequence segmentation along the time on stacked waveform recognition, with the“bigger first” principle, MUAP was automatically extracted from the whole sEMGsignal one by one based on the “size principle” of MU recruitment. It is noticing that,there is no adding assumption of MUAP firing rule during the MUAP waveformrecognition and just limited the scope of MUAP firing frequency. In addition, theproposed sEMG signal decomposition algorithm can decompose the single channelsignal independently. It overcomes the dependence on the information of the otherchannels.In order to obtaining high SNR ratio sEMG signal for decomposition, three stepsof preprocessing of raw sEMG signals were adopted. Firstly, the elliptical bandpassdigital filter was adopted to remove the low and high frequency noise outside of thesEMG signal main frequency zone. Second, the power frequency interference wasseparated based on the fast independent component analysis (FastICA) algorithm. Lastly,coif2 mother wavelet that has the advantage of orthogonality, compactly supportting andIII 重庆大学硕士学位论文英文摘要approximate symmetry was used for wavelet packet denoising. The decompositionresult of real sEMG signals shows that, the proposed preprocessing algorithm can notonly remove the power frequency interference and other noise, but also keep thesharpness of MUAP waveform.Due to the lack of priori knowledge of constituent motor unit action potential trainsin sEMG signals, specially designed accuracy validation protocol is often needed.Therefore, simple sEMG signal model is constructed for the proposed algorithmaccuracyvalidation.SimulatedsEMGsignals(5s)withdifferentSNR(5dB/10dB/15dB/20dB) and different degree of superposition (0%/10%/20%/30%)were decomposed, 20 groups for each case. The decomposition result of SimulatedsEMG signals shows that, the accuracy (90.94%1.27%) of this decompositionalgorithm is higher when the SNR is higher (SNR=20dB) and degree of superposition islower (10%).In order to verify if there is some relation between the MUAPTs and thecorresponding neuromuscular activity, the proposed algorithm was applied to thedecompose in the multi-channel (12 channels)sEMG recorded from 8 subjects (3 groupsper subject) flexor digitorum superficialis muscle in different fingers activity mode. Thedecomposition result of single channel sEMG indicate that ,the main MUAPts of highforce level sEMG can be effectively distinguished and classified; Statistics show that，with the increasing of force level，the number of MUAP increases; The proportion ofdifferent size of MUAP varies significantly with active finger and force level. Theexperimental results preliminary proved that it is feasible to gradually extract MUAPsfrom sEMG signals using prior templates，it also provides a new method for sEMGdecomposition and the evaluation of MU firing patterns.Keywords：Surface electromyography, Decomposition, Prior templates,Motor unit action potentialsIV 重庆大学硕士学位论文目录目录中文摘要.I英文摘要. III1绪论. 11.1课题背景及意义 . 11.2肌电信号的研究现状 . 11.2.1肌内肌电信号的研究现状 . 11.2.2表面肌电信号的研究现状 . 31.3研究目的及内容 . 61.4本文内容安排 . 72表面肌电信号运动单元动作电位检测 . 82.1肌电信号的形成机制及检测方法 . 82.1.1单纤维动作电位 . 82.1.2运动单元动作电位 . 92.1.3运动单元动作电位序列 . 102.1.4肌电信号检测 . 102.2 MUAP的特征与检测分析 . 112.3基于 sEMG信号的 MUAP提取与分析 . 122.3.1 sEMG信号分解算法. 122.3.2分解算法有效性评价 . 172.4本章小节 . 183基于经验模板的 sEMG渐进分解算法设计. 203.1表面肌电信号预处理 . 213.1.1带通数字滤波 . 213.1.2 FastICA去工频干扰. 223.1.3小波包去噪 . 263.2基于经验模板的 sEMG渐进分解算法设计 . 293.2.1模板建立 . 293.2.2渐进分解 . 303.3算法的仿真信号测试 . 323.3.1仿真 sEMG信号构建. 323.3.2仿真信号分解结果 . 33V 重庆大学硕士学位论文目录3.4实测表面肌电信号分解 . 353.5本章小结 . 364手指活动模式对 MUAP发放的影响. 384.1实验数据采集 . 394.2基于 sEMG分解的 MUAP序列提取 . 414.3手指活动模式对 MUAP发放数目的影响 . 424.3.1 MUAP发放总数随力量水平的变化趋势. 424.3.2不同手指活动模式下不同大小 MUAP发放比重变化. 444.4本章小结 . 475总结与展望 . 485.1总结. 485.2展望. 49致谢. 51参考文献. 52附录. 57A.作者在攻读学位期间发表的论文 . 57B.作者在攻读学位期间参与的科研项目 . 57VI 重庆大学硕士学位论文1绪论1 绪论1.1课题背景及意义运动是人体的基本功能之一1。人体完成任何一项运动都伴随着肌肉的收缩或伸展1-5态与身体机能息息相关，研究肌肉的功能状态具有极其重要的生理意义2, 3。当人体处于不同的运动状态时，肌肉的活动状态也不同；肌肉的功能状。骨骼肌肌肉由肌外膜构成的结缔组织包裹，肌肉组织中包含有许多肌束，各个肌束由肌束膜包裹着多个肌纤维组成4。人体的骨骼肌的活动受到运动神经元的控制，占主导作用的 a运动神经元，a运动神经元与其所支配的所有肌纤维构成一个运动单位（Motor Unit，MU），肌肉自主收缩或受到外部刺激时，神经肌肉系统中的运动神经元被激活，所募集 MU发放的运动单位动作电位（Motor Unit ActionPotential，MUAP）经过由肌肉、脂肪及皮肤等组织构成的容积导体后在检测电极处与各种噪声信息综合叠加就形成了肌电信号（Electromyography，EMG）2-4, 6EMG信号与运动中枢对骨骼肌的控制的密切相关，包含了丰富的生理信息，研究EMG信号具有重要的实用价值2-4, 6近年来，随着 EMG信号检测技术的进步和 EMG信号处理方法的快速发展，EMG已广泛应用于神经肌肉控制相关的基础理论研究6-8和临床医学9, 10。、假肢控制11, 12、康复医学13、运动医学等应用领域，并且展现了良好的应用前景6。传统的时域分析、频域分析和时频域的线性分析方法，虽然在一定程度上实现了 EMG信号的特征信息提取，但所提取的信息量有限，只能作为神经对肌肉控制信息的一种简单表征。同样，不断兴起的非线性分析方法虽然能够从 EMG信号中提取出许多传统分析方法不能提取的隐含信息，但是所提取的信息仍然是运动单元活动信息的间接形式。因此，要想深入研究中枢神经系统对运动单元的募集模式，最好的方式还是分解 EMG还原为构成它的主体 MUAP序列。1.2肌电信号的研究现状根据检测方式的不同，可以将 EMG信号分为两种，一种是通过直接插入到肌肉组织内部的针电极或线电极记录的肌内肌电信号 (intramuscular EMG signal，iEMG)，另一种是通过置于皮肤表面的表面电极检测的表面肌电信号(surface EMGsignal，sEMG)4。1.2.1肌内肌电信号的研究现状由于针电极尺寸小，能够与肌纤维直接进行接触，采集到的 iEMG信号信噪比高，可以比较容易地从 iEMG信号中检测到各运动单元发放的 spike信息。因此，1 重庆大学硕士学位论文1绪论早期的研究主要通过提取 iEMG信号的峰峰值、过零点数、中值频率和小波变换系数等时域、频域或时频域的特征参数来研究运动单元的活动情况14。随着高采样率 iEMG信号检测技术的出现，iEMG信号中 MUAP波形细节更易于分辨，区分不同 MU发放的 MUAP波形逐渐成为可能。在此基础上，De Luca研究团队于 1982年首先发表了使用信号处理和模式识别的算法分解 EMG信号的报道15。De Luca等指出 EMG分解就是从 EMG信号中还原出其主体运动单位动作电位序列（Motor Unit Action Potential Trains，MUAPTs）的过程14。图 1. 1所示为 EMG分解示意图，该图概念性地显示了 EMG分解的过程，描述了分解后的EMG信号与单个 MU活动之间的联系。图 1. 1EMG信号分解示意图14Fig.1. 1 schematic diagram of EMG signal decomposition141982年以来，以波士顿大学的De Luca和Le Fever研究团队为代表的众多科研工作者做了大量的EMG分解工作。经过30多年的发展，处理技术己经日趋成熟，许多研究小组都开发了比较完善的 EMG分解系统，如：De Luca研究团队研发的Precision Decomposition II（PD II）系统，Stanford大学的K.C. McGill研究团队研发的EMGLAB软件等4。目前，各种商业的EMG分解系统软件也已广泛应用到科学研究与临床诊断中。总的来说，iEMG信号检测与研究工作的成果是相当显著的，但由于 iEMG自身的局限性，限制了它的应用范围：用针电极提取EMG信号需将电极插入肌肉内部，会给被试造成一定程度的损伤甚至感染；采集iEMG信号需要专业的操作技术和被试的积极配合，在临床应用的过程中，如小儿麻痹症患者和帕金森综合症的检测时，信号的采集难度大；使用针电极进行EMG信号采集的可重复性差；采集iEMG信号一般在低肌肉活动水平下进行，不适用于高力量水平下的2 重庆大学硕士学位论文1绪论肌肉活动情况的研究；针电极的尺寸小，所能反映的肌肉活动区域有限，无法反映整块肌肉的放电情况4, 6因此，具有无创性、便宜性和可重复性等优点的sEMG信号逐渐成为了研究的。焦点。1.2.2表面肌电信号的研究现状sEMG信号是通过粘贴在皮肤表面的表面电极所记录的弱电信号，因此 sEMG信号比 iEMG信号更容易受工频干扰等环境噪声和设备固有噪声的干扰，加上人体组织的容积导体滤波作用，使得 sEMG信号的相关研究的面临着比 iEMG信号更大的挑战4, 6, 16, 17。尽管如此，得益于不断提高的 sEMG信号检测技术和信号处理分析方法，以及不断投入的大量人力与物力，sEMG信号的研究在近年来有了较快的发展，并逐渐在许多基础理论研究和临床诊断方面取代 iEMG信号成为相关研究的重要手段。根据 sEMG信号研究的目的、应用场合以及需要提取的信息量不同，各种信号处理方法都在 sEMG研究中得到了不同程度的应用。目前，sEMG的分析方法主要有传统的线性分析方法（主要有时域分析、频域分析和时频域分析）、非线性分析方法（包括分形分维、李雅普诺夫指数和近似熵）和基于 sEMG分解的 MUAP的特征分析。传统线性分析方法1）时域分析时域分析方法是把 sEMG信号看成是一个时间的函数，通过提取时域特征参数来表征肌肉活动的变化 18，主要用于肌电假肢控制和肌肉疲劳检测等研究。通常应用于肌电假肢控制的时域特征参数主要有均方值、过零点数、方差、绝对值积分、直方图等12；AR模型、过零率、平均整流值和均方根值等19则是肌肉疲劳分析常用的时域特征参数。本人所在实验室也采用均方根值进行了假肢控制的研究11。但是，由于 sEMG信号本身比较微弱，人体组织对其又有衰减作用，再加上环境噪声、设备固有噪声等对信号质量也有很大的影响，因而时域特征分析所提取出的特征参数并不能十分准确反映肌肉的活动状态，从而使得其应用也相应的受到了限制。2）频域分析频域分析主要是对信号进行功率谱分析，对比分析不同肌肉活动模式下信号频谱成分的变化。Ronager等人通过对比分析了 sEMG信号功率谱之间的差异来区分神经肌肉疾病患者与正常人，取得了比较好的区分效果 20。sEMG信号的频谱分析也常用于肌肉疲劳等研究21。值得注意的是，sEMG本身是一种非平稳的时变信号，频谱分析对于肌肉活动随时间变化很快的研究并不适用。3 重庆大学硕士学位论文1绪论3）时频分析针对 sEMG信号的非平稳时变特性，近年来 sEMG信号的时间-频率分析方法引起了人们的关注。其中短时傅立叶变换（ STFT）、维格纳分布（WVD）、乔伊-威廉斯分布（CWD）、小波变换（CWT）和希尔伯特-黄变换（HHT）等时频分析方法已广泛应用于基于 sEMG信号的动作识别、假肢控制和疲劳分析等研究领域21。小波分析在 sEMG信号分析中的应用相对较广。K. Englehart等将小波变换应用到了手势识别与假肢控制的研究，实现了 6种手势的实时识别22，刘晨曦等将小波变换应用于静态收缩时前臂肌肉疲劳的研究 23也得到了有效地应用21。此外，小波变换在疾病诊断方面诚然，传统的线性分析方法在 sEMG研究中扮演了重要的角色，但是由于人体本身是一个复杂的系统，sEMG信号中包含了肌肉活动的复杂信息，单纯的线性分析方法对 sEMG信息提取毕竟是有限。因此，对 sEMG进行复杂性分析成为了一种必然的趋势21。非线性分析继脑电和心电等生理电信号的非线性特征得到了充分的证实之后，一些常用的复杂性分析方法及非线性时间序列分析方法也开始应用到 EMG信号的分析中。目前，在 sEMG信号处理中，非线性分析方法的典型应用包括：1）分形分维定义一个分形体特性就是自相似性，大部分的生理信号呈现统计自相似性。最近的研究发现 sEMG信号中包含许多传统分析方法不能提取的隐含信息。比如，王人成等的研究指出人体上臂肌肉在不同运动模式下进行不同程度的收缩，所采集到 sEMG信号具有不同的分形维24。Gitter等研究表明分形维和肌力有高度的相关性，可用于 MU的募集模式的量化25。Nieminen和 Takala的研究结果显示，随着肌肉疲劳程度的增加，sEMG信号的嵌入维数也随着下降262）李雅普诺夫（Lyapunov）指数。Lyapunov指数是定量评价系统动力学特性的一个重要指标，它表征了系统在相空间中相邻轨道间收敛或发散的平均指数率27。当其应用于 sEMG信号分析时，Lyapunov指数主要用于反映肌肉活动随时间的变化情况。邹晓阳将 Lyapunov指数应用到人体前臂 sEMG的动作模式识别中，对前臂 6类动作进行了有效识别27, 28Dingwell和 Cusumano利用 Lyapunov指数来量化人类行走运动学中的局部动态稳定性21, 29。3）近似熵（ApEn）近似熵（approximate entropy，ApEn）是 Pincus为了定量检测时间序列的复杂性于 1991年提出的一个概念21。由于估算出合理的 ApEn值只需要较短的数据，4 重庆大学硕士学位论文1绪论克服了样本熵、K-S熵、分形维和 Lyapunov指数等参数的局限性，且 ApEn值对噪声的敏感度较低。因此，在 EMG信号与许多生理电信号的非线性分析中普遍以ApEn值作为特征参数。陈伟婷等通过计算 ApEn等熵值来对肌肉疲劳进行了评价21。熊安斌等通过 ApEn值对面瘫患者健康侧与患侧 sEMG信号的不确定性进行了量化与分析，结果显示 ApEn值等表现出显著的统计学差异30另外，其他一些团队还结合多种非线性指标来对 sEMG信号进行处理与分析也取得了比较好的效果31。综上，就目前非线性分析的方法在 sEMG信号研究的现状而言，这类分析方法在疲劳分析、临床疾病初步诊断等方面确实得到了比较好的应用。但是，由于这方面的研究尚浅，能够获得信息量毕竟有限，非线性特征参数通常只能作为一个肌肉功能特性的评价指标，通过进一步与其他分析方法结合或许能够挖掘到更多中枢神经系统对肌肉的控制信息。基于 sEMG分解的 MUAP特征分析sEMG信号是肌肉收缩时所有募集的 MU发放出来的动作电位在表面检测电极处时间和空间上的综合叠加，隐含着不同神经肌肉单元 MU的募集特征。将sEMG信号分解为 MUAPTs，国内外研究团队作了大量的理论研究和实验探索工作。De Luca实验团队从事 EMG信号分解的研究工作的时间最长，后续研究多以其研究方法、结果为基础6。该实验室在 iEMG信号分解的基础上对算法进行了使得的改进以适应 sEMG信号分解15力量水平（最高 100%MVC）的 sEMG分解，所分解的 MU类型数高达 40，并且准确度也高达 97%32，并将其应用于帕金森综合症的研究33, 34认为其准确性仍需改进35。最新报道证实该算法可用于不同肌肉、不同。但是，也有学者。除此之外，还有多种先进信号测量与处理方法应用于sEMG信号的分解研究，典型的算法有人工智能算法、统计决策理论和盲源分离等4, 6, 36。上述算法在特定的条件下都取得了一定的进展。但也各有不足。就人工智能算法和统计决策理论而言，只有对 MUAP具体发放