神经肌肉电刺激诱发的下肢运动疲劳信息检测与处理技术研究

张 希 明 东 李林枫 张雪君

1(天津大学精密仪器与光电子工程学院，天津 300072)

2(天津医科大学医学影像学院，天津 300072)

引言

目前，因疾病、交通事故、地震、战争及意外伤害而造成的下肢截肢、偏瘫和截瘫的患者数目正以成倍的速度增长。神经肌肉电刺激(neuromuscular electrical stimulation，NMES)作为一种安全、有效的治疗方式，能够恢复患者损伤下肢的运动功能，使他们恢复生活自理能力和参加力所能及的活动。NMES是运用电刺激的手段、用精确的刺激顺序和强度激活瘫痪肌群来帮助患者实现自主运动功能的方法[1]，其应用最早出现在20世纪60年代，用于截瘫患者的站立、行走，还有偏瘫患者的肌力恢复，目前正逐渐应用于临床治疗多种神经肌肉疾病。NMES的应用前景广阔，其控制方式是康复工程领域的重要发展方向[2-3]。

限制NMES广泛应用的一个主要因素就是由其诱发的肌疲劳[4]。作为康复训练中不可避免的一种生理现象，肌疲劳会直接影响刺激的康复效果。通过研究如何探测及评价肌疲劳程度，进一步克服肌疲劳带来的负面影响，对于神经肌肉系统认知、残障康复工程和临床理疗评价等方面，均具有重要的方法学意义。

使用者感觉的反馈对NMES系统的有效控制是十分必要的[5]。由于由于表面肌电图(surface electromyography，sEMG)信号是神经肌肉系统活动时生物电活动的总和，因此探讨sEMG信号的变化，也有助于从理论上了解该系统的基本活动规律，为揭示肌肉活动的神经控制机制提供科学依据。

本研究设计了一套NMES诱发下肢运动疲劳的信息检测系统，采集NMES诱发伸膝运动条件下的膝关节角度和疲劳肌电。从sEMG进行预处理、频谱特征提取等环节，结合膝关节角度变化，对NMES诱发肌疲劳做出正确评价。

1 实验对象与方法

1.1 实验对象

本研究的NMES诱发下肢运动疲劳信息检测实验征集了10名受试者，包括7名男性大学生和3名女性大学生，年龄 22～25 岁(23.5±1.0)岁，体重45 ～75 kg，(64.8 ±19.2)kg，身高 162 ～ 180 cm，(172.2±10.2)cm。所有受试者均自愿参加此项目的测试，健康且无肌肉骨骼病史。受试前24 h内未做剧烈运动，以排除过度运动残余疲劳的影响。

1.2 实验系统

实验系统框图如图1所示。首先，NMES作用于受试者下肢肌群，然后分两路通道进行信号采集：一路由表面电极采集肌电信号，传送至肌电放大器进行放大，并通过AD采集卡进行模数转换，最终将数字信号传送入LabVIEW8.2软件平台进行显示和记录;另一路由角度探测器采集膝关节角度信号，并通过角度传感器和USB接口传送至计算机显示和存储。

NMES仪器选用美国芝加哥伊利诺伊大学与Sigmedics公司合作开发的Parastep-1型。北京祥云计算机技术公司生产的隔离肌电放大器EMG-200 μV 的分辨率 ＜1.5 μV，带宽为 0.5～2000 Hz，增益为5000倍，共模抑制比CMRR＞120 dB。角度采集选用了PASCO测角仪，包括角度探测器、绑带、角度传感器、USB接口和角度采集软件。

1.3 实验方法

图1 神经肌肉电刺激诱发疲劳实验系统框图Fig.1 NMES inducing fatigue experiment system

实验对象坐在实验椅上，髋关节与椅子成90°，小腿自然下垂，股四头肌起到伸膝的主要作用。将NMES作用于股四头肌，正负电极距离6 cm。表面电极采用三点式差动输入电极，正负电极置于股直肌肌腹处，电极与肌纤维平行放置，间距20 mm。参考电极放置在膝关节骨骼位置处。最后用测角仪的弹性绑带将角度探测器绑在腿的侧面，以膝关节处为中心。实验场景如图2所示。

图2 实验场景Fig.2 Experiment scene

为了将NMES诱发疲劳和自主收缩疲劳进行比较，本研究10名受试者均参加刺激组和对照组的实验。刺激组首先对受试者施加一次电刺激来确定刺激级数，即膝关节伸展达到水平时的刺激强度。整个实验过程保持刺激强度为这个级数不变。采用每间隔10 s进行一次持续20 s的电刺激，刺激电流波宽3 ms、频率40 Hz，并在20 s的电刺激后立即进行3 s的自主小腿伸直动作，记录下3 s的sEMG信号。为了使刺激组和对照组实验过程保持一致，对照组每间隔10 s进行一次持续20 s的自主小腿伸直动作，并记录下其后3 s动作的肌电信号。采样率1 kHZ，由Labview软件保存。用 PASCO角度采集软件，记录角度变化数据。

2 数据处理

受试者的sEMG信号和角度信号由软件存储，再利用Matlab软件离线处理。将记录的原始表面肌电信号进行预处理，然后进行频谱分析。肌电信号本质上是一种均值为零的类随机信号，绝大部分能量集中在50～150 Hz之间。因此，肌电信号要通过20 Hz高通滤波器和自适应陷波器来分别滤除低频漂移和50 Hz工频干扰[6]。

sEMG信号的功率谱分析广泛应用于肌肉疾病诊断和肌疲劳检测。传统的谱分析方法是通过傅里叶变换FT(Fourier transform，FT)将时域信号转换为频域信号进行频谱分析，这里使用快速傅里叶变换FFT(fast Fourier transform FFT)，再对功率谱密度提取平均频率MNF(mean frequency，MNF)和中值频率 MDF(median frequency，MDF)指标[7]。

MNF的定义式为

MDF是将功率谱分成上、下两个相等面积区域的频率，其定义式为

式中，f为肌电信号的频率，P(f)为其功率密度谱，f0是功率谱密度的上限频率，即为采样频率的一半。

由于表面肌电信号本质上是非平稳的时变信号，采用时变参数模型可将时变参数用一些基函数的加权和来近似，从而就将线性非平稳问题转化为了线性平稳时不变问题[8]。AR(autoregressive)参数模型法是现代谱分析法中重要的内容，也是肌电信号分析的一个重要方法，肌电信号的性质可以用AR模型的若干系数来表征。零均值、n阶时变参数的AR模型为

式中，et是平稳白噪声过程，零均值，方差为σ2。

3 实验结果

3.1 膝关节角度变化曲线

图3所示为10名受试者NMES诱发伸膝运动的膝关节角度变化数据。可以看到，10名受试者在NMES诱发伸膝动作条件下，膝关节的伸展角度不断减小，从初始的水平90°逐次降低。这说明周期性的刺激诱发了肌疲劳，导致肌力的减退，从而刺激级数达到最大级时膝关节伸展的幅度不能维持最大的伸展水平。

图3 10名受试者NMES诱发伸膝运动的膝关节角度变化及线性回归Fig.3 The knee joint change induced by NMES and linear regression curve of ten subjects

3.2 sEMG信号谱分析

3.2.1 经典谱分析

本研究将NMES诱发疲劳和自主收缩疲劳的频谱与AR模型参量进行对比，所有受试者显示了相同的变化规律。图4显示了其中一名受试者的参加刺激组和对照组实验的MNF和MDF及其线性回归曲线，线性回归的通式为y=b1x+b0。可以看出，对信号做频谱分析，提取的频率特征逐次降低，表明疲劳过程伴随着频谱左移，并且刺激组的下降斜率要远大于对照组，说明了在相同条件下，NMES较自主收缩更容易诱发肌疲劳。

根据前人的研究经验，式(3)中AR模型的阶数n取为4时，对信号的分析和识别性能是最好的[9]。此时式(3)计算得到的4组时变参数，一阶参数 a1是最重要的，是肌肉状态随时间变化的最直接的量[10]。其中一名受试者的a1值的线性回归曲线的图5所示，参数a1随着肌肉的疲劳有逐次变小的趋势。并且，刺激组的下降斜率要远大于对照组，说明了相同条件下NMES较自主收缩更容易诱发肌疲劳。

图4 一受试者NMES诱发伸膝运动的频率变化曲线。(a)平均频率;(b)中值频率Fig.4 The frequency change curve of one subject knee stretch motion induced by NMES.(a)mean frequency;(b)median frequency

图5 一名受试者NMES诱发伸膝运动的a1参数变化曲线Fig.5 The a1 value curve of one subject knee stretch motion induced by NMES

图6和表1显示了实验结束时10名受试者的3种肌电特征参量的变化率，即结束时相对于初始值的下降百分比。可以看出，MNF和MDF变化情况比较接近，受试者中最大分别下降了 43.82%和55.49%;而AR模型的 a1参数下降百分比要明显大于前两者，受试者中最大下降了90.23%。这说明，用时变参数跟踪肌疲劳与经典的频谱参数相比，具有显著性的高分辨灵敏度。

图6 10名受试者的3种肌电特征参量的变化率比较Fig.6 The compare of three sEMG character parameter change of ten subjects

表1 10名受试者的3种肌电特征参量的变化率比较Tab.1 The compare of three sEMG character parameter change of ten subjects

4 讨论

NMES诱发肌疲劳是影响NMES作用的一个主要限制因素，刺激强度是影响肌疲劳的一个重要因素。为了采集股直肌最大收缩时的信号，刺激组首先对受试者施加一次电刺激来确定级数。由于不同个体自身的最大收缩力不同，所以达到膝关节最大伸展时的级数也不同，并且每周期刺激保持这个级数不变。而对照组是受试者自主最大收缩(即小腿伸直)的状态。这样，实验组和对照组采集的都是肌肉最大收缩时的sEMG信号。

探究表征肌疲劳的因素，从而实现对NMES作用过程中肌疲劳的检测，是一个研究热点。随着刺激过程疲劳的增加，直接导致肌力下降[11]。由于sEMG信号中时域和频域指标的变化与肌肉活动及功能状态之间有较好的关联性，故sEMG信号常被用于肌肉活动状态与疲劳状态的分析。持续性最大力量肌肉收缩至肌疲劳后，sEMG幅值呈明显降低[12]。近年来，平均频率和中值频率通常被认为是稳定的肌疲劳指示器，肌疲劳时频谱曲线左移[13]，本研究也发现了这一现象。这一现象与以下因素有关：在持续运动时，ATP耗竭，快运动单位(motor units，MUS)很快疲劳，而以慢运动单位替代;在疲劳时，为维持肌肉张力，加强了运动单位兴奋的同步化;在肌疲劳时，肌内压升高，使血流受阻 、肌肉乳酸积累，引起肌膜兴奋性降低，导致肌纤维传导速度(muscle fiber conduction velocity，MFCV)降低。本研究结果显示，随NMES刺激时间增加，基于傅里叶变换的频域指标 MDF下降率比MNF更大，其中最大的一组为 MDF(55.49%)和 MNF(43.82%)，表明反映肌疲劳方面MDF具有更高的敏感性。扬丹[14]等在研究长负荷诱发肱二头肌疲劳过程中也发现sEMG变化时，MDF随负荷时间延长，MDF下降斜率明显大于 MNF，同样认为 MDF对于反映肌疲劳较MNF更具敏感性。

然而，傅里叶变换是对线性时不变信号的频域分析，sEMG信号具有典型的非稳态信号的基本特征，而且时域分辨率与频域分辨率的局部化矛盾也限制了傅里叶变换的使用。现代谱估计中的参数模型法从最小均方误差拟合角度，对随机信号进行分析，克服了经典谱估计的频率分辨率低和方差性能不好的缺点。与傅里叶变换相比，具有对时间窗口宽度和背景噪声敏感性低的优点，因此已被一些学者用于表面肌电信号的分析中。采用时变参数模型，可将时变参数用一些基函数的加权和来近似，从而就把线性非平稳问题转化为了线性平稳时不变问题。曹玉珍等人应用时变AR模型研究仰卧起坐诱发肌疲劳的 a1参数，它比传统的 MDF对疲劳反应的灵敏度高[10]。本研究证明，在NMES诱发肌疲劳条件下，肌电的a1参量比传统的频谱指标对疲劳反应的灵敏度更高。此方法充分考虑了表面肌电信号的非平稳性，对于数据长度要求不高，可用于在线实时分析，且提取的特征稳定且易于识别。

本研究进一步对NMES诱发收缩和自主收缩疲劳进行比较，探究NMES诱发疲劳的机理。在第二个周期时，刺激组的 sEMG频谱和AR模型参量已表现出明显的下降，而对照组则变化缓慢。比起意志性收缩，NMES作用于骨骼肌所产生疲劳的速度更快。3个因素作用于这一现象：一是在意志性收缩中，运动单位的募集顺序允许在低强度收缩时选择性募集耐疲劳的运动单位[16]，而电激发的收缩即使是低强度刺激，也会使许多快速易疲劳的运动单位被募集起来;二是在意志性收缩中，产生接近最大肌力所需频率似乎比电激发收缩要低得多，较高的频率可引起更快的疲劳;三是进行意志性收缩时，中枢神经系统(central nervous system，CNS)能变换使用运动单位，并且调整它们的放电率，以帮助维持目标水平的用力。但是，这种机制未曾在NMES诱发肌肉收缩中出现。

目前，对于肌疲劳的研究，大多在肌肉静态等长收缩条件下，本研究在NMES诱发下肢动态收缩的条件下，检测了膝关节角度的变化。实验表明，随着周期性的刺激，在达到初始设定的最大刺激强度时，膝关节的角度逐渐减小，说明随着疲劳增加，肌力下降，膝关节不足以维持最大伸展水平。

由于神经肌肉系统的高度复杂性和精细性，以及肌电信号本身具有非稳态性和非线性等特点，故采用传统的线性变化尚不足以反映信号变化的全部特征，因此还将对电刺激致肌疲劳后的肌电信号变化进行非线性分析。今后将疲劳过程中sEMG信号的变化实时地反馈到NMES系统，建立sEMG和肌疲劳关系的模型;通过对肌疲劳的检测，相应地调节NMES的作用参数，使NMES系统获得最佳刺激效果，从而实现最优控制。

5 结论

本研究设计了NMES诱发下肢运动条件下的肌疲劳的检测系统。在NMES诱发伸膝运动疲劳过程中，实时采集膝关节角度和表面肌电信号，应用频谱分析对sEMG信号进行特征提取，表明中值频率、平均频率和AR模型指标可以作为NMES诱发疲劳的检测指标。并且，AR模型法充分考虑了表面肌电信号的非平稳性，在NMES诱发肌疲劳的判定等领域具有很大的应用潜力。

[1]明东，万柏坤.功能性电刺激技术在截瘫行走中的应用研究进展[J].生物医学工程学杂志，2007，24(4)：932-936.

[2]Tepavac D，Schwirtlich L.Detection and Prediction of FES-induced Fatigue[J].Electromyogr Kinesiol，1997，7(1)：39-50.

[3]金德闻，张济川.康复工程学的研究与发展[J].现代康复，2000，4(5)：643-646.

[4]Joseph Mizrahi，Oron Levin，Aviv Aviram，et al.Muscle fatigue in interrupted stimulation：Effect of partial recovery on force and EMG dynamics[J].Electromyogr Kinesiol，1997，7(1)：51-65.

[5]Patrick E. Sensors for use with functional neuromuscular stimulation[J].IEEE Trans on BME，1986，33(1)：256-268.

[6]吴丽萍.面向智能假肢的表面肌电信号采集与处理[D].南京：东南大学，2006.

[7]Merletti R，Knaflitz M，De Luca CJ.Myoelectric manifestations of fatigue in voluntary and electrically elicited contractions[J].J ADDI Phvsiol，1990，3(69)：1810-1820.

[8]Karlsson SB，Gerdle MA.Analyzing surface signals recorded during isokinetic contractions[J].IEEE Eng Med Bio Mag，2001，20(6)：97-105.

[9]罗志增，杨广映.表面肌电信号的 AR参数模型分析方法[J].传感技术学报，2003，12(4)：384-387.

[10]刘洪涛，曹玉珍，谢小波.表面肌电信号的时变 AR模型参数评估肌疲劳程度的研究[J].中国生物医学工程学报，2007，26(4)：493-497.

[11]Mizrahi J.，Levy M.，Ring H.，et al.EMG as an indicator of fatigue in isometrically fes-activated paralyzed muscles[J].IEEE Trans Neural Syst Rehabil Eng，1994，10(2)：57-65.

[12]Tepavac D，Schwirtlich L.Detection and prediction of FES-induced fatigue[J].J Electromyogr Kinesiol.1997，7(1)：39-50.

[13]Merletti R， Conte LR.Advances in processing of surface myoeletric signals：Part l[J].Med Biol Eng Comput，1995，33(3)：362-372.

[14]扬丹.等长负荷诱发肱二头肌疲劳过程中 sEMG信号变化[J].体育与科学，2000，21(5)：27-28.

[15]Fitts RH，Metzger JM，Mechanisms of muscular fatigue[J].Principles of Exercise Biochemistry，1988，5(10)：212-229.