目前，上肢運動的疲勞狀態監測，一般單純依賴表面肌電信號（sEMG）對疲勞進行識別和分類，導致結果不穩定，存在一定局限。為此，本文將sEMG信號識別與動作捕捉技術引入到疲勞狀態監測過程中，提出了一種融合改進的肌電疲勞閾值算法與生物力學分析的疲勞分析方法。本研究通過右上肢負載屈肘試驗，同步采集肱二頭肌sEMG信號與上肢動作捕捉數據，并同時運用柏格（Borg）疲勞度主觀自覺量表記錄受試者疲勞感受。然后，將融合改進的肌電疲勞閾值算法和生物力學分析的疲勞分析方法與平均功率頻率（MPF）、譜矩比（SMR）、模糊近似熵（fApEn）、Lempel-Ziv復雜度（LZC）四種單一評價指標疲勞評價方法的試驗結果進行對比。試驗結果表明，本文方法對總體疲勞狀態識別率結果達到98.6%，對輕松、過渡、疲勞三種狀態的識別率分別達到97%、100%、99%，較其他方法更有優勢。本文研究結果證明，本文方法在上肢運動過程中能夠有效預防過度訓練引起的二次損傷，對于疲勞監護具有重要意義。
引用本文： 徐兆, 呂健, 潘偉杰, 何愷倫. 基于表面肌電信號和動作捕捉的上肢運動疲勞分析. 生物醫學工程學雜志, 2022, 39(1): 92-102. doi: 10.7507/1001-5515.202108026 復制
版權信息： ?四川大學華西醫院華西期刊社《生物醫學工程學雜志》版權所有，未經授權不得轉載、改編
引言
人體上肢運動是骨骼肌在中樞神經的刺激下收縮并帶動骨骼以關節為支點所產生的定向運動[1]，持續性的負載上肢運動會使上肢疲勞產生一定程度的積累，導致肌肉骨骼勞損，造成肌肉骨骼失調[2]。上肢疲勞分析在上肢康復訓練、上肢體能訓練、上肢運動輔助設備的設計和優化等方面被廣泛運用[3]。有效的上肢疲勞監測能夠為上肢運動疲勞提供量化參考和評價，能將疲勞特征作為訓練策略制定和調整的重要依據，從而降低上肢勞損的風險。
目前，國內外各研究機構利用表面肌電信號（surface electromyography，sEMG）特征、生物力學分析等手段，在疲勞監測以及分析等方面取得了一定成果[4-7]。例如，Matsumoto等[4]首次提出了應用sEMG信號特征值對骨骼肌負荷進行量化判斷的方法，稱為肌電疲勞閾值（electromyogram fatigue threshold，EMGFT）算法。Kang等[5]對69名普通男性大學生進行負載遞增騎行試驗，分析了股外側肌的EMGFT，并對EMGFT和通氣閾值（ventilatory threshold，VT）進行了相關性分析，結果顯示EMGFT和VT均與疲勞具有高度的相關性。王琳等[6]結合生物力學理論分析了sEMG信號對頸部肌肉疲勞狀態的反映效果。劉曉光等[7]通過對10名健康受試者在負載可調的功率自行車上進行負載遞增騎行試驗，驗證了基于EMGFT來分析肌肉疲勞的有效性。然而，疲勞是一類多維驅動的復雜生理現象，但目前的研究中大多數將sEMG信號作為疲勞判斷的唯一標準，隨著上肢訓練時長的增加，生物電信號信噪比降低[8]，單純依賴sEMG信號對疲勞進行識別和分類會使結果不穩定，存在一定局限性。人的肢體活動是由骨骼肌以及關節共同協作完成的多自由度復雜運動[9-12]，sEMG信號主要用于量化所測肌肉的活動狀態，常用于身體局部肌肉負荷以及疲勞程度的評估，并不能對關節的疲勞程度進行評估。
針對以上問題，本文通過sEMG信號和動作捕捉技術，提出了一種融合EMGFT算法和生物力學分析的疲勞分析方法。通過對15名健康受試者進行右上肢負載曲肘試驗，同步采集肱二頭肌sEMG信號、右上肢動作捕捉數據，并通過柏格（Borg）疲勞度主觀自覺量表記錄受試者主觀疲勞感受，分析肱二頭肌的EMGFT，建立肘關節生物力學模型、計算肘關節切向力方差。本文通過綜合評價上肢骨骼肌以及關節疲勞程度，并與平均功率頻率（mean power frequency，MPF）、譜矩比（spectral moments ratio，SMR）、模糊近似熵（fuzzy approximate entropy，fApEn）、Lempel-Ziv復雜度（Lempel-Ziv complxity，LZC）四種單一評價指標疲勞評價方法進行對比，以驗證本文方法的有效性，本文研究方法如圖1所示。
圖1
基于sEMG信號和動作捕捉的上肢疲勞分析研究方法
Figure1.
Upper limb fatigue analysis method based on sEMG signal and motion capture
1 基于改進的EMGFT算法的肌肉疲勞分析方法
1.1 sEMG原始信號預處理
肌電信號（electromyography，EMG）是神經肌肉系統活動時的生物電信號變化經電極引導、放大、顯示和記錄所獲得的一維電壓時間序列信號，而sEMG信號則是淺層肌肉EMG信號和神經干上電活動在皮膚表面形成的具有非線性特征的綜合生物電效應[13]，sEMG信號通過電極在皮膚表面進行EMG信號的采集，具有非損傷性、實時性、多靶點測量等特點，常用于身體局部肌肉負荷以及疲勞程度的評估[14]。
在進行EMGFT的計算之前，需要對采集到的包含一定噪聲的原始sEMG信號進行預處理，首先設計0～500 Hz低通濾波器對經由sEMG信號采集模塊（EMGduino8通道肌電采集模塊，航奕生物，中國）采集到的原始sEMG信號進行低通濾波，去除高頻干擾；然后采用49.5～50.0 Hz自適應陷波器濾除信號中的工頻和諧波干擾；最后運用小波閾值去噪。由于小波變換中的時與頻的定位特性，可以進行信號的時變譜分析，能夠在各類細節上分析sEMG信號。因此小波閾值去噪是sEMG信號分析的有力工具[15]。本文通過觀察sEMG信號的頻域情況，發現能量在300 Hz以下的部分占全部能量的80%，主要集中在50～220 Hz之間。因此，本研究分解次數選擇4層，如圖2所示為sEMG信號預處理過程的時域圖。
圖2
sEMG信號預處理過程時域圖
Figure2.
Time domain diagram of the process of sEMG signal preprocessing
1.2 改進的EMGFT算法及其原理
EMGFT是指利用sEMG信號及其分析技術來確定運動肌疲勞閾值的方法。與傳統的反映疲勞閾值的生理指標相比，sEMG信號變化特征與肌肉活動狀態及其功能特征之間存在較好的關聯性[16-17]。目前，廣泛認同肌肉運動達到無氧閾值(anaerobic threshold，AT)時，骨骼肌進入疲勞狀態[18-21]。通過對EMGFT進行計算得到的結果，能夠用來測定AT和對應骨骼肌負荷強度。
劉曉光等[7]根據Matsumoto等[4]首次提出的EMGFT算法理論，進行固定時長的負載遞增騎行試驗，將sEMG信號的均方根（root mean square，RMS）值作為特征值，結合最小二乘法分析了股外側肌的EMGFT，驗證了其EMGFT算法的有效性。但是由于劉曉光等[7]試驗設計的特殊性，固定了試驗時長，試驗中采用了固定的時間窗長度（以符號tim表示）與移動窗長度（以符號mov表示），若受試者個體差異導致試驗時長不穩定，固定的tim與mov會造成一定程度的數據損失。為提高肌肉疲勞程度的識別率，本文針對以上問題，本文設置自適應tim與mov對EMGFT算法進行改進。
1.2.1 自適應計算時間窗和移動窗
由于受試者存在個體差異導致疲勞時間不穩定，即sEMG數據集中收集到的數據點數量不穩定。若固定tim與mov會使得時間窗與移動窗超出數據集邊緣，導致計算不完全，造成數據損失，影響試驗結果，因此本文依據具體試驗情況及需求，設置自適應窗計算RMS值。
為提升EMGFT算法疲勞狀態識別準確率，本文每組試驗需要通過自適應窗獲取RMS值50個，即mov需要移動50次，其中設置tim為mov的兩倍長度，則計算mov的計算公式如式(1)所示，而tim的計算公式如式(2)所示：
 |
 |
其中，n為sEMG數據集數據點數量，f為sEMG信號采集模塊采樣頻率。
1.2.2 改進的EMGFT算法流程
改進的EMGFT算法流程如圖3所示，具體算法步驟如下：
圖3
改進的EMGFT算法流程
Figure3.
Flowchart of improved EMGFT algorithm
步驟（1）：使用1.1小節中提到的sEMG信號預處理方法，對通過sEMG信號采集模塊采集到的右上肢肱二頭肌sEMG原始信號進行預處理，建立sEMG數據集。
步驟（2）：根據1.2.1小節中方法自適應計算tim、mov。
步驟（3）：由式(3)以及步驟（2）中計算得到的tim和mov，加窗計算sEMG數據集RMS值(以符號RMS表示)，建立RMS數據集，如式（3）所示：
 |
式（3）中，Nt為tim內sEMG 信號數據點數量，i為sEMG信號數據點序號，Ei為sEMG數據集中第i個sEMG信號。
步驟（4）：將步驟（3）中得到的RMS數據集分為X、Y兩個部分，前5個RMS數據點為X1，剩余RMS數據點為Y1，X1和Y1組成擬合組1。再將前6個RMS數據點作為X2，剩余RMS數據點為Y2，X2和Y2組成擬合組2。直到Yn的RMS數據點數量等于X1的RMS數據點數量時，Xn和Yn組成擬合組n。分別對每組進行一階最小二乘擬合，每組得到兩條擬合直線，然后計算每組兩條擬合直線的斜率乘積。
步驟（5）：計算斜率乘積kXn·kYn最大的一組數據所擬合直線的交點，即為EMGFT。
2 融合改進的EMGFT算法和生物力學分析的疲勞分析流程
2.1 肘關節疲勞分析方法
上肢屈肘動作主要由肱二頭肌群以及肘關節協作完成[22]，僅利用sEMG信號的EMGFT算法只能對上肢屈肘運動中的肌肉疲勞程度進行評價。因此本文在改進的EMGFT算法的基礎上，通過動作捕捉技術采集三維人體肢體運動的位置數據，獲得各關節點絕對坐標，建立上肢球棍模型、計算肘關節空間位置關系，通過角速度與角加速度、切向力變化及其方差對關節疲勞程度進行評價，對照EMGFT算法輸出結果共同評價上肢負載屈肘過程中骨骼肌與關節的疲勞程度，建立基于sEMG信號和動作捕捉技術的上肢康復疲勞評價模型。
2.1.1 建立上肢生物力學模型
本文采用多剛體節段和鉸鏈關節的球棍結構建立人體上肢球棍模型，在球棍模型中對上肢關節角度進行定義。人體額狀面上肢球棍模型如圖4所示。
圖4
人體額狀面上肢球棍模型
Figure4.
Ball and stick model of the upper limb on the frontal plane of the human body
圖4中，點a為肩關節中心點，點b為肘關節中心點，點c為腕關節中心點。c1為屈肘運動過程中腕關節某一點位置。上臂與前臂的夾角θ為肘關節角度。屈肘運動可以近似理解為前臂繞肘關節點（點b）作非勻速運動。
根據生物力學分析以及人體測量學中肢體參數計算公式，對上肢屈肘運動過程進行受力分析，上肢屈肘運動過程受力分析如圖5所示。
圖5
上肢屈肘運動過程受力分析
Figure5.
Stress analysis of elbow flexion in upper extremity
圖5中，P為前臂和手部質心，H為肘關節中心，L1為前臂和手的質心到肘關節中心點的距離，θ為上臂與前臂的夾角，G1為前臂和手部的重力，G2是重物的重力，Fg是手部對重物的支持力，FX為前臂和手繞肘關節點作非勻速運動時的離心力，Fn為前臂和手繞肘關節點作非勻速運動時的切向力。
2.1.2 肘關節疲勞分析方法
為評價上肢負載屈肘運動中的肘關節疲勞程度，本文運用上肢位置信息絕對坐標，計算肘關節角度及其變化情況，計算肘關節切向力方差，以此反映肘關節切向力的穩定性。肘關節疲勞分析方法具體步驟如下：
步驟（1）：三維人體肢體運動的位置數據獲取
根據人體額狀面上肢球棍模型，如圖4所示，于a、ｂ、ｃ三點設置動作捕捉慣性傳感器，記錄a、ｂ、ｃ三點運動過程中的絕對坐標數據。
步驟（2）：肘關節角度計算
根據動作捕捉慣性傳感器獲得的絕對坐標數據，定義 
代表受試者m身上第n個點在t時刻的三維坐標值。如圖4所示，a、b、c是空間中3個點，θ代表肘關節角度，θ計算方法如式（4）所示:
 |
步驟（3）：肘關節角速度計算
由于上肢屈肘運動中肘關節角度θ隨時間t非勻速變化，根據本小節步驟（2）中計算得到的肘關節角度θ加窗計算角速度
，計算方法如式(5)所示：
 |
步驟（4）：肘關節角加速度計算
由于人體上肢屈肘運動并不是完全的勻速運動，因此上肢屈肘運動過程中角速度處于變化中，根據步驟（3）中得到的角速度
加窗計算角加速度a，計算方法如式（6）所示：
 |
步驟（5）：肘關節切向力計算
肘關節切向力反映肘關節負載運動過程中的受力情況。如圖5所示，人體上肢屈肘運動時肘關節切向力Fn處于不規則變化中，根據身體各部重量占體重比例[9]計算前臂和手部的重力G1，G2為重物的重力，由此計算肘關節負載m，計算方法如式（7）所示，通過運動學原理以及本小節步驟（4）中計算得到的角加速度a，計算切向力Fn，計算方法如式（8）所示：
 |
 |
步驟（6）：肘關節切向力方差計算
方差是衡量一組隨機變量離散程度的度量，常用來衡量一批數據的波動大小，反映數據的穩定性[23-24]，由于上肢肘關節處于疲勞狀態時肘關節受力會出現波動，為評價肘關節受力穩定性，根據本小節步驟（5）中得到的肘關節切向力Fn，加窗計算Fn方差s2，計算方法如式（9）所示：
 |
式(9)中，n為肘關節切向力Fn數量，
為肘關節切向力Fn均值。
2.2 評價流程
由于單一sEMG信號無法作為評價關節疲勞的指標，本文融合改進的EMGFT算法和生物力學分析，對上肢負載屈肘運動進行疲勞分析。在上肢屈肘運動過程中，綜合評價肱二頭肌及肘關節的疲勞程度，提出融合EMGFT算法和生物力學分析的疲勞評價流程，在時間序列中，定義達到EMGFT閾值時間點為肌肉進入疲勞狀態，達到關節切向力方差峰值時間點為關節進入疲勞狀態。具體步驟如圖6所示。該評價流程輸出四類結果分別為：P1（肌肉關節均疲勞）、P2（僅肌肉疲勞）、P3（僅關節疲勞）、P4（肌肉關節均未疲勞）。
圖6
融合EMGFT算法和生物力學分析的疲勞分析方法
Figure6.
Fatigue analysis method integrating EMGFT and biomechanics
3 試驗及結果分析
3.1 試驗
為驗證本文所提出的融合改進的EMGFT算法和生物力學分析的疲勞分析方法的有效性，本試驗招募15名健康男性志愿者參與研究，年齡（22.4 ± 1.5）周歲，身高（173.2 ± 5.8）cm，體重（65.3 ± 7.2）kg，體質指數（body mass index，BMI）為22.7 ± 1.8。受試者均身體健康，無肢體損傷等疾病和不良嗜好。數據采集地為現代制造技術教育部重點實驗室，受試者在試驗前均已熟悉本試驗方案和流程且簽署志愿者知情同意書。整個試驗已通過貴州大學人體醫學實驗倫理分委員會審查并獲得倫理審查通知書。受試者在試驗前3月內每周鍛煉次數不超過2次且在試驗前48 h內未進行劇烈運動。sEMG信號采集設備采用sEMG信號采集模塊（EMGduino8通道肌電采集模塊，航奕生物，中國）。動作捕捉設備采用Smartsuit（pro，Rokoko，丹麥）。根據人體生理學結構，同步采集每位受試者右上肢肱二頭肌sEMG信號和右上肢動作捕捉數據，sEMG信號采樣頻率為1 kHz，動作捕捉數據采樣頻率為30 Hz，傳感器位置如圖7所示，S1為右上肢肱二頭肌處的sEMG信號傳感器，S2、S3為sEMG信號參考電極。C1、C2、C3、C4分別為右上肢肩、上臂、前臂以及手處的9自由度慣性測量單元（inertial measurement unit，IMU）運動跟蹤傳感器（Smartsuit pro，Rokoko，丹麥）。
圖7
上肢試驗過程及傳感器放置
Figure7.
Upper limb experiment procedure and sensor placement
受試者測試前刮除右上肢肱二頭肌周圍體毛，用75%酒精擦拭傳感器處皮膚，降低汗液導致電極移動引起的低頻噪聲。受試者保持坐立，右上肢進行負載為7.5 kg的屈肘運動，定義受試者上肢由圖7位置a進行負載屈肘運動到圖7位置b再回到圖7位置a為一次完整屈肘運動。每位受試者進行15次完整的屈肘運動，并在負載屈肘運動期間在測試者幫助下完成Borg量表評分。在Borg量表中將受試者主觀疲勞程度分為1～12評分區間，分別選取Borg量表中 1～4、5～8、9～12 分數段對應輕松、過渡和疲勞 3 種狀態進行相關問題的分析和討論。
3.2 試驗結果分析
利用1.2節所述改進的EMGFT算法分析受試者右上肢肱二頭肌的疲勞時間，如圖8所示為隨機抽取的其中一名受試者RMS值分組擬合過程，X-RMS為Xn擬合直線，Y-RMS為Yn擬合直線。選取所擬合直線斜率乘積最大的一組，獲取該受試者在進行上肢負載屈肘運動時肱二頭肌的EMGFT。
圖8
RMS一階最小二乘擬合過程
Figure8.
RMS first order least squares fitting process
利用2.1.2節中所述肘關節疲勞分析方法，計算分析受試者右上肢肘關節疲勞情況。其中一名受試者進行上肢負載屈肘運動時，手與前臂繞肘關節中心點作非勻速運動時的肘關節角度變化情況，以及該過程中切向力變化情況如圖9所示。
圖9
肘關節角度和切向力變化情況
Figure9.
Changes in elbow angle and tangential force
表1為融合改進的EMGFT算法和生物力學分析的疲勞分析方法輸出結果與Borg量表評分區間對照。如圖10所示，為受試者中隨機3名受試者進行上肢負載屈肘運動時，EMGFT和肘關節切向力方差變化對照。由圖10可知，3名受試者先后出現肱二頭肌EMGFT以及肘關節切向力方差峰值。上肢負載屈肘運動中，肱二頭肌EMGFT與肘關節切向力方差峰值存在一定相關性，且肘關節切向力方差峰值出現早于肱二頭肌EMGFT。
表1
評價流程輸出結果與Borg量表評分對照關系
Table1.
Comparison between the output results of the evalua tion process and Borg scale scores
圖10
EMGFT和肘關節切向力方差變化對照
Figure10.
Comparison of EMGFT and variance of tangential force at elbow joint
如表2所示，為15名受試者進行上肢負載屈肘運動時，疲勞分析評價流程與Borg主觀量表輸出結果對比，反映了受試者疲勞情況與主觀疲勞感受的關系。由表2可知，上肢屈肘運動中肘關節平均疲勞時間為17.10 s，肱二頭肌平均疲勞時間為19.92 s，且肘關節先于肱二頭肌出現疲勞狀態。平均時間差為1.92 s（過渡狀態時長平均值）。評價流程輸出結果與Borg量表評分區間結果一致，在肘關節疲勞前，Borg量表評分區間結果為輕松；在肘關節疲勞后肱二頭肌疲勞前，Borg量表評分區間結果為過渡；在肱二頭肌疲勞后，Borg量表評分區間結果為疲勞，驗證了本文提出模型的有效性。
表2
疲勞分析評價流程與Borg主觀量表輸出結果對比
Table2.
Comparison of fatigue analysis and evaluation process with Borg subjective scale output results
3.3 不同疲勞評價方法比較
基于單一sEMG信號特征值對疲勞進行評價的算法主要分為兩類：第一類是基于sEMG信號頻譜分析的疲勞評價算法[25]。例如： MPF、SMR。第二類是從sEMG信號的非線性特征對肌肉疲勞進行評價的算法。例如：fApEn、LZC。為驗證本文方法相較于單一sEMG信號特征值對疲勞狀態識別的優勢，采用試驗中得到的RMS值作為輸入值，將本文方法輸出結果與MPF?SMR?fApEn?LZC四種單一評價指標疲勞評價方法的試驗結果進行對比。
（1）MPF是以功率譜密度為基礎計算得出的特征量。大量研究成果表明，隨著疲勞狀態的推移，sEMG信號功率譜中的主頻率逐漸向較低頻率轉換，即MPF逐漸下降[26]。
（2）SMR是基于譜矩概念提出的一種分析肌肉疲勞的評價算法[27]，隨著疲勞的推移，SMR呈線性上升趨勢。
（3）fApEn是基于近似熵提出的一種改進肌肉疲勞監測算法[28]，隨著疲勞的推移，fApEn逐漸降低。相較于原始近似熵，fApEn呈現出更好的一致性及穩定性。
（4）LZC復雜度是從非線性分析角度來度量sEMG信號疲勞特征的指標[29-31]，LZC反映了一組信號時間序列隨時間增加出現新模式的速率，對于分析非線性的信號具有較好的效果。
為驗證本文方法優勢，將本文輸出結果分別與MPF、SMR、fApEn、LZC 四種疲勞評價算法進行比較，對比結果如圖11所示。
圖11
結果對比柱狀圖
Figure11.
Bar chart of results comparison
由圖11可知，本文提出的方法對總體疲勞狀態識別率結果達到了98.6%。相較于MPF、SMR、fApEn、LZC四種方法有明顯優勢，本文方法對過渡狀態的識別率達到100%，明顯優于其他方法。但在輕松狀態的識別上識別率略低于fApEn（98%），達到了97%。對疲勞狀態的識別率達到99%，與fApEn相同。本文方法對過渡及疲勞狀態的識別能保證較高敏感度和準確率，可以有效降低依賴單一sEMG信號特征值的疲勞評價的誤判率。本文方法識別率方差與fApEn相同，低于MPF和SMR且高于LZC，本文方法對于三種疲勞狀態的識別穩定性有待提高。
4 結論
在上肢運動過程中，疲勞是一類多維驅動的復雜生理現象。通過集成多源生物信號、生理信號的疲勞評價方法，能有效提高疲勞分析的準確率。本文針對目前單純依賴sEMG信號對疲勞狀態進行監測所導致的局限性問題，以上肢訓練中負載屈肘運動的疲勞狀態為分析對象，運用sEMG信號和動作捕捉技術手段，提出了一種融合改進的EMGFT算法與生物力學分析的疲勞分析方法，并將其輸出結果與Borg量表進行映射。通過對15名受試者進行上肢負載屈肘運動試驗，驗證了本文方法的有效性。本文方法對綜合疲勞狀態的識別率達到了98.6%，明顯優于其他基于單一sEMG信號特征值疲勞評價算法。試驗證明，本文提出的基于sEMG信號與動作捕捉的上肢疲勞度評價模型能為上肢運動及訓練提供有效的疲勞監護，同時為上肢運動輔助設備的效能評估和測試提供理論技術支持。
但本文方法一方面僅通過單一肌肉與單一關節的特征值對上肢負載屈肘運動的疲勞狀態進行分析，且在融合判斷的過程中沒有引入權重參考。因此本文方法對于上肢運動過程中的疲勞判斷的準確率有待提高；另一方面，本研究試驗受試者均為（22.4±1.5）歲的健康男性，因此在接下來的相關研究中，應將性別以及年齡等差異所導致的疲勞狀態差異納入研究考慮范圍。
重要聲明
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
作者貢獻聲明：徐兆負責論文及試驗思路構思及分析撰寫，呂健負責試驗設計的核實和論文修改，潘偉杰負責論文修改及排版，何愷倫負責試驗實施和數據整理。
倫理聲明：本研究通過了貴州大學人體醫學實驗倫理分委員會的審批（批文編號：HMEE-GZU-2021-T005）。
引言
人體上肢運動是骨骼肌在中樞神經的刺激下收縮并帶動骨骼以關節為支點所產生的定向運動[1]，持續性的負載上肢運動會使上肢疲勞產生一定程度的積累，導致肌肉骨骼勞損，造成肌肉骨骼失調[2]。上肢疲勞分析在上肢康復訓練、上肢體能訓練、上肢運動輔助設備的設計和優化等方面被廣泛運用[3]。有效的上肢疲勞監測能夠為上肢運動疲勞提供量化參考和評價，能將疲勞特征作為訓練策略制定和調整的重要依據，從而降低上肢勞損的風險。
目前，國內外各研究機構利用表面肌電信號（surface electromyography，sEMG）特征、生物力學分析等手段，在疲勞監測以及分析等方面取得了一定成果[4-7]。例如，Matsumoto等[4]首次提出了應用sEMG信號特征值對骨骼肌負荷進行量化判斷的方法，稱為肌電疲勞閾值（electromyogram fatigue threshold，EMGFT）算法。Kang等[5]對69名普通男性大學生進行負載遞增騎行試驗，分析了股外側肌的EMGFT，并對EMGFT和通氣閾值（ventilatory threshold，VT）進行了相關性分析，結果顯示EMGFT和VT均與疲勞具有高度的相關性。王琳等[6]結合生物力學理論分析了sEMG信號對頸部肌肉疲勞狀態的反映效果。劉曉光等[7]通過對10名健康受試者在負載可調的功率自行車上進行負載遞增騎行試驗，驗證了基于EMGFT來分析肌肉疲勞的有效性。然而，疲勞是一類多維驅動的復雜生理現象，但目前的研究中大多數將sEMG信號作為疲勞判斷的唯一標準，隨著上肢訓練時長的增加，生物電信號信噪比降低[8]，單純依賴sEMG信號對疲勞進行識別和分類會使結果不穩定，存在一定局限性。人的肢體活動是由骨骼肌以及關節共同協作完成的多自由度復雜運動[9-12]，sEMG信號主要用于量化所測肌肉的活動狀態，常用于身體局部肌肉負荷以及疲勞程度的評估，并不能對關節的疲勞程度進行評估。
針對以上問題，本文通過sEMG信號和動作捕捉技術，提出了一種融合EMGFT算法和生物力學分析的疲勞分析方法。通過對15名健康受試者進行右上肢負載曲肘試驗，同步采集肱二頭肌sEMG信號、右上肢動作捕捉數據，并通過柏格（Borg）疲勞度主觀自覺量表記錄受試者主觀疲勞感受，分析肱二頭肌的EMGFT，建立肘關節生物力學模型、計算肘關節切向力方差。本文通過綜合評價上肢骨骼肌以及關節疲勞程度，并與平均功率頻率（mean power frequency，MPF）、譜矩比（spectral moments ratio，SMR）、模糊近似熵（fuzzy approximate entropy，fApEn）、Lempel-Ziv復雜度（Lempel-Ziv complxity，LZC）四種單一評價指標疲勞評價方法進行對比，以驗證本文方法的有效性，本文研究方法如圖1所示。
圖1
基于sEMG信號和動作捕捉的上肢疲勞分析研究方法
Figure1.
Upper limb fatigue analysis method based on sEMG signal and motion capture
1 基于改進的EMGFT算法的肌肉疲勞分析方法
1.1 sEMG原始信號預處理
肌電信號（electromyography，EMG）是神經肌肉系統活動時的生物電信號變化經電極引導、放大、顯示和記錄所獲得的一維電壓時間序列信號，而sEMG信號則是淺層肌肉EMG信號和神經干上電活動在皮膚表面形成的具有非線性特征的綜合生物電效應[13]，sEMG信號通過電極在皮膚表面進行EMG信號的采集，具有非損傷性、實時性、多靶點測量等特點，常用于身體局部肌肉負荷以及疲勞程度的評估[14]。
在進行EMGFT的計算之前，需要對采集到的包含一定噪聲的原始sEMG信號進行預處理，首先設計0～500 Hz低通濾波器對經由sEMG信號采集模塊（EMGduino8通道肌電采集模塊，航奕生物，中國）采集到的原始sEMG信號進行低通濾波，去除高頻干擾；然后采用49.5～50.0 Hz自適應陷波器濾除信號中的工頻和諧波干擾；最后運用小波閾值去噪。由于小波變換中的時與頻的定位特性，可以進行信號的時變譜分析，能夠在各類細節上分析sEMG信號。因此小波閾值去噪是sEMG信號分析的有力工具[15]。本文通過觀察sEMG信號的頻域情況，發現能量在300 Hz以下的部分占全部能量的80%，主要集中在50～220 Hz之間。因此，本研究分解次數選擇4層，如圖2所示為sEMG信號預處理過程的時域圖。
圖2
sEMG信號預處理過程時域圖
Figure2.
Time domain diagram of the process of sEMG signal preprocessing
1.2 改進的EMGFT算法及其原理
EMGFT是指利用sEMG信號及其分析技術來確定運動肌疲勞閾值的方法。與傳統的反映疲勞閾值的生理指標相比，sEMG信號變化特征與肌肉活動狀態及其功能特征之間存在較好的關聯性[16-17]。目前，廣泛認同肌肉運動達到無氧閾值(anaerobic threshold，AT)時，骨骼肌進入疲勞狀態[18-21]。通過對EMGFT進行計算得到的結果，能夠用來測定AT和對應骨骼肌負荷強度。
劉曉光等[7]根據Matsumoto等[4]首次提出的EMGFT算法理論，進行固定時長的負載遞增騎行試驗，將sEMG信號的均方根（root mean square，RMS）值作為特征值，結合最小二乘法分析了股外側肌的EMGFT，驗證了其EMGFT算法的有效性。但是由于劉曉光等[7]試驗設計的特殊性，固定了試驗時長，試驗中采用了固定的時間窗長度（以符號tim表示）與移動窗長度（以符號mov表示），若受試者個體差異導致試驗時長不穩定，固定的tim與mov會造成一定程度的數據損失。為提高肌肉疲勞程度的識別率，本文針對以上問題，本文設置自適應tim與mov對EMGFT算法進行改進。
1.2.1 自適應計算時間窗和移動窗
由于受試者存在個體差異導致疲勞時間不穩定，即sEMG數據集中收集到的數據點數量不穩定。若固定tim與mov會使得時間窗與移動窗超出數據集邊緣，導致計算不完全，造成數據損失，影響試驗結果，因此本文依據具體試驗情況及需求，設置自適應窗計算RMS值。
為提升EMGFT算法疲勞狀態識別準確率，本文每組試驗需要通過自適應窗獲取RMS值50個，即mov需要移動50次，其中設置tim為mov的兩倍長度，則計算mov的計算公式如式(1)所示，而tim的計算公式如式(2)所示：
|
|
其中，n為sEMG數據集數據點數量，f為sEMG信號采集模塊采樣頻率。
1.2.2 改進的EMGFT算法流程
改進的EMGFT算法流程如圖3所示，具體算法步驟如下：
圖3
改進的EMGFT算法流程
Figure3.
Flowchart of improved EMGFT algorithm
步驟（1）：使用1.1小節中提到的sEMG信號預處理方法，對通過sEMG信號采集模塊采集到的右上肢肱二頭肌sEMG原始信號進行預處理，建立sEMG數據集。
步驟（2）：根據1.2.1小節中方法自適應計算tim、mov。
步驟（3）：由式(3)以及步驟（2）中計算得到的tim和mov，加窗計算sEMG數據集RMS值(以符號RMS表示)，建立RMS數據集，如式（3）所示：
|
式（3）中，Nt為tim內sEMG 信號數據點數量，i為sEMG信號數據點序號，Ei為sEMG數據集中第i個sEMG信號。
步驟（4）：將步驟（3）中得到的RMS數據集分為X、Y兩個部分，前5個RMS數據點為X1，剩余RMS數據點為Y1，X1和Y1組成擬合組1。再將前6個RMS數據點作為X2，剩余RMS數據點為Y2，X2和Y2組成擬合組2。直到Yn的RMS數據點數量等于X1的RMS數據點數量時，Xn和Yn組成擬合組n。分別對每組進行一階最小二乘擬合，每組得到兩條擬合直線，然后計算每組兩條擬合直線的斜率乘積。
步驟（5）：計算斜率乘積kXn·kYn最大的一組數據所擬合直線的交點，即為EMGFT。
2 融合改進的EMGFT算法和生物力學分析的疲勞分析流程
2.1 肘關節疲勞分析方法
上肢屈肘動作主要由肱二頭肌群以及肘關節協作完成[22]，僅利用sEMG信號的EMGFT算法只能對上肢屈肘運動中的肌肉疲勞程度進行評價。因此本文在改進的EMGFT算法的基礎上，通過動作捕捉技術采集三維人體肢體運動的位置數據，獲得各關節點絕對坐標，建立上肢球棍模型、計算肘關節空間位置關系，通過角速度與角加速度、切向力變化及其方差對關節疲勞程度進行評價，對照EMGFT算法輸出結果共同評價上肢負載屈肘過程中骨骼肌與關節的疲勞程度，建立基于sEMG信號和動作捕捉技術的上肢康復疲勞評價模型。
2.1.1 建立上肢生物力學模型
本文采用多剛體節段和鉸鏈關節的球棍結構建立人體上肢球棍模型，在球棍模型中對上肢關節角度進行定義。人體額狀面上肢球棍模型如圖4所示。
圖4
人體額狀面上肢球棍模型
Figure4.
Ball and stick model of the upper limb on the frontal plane of the human body
圖4中，點a為肩關節中心點，點b為肘關節中心點，點c為腕關節中心點。c1為屈肘運動過程中腕關節某一點位置。上臂與前臂的夾角θ為肘關節角度。屈肘運動可以近似理解為前臂繞肘關節點（點b）作非勻速運動。
根據生物力學分析以及人體測量學中肢體參數計算公式，對上肢屈肘運動過程進行受力分析，上肢屈肘運動過程受力分析如圖5所示。
圖5
上肢屈肘運動過程受力分析
Figure5.
Stress analysis of elbow flexion in upper extremity
圖5中，P為前臂和手部質心，H為肘關節中心，L1為前臂和手的質心到肘關節中心點的距離，θ為上臂與前臂的夾角，G1為前臂和手部的重力，G2是重物的重力，Fg是手部對重物的支持力，FX為前臂和手繞肘關節點作非勻速運動時的離心力，Fn為前臂和手繞肘關節點作非勻速運動時的切向力。
2.1.2 肘關節疲勞分析方法
為評價上肢負載屈肘運動中的肘關節疲勞程度，本文運用上肢位置信息絕對坐標，計算肘關節角度及其變化情況，計算肘關節切向力方差，以此反映肘關節切向力的穩定性。肘關節疲勞分析方法具體步驟如下：
步驟（1）：三維人體肢體運動的位置數據獲取
根據人體額狀面上肢球棍模型，如圖4所示，于a、ｂ、ｃ三點設置動作捕捉慣性傳感器，記錄a、ｂ、ｃ三點運動過程中的絕對坐標數據。
步驟（2）：肘關節角度計算
根據動作捕捉慣性傳感器獲得的絕對坐標數據，定義 代表受試者m身上第n個點在t時刻的三維坐標值。如圖4所示，a、b、c是空間中3個點，θ代表肘關節角度，θ計算方法如式（4）所示:
|
步驟（3）：肘關節角速度計算
由于上肢屈肘運動中肘關節角度θ隨時間t非勻速變化，根據本小節步驟（2）中計算得到的肘關節角度θ加窗計算角速度，計算方法如式(5)所示：
|
步驟（4）：肘關節角加速度計算
由于人體上肢屈肘運動并不是完全的勻速運動，因此上肢屈肘運動過程中角速度處于變化中，根據步驟（3）中得到的角速度加窗計算角加速度a，計算方法如式（6）所示：
|
步驟（5）：肘關節切向力計算
肘關節切向力反映肘關節負載運動過程中的受力情況。如圖5所示，人體上肢屈肘運動時肘關節切向力Fn處于不規則變化中，根據身體各部重量占體重比例[9]計算前臂和手部的重力G1，G2為重物的重力，由此計算肘關節負載m，計算方法如式（7）所示，通過運動學原理以及本小節步驟（4）中計算得到的角加速度a，計算切向力Fn，計算方法如式（8）所示：
|
|
步驟（6）：肘關節切向力方差計算
方差是衡量一組隨機變量離散程度的度量，常用來衡量一批數據的波動大小，反映數據的穩定性[23-24]，由于上肢肘關節處于疲勞狀態時肘關節受力會出現波動，為評價肘關節受力穩定性，根據本小節步驟（5）中得到的肘關節切向力Fn，加窗計算Fn方差s2，計算方法如式（9）所示：
|
式(9)中，n為肘關節切向力Fn數量，為肘關節切向力Fn均值。
2.2 評價流程
由于單一sEMG信號無法作為評價關節疲勞的指標，本文融合改進的EMGFT算法和生物力學分析，對上肢負載屈肘運動進行疲勞分析。在上肢屈肘運動過程中，綜合評價肱二頭肌及肘關節的疲勞程度，提出融合EMGFT算法和生物力學分析的疲勞評價流程，在時間序列中，定義達到EMGFT閾值時間點為肌肉進入疲勞狀態，達到關節切向力方差峰值時間點為關節進入疲勞狀態。具體步驟如圖6所示。該評價流程輸出四類結果分別為：P1（肌肉關節均疲勞）、P2（僅肌肉疲勞）、P3（僅關節疲勞）、P4（肌肉關節均未疲勞）。
圖6
融合EMGFT算法和生物力學分析的疲勞分析方法
Figure6.
Fatigue analysis method integrating EMGFT and biomechanics
3 試驗及結果分析
3.1 試驗
為驗證本文所提出的融合改進的EMGFT算法和生物力學分析的疲勞分析方法的有效性，本試驗招募15名健康男性志愿者參與研究，年齡（22.4 ± 1.5）周歲，身高（173.2 ± 5.8）cm，體重（65.3 ± 7.2）kg，體質指數（body mass index，BMI）為22.7 ± 1.8。受試者均身體健康，無肢體損傷等疾病和不良嗜好。數據采集地為現代制造技術教育部重點實驗室，受試者在試驗前均已熟悉本試驗方案和流程且簽署志愿者知情同意書。整個試驗已通過貴州大學人體醫學實驗倫理分委員會審查并獲得倫理審查通知書。受試者在試驗前3月內每周鍛煉次數不超過2次且在試驗前48 h內未進行劇烈運動。sEMG信號采集設備采用sEMG信號采集模塊（EMGduino8通道肌電采集模塊，航奕生物，中國）。動作捕捉設備采用Smartsuit（pro，Rokoko，丹麥）。根據人體生理學結構，同步采集每位受試者右上肢肱二頭肌sEMG信號和右上肢動作捕捉數據，sEMG信號采樣頻率為1 kHz，動作捕捉數據采樣頻率為30 Hz，傳感器位置如圖7所示，S1為右上肢肱二頭肌處的sEMG信號傳感器，S2、S3為sEMG信號參考電極。C1、C2、C3、C4分別為右上肢肩、上臂、前臂以及手處的9自由度慣性測量單元（inertial measurement unit，IMU）運動跟蹤傳感器（Smartsuit pro，Rokoko，丹麥）。
圖7
上肢試驗過程及傳感器放置
Figure7.
Upper limb experiment procedure and sensor placement
受試者測試前刮除右上肢肱二頭肌周圍體毛，用75%酒精擦拭傳感器處皮膚，降低汗液導致電極移動引起的低頻噪聲。受試者保持坐立，右上肢進行負載為7.5 kg的屈肘運動，定義受試者上肢由圖7位置a進行負載屈肘運動到圖7位置b再回到圖7位置a為一次完整屈肘運動。每位受試者進行15次完整的屈肘運動，并在負載屈肘運動期間在測試者幫助下完成Borg量表評分。在Borg量表中將受試者主觀疲勞程度分為1～12評分區間，分別選取Borg量表中 1～4、5～8、9～12 分數段對應輕松、過渡和疲勞 3 種狀態進行相關問題的分析和討論。
3.2 試驗結果分析
利用1.2節所述改進的EMGFT算法分析受試者右上肢肱二頭肌的疲勞時間，如圖8所示為隨機抽取的其中一名受試者RMS值分組擬合過程，X-RMS為Xn擬合直線，Y-RMS為Yn擬合直線。選取所擬合直線斜率乘積最大的一組，獲取該受試者在進行上肢負載屈肘運動時肱二頭肌的EMGFT。
圖8
RMS一階最小二乘擬合過程
Figure8.
RMS first order least squares fitting process
利用2.1.2節中所述肘關節疲勞分析方法，計算分析受試者右上肢肘關節疲勞情況。其中一名受試者進行上肢負載屈肘運動時，手與前臂繞肘關節中心點作非勻速運動時的肘關節角度變化情況，以及該過程中切向力變化情況如圖9所示。
圖9
肘關節角度和切向力變化情況
Figure9.
Changes in elbow angle and tangential force
表1為融合改進的EMGFT算法和生物力學分析的疲勞分析方法輸出結果與Borg量表評分區間對照。如圖10所示，為受試者中隨機3名受試者進行上肢負載屈肘運動時，EMGFT和肘關節切向力方差變化對照。由圖10可知，3名受試者先后出現肱二頭肌EMGFT以及肘關節切向力方差峰值。上肢負載屈肘運動中，肱二頭肌EMGFT與肘關節切向力方差峰值存在一定相關性，且肘關節切向力方差峰值出現早于肱二頭肌EMGFT。
表1
評價流程輸出結果與Borg量表評分對照關系
Table1.
Comparison between the output results of the evalua tion process and Borg scale scores
圖10
EMGFT和肘關節切向力方差變化對照
Figure10.
Comparison of EMGFT and variance of tangential force at elbow joint
如表2所示，為15名受試者進行上肢負載屈肘運動時，疲勞分析評價流程與Borg主觀量表輸出結果對比，反映了受試者疲勞情況與主觀疲勞感受的關系。由表2可知，上肢屈肘運動中肘關節平均疲勞時間為17.10 s，肱二頭肌平均疲勞時間為19.92 s，且肘關節先于肱二頭肌出現疲勞狀態。平均時間差為1.92 s（過渡狀態時長平均值）。評價流程輸出結果與Borg量表評分區間結果一致，在肘關節疲勞前，Borg量表評分區間結果為輕松；在肘關節疲勞后肱二頭肌疲勞前，Borg量表評分區間結果為過渡；在肱二頭肌疲勞后，Borg量表評分區間結果為疲勞，驗證了本文提出模型的有效性。
表2
疲勞分析評價流程與Borg主觀量表輸出結果對比
Table2.
Comparison of fatigue analysis and evaluation process with Borg subjective scale output results
3.3 不同疲勞評價方法比較
基于單一sEMG信號特征值對疲勞進行評價的算法主要分為兩類：第一類是基于sEMG信號頻譜分析的疲勞評價算法[25]。例如： MPF、SMR。第二類是從sEMG信號的非線性特征對肌肉疲勞進行評價的算法。例如：fApEn、LZC。為驗證本文方法相較于單一sEMG信號特征值對疲勞狀態識別的優勢，采用試驗中得到的RMS值作為輸入值，將本文方法輸出結果與MPF?SMR?fApEn?LZC四種單一評價指標疲勞評價方法的試驗結果進行對比。
（1）MPF是以功率譜密度為基礎計算得出的特征量。大量研究成果表明，隨著疲勞狀態的推移，sEMG信號功率譜中的主頻率逐漸向較低頻率轉換，即MPF逐漸下降[26]。
（2）SMR是基于譜矩概念提出的一種分析肌肉疲勞的評價算法[27]，隨著疲勞的推移，SMR呈線性上升趨勢。
（3）fApEn是基于近似熵提出的一種改進肌肉疲勞監測算法[28]，隨著疲勞的推移，fApEn逐漸降低。相較于原始近似熵，fApEn呈現出更好的一致性及穩定性。
（4）LZC復雜度是從非線性分析角度來度量sEMG信號疲勞特征的指標[29-31]，LZC反映了一組信號時間序列隨時間增加出現新模式的速率，對于分析非線性的信號具有較好的效果。
為驗證本文方法優勢，將本文輸出結果分別與MPF、SMR、fApEn、LZC 四種疲勞評價算法進行比較，對比結果如圖11所示。
圖11
結果對比柱狀圖
Figure11.
Bar chart of results comparison
由圖11可知，本文提出的方法對總體疲勞狀態識別率結果達到了98.6%。相較于MPF、SMR、fApEn、LZC四種方法有明顯優勢，本文方法對過渡狀態的識別率達到100%，明顯優于其他方法。但在輕松狀態的識別上識別率略低于fApEn（98%），達到了97%。對疲勞狀態的識別率達到99%，與fApEn相同。本文方法對過渡及疲勞狀態的識別能保證較高敏感度和準確率，可以有效降低依賴單一sEMG信號特征值的疲勞評價的誤判率。本文方法識別率方差與fApEn相同，低于MPF和SMR且高于LZC，本文方法對于三種疲勞狀態的識別穩定性有待提高。
4 結論
在上肢運動過程中，疲勞是一類多維驅動的復雜生理現象。通過集成多源生物信號、生理信號的疲勞評價方法，能有效提高疲勞分析的準確率。本文針對目前單純依賴sEMG信號對疲勞狀態進行監測所導致的局限性問題，以上肢訓練中負載屈肘運動的疲勞狀態為分析對象，運用sEMG信號和動作捕捉技術手段，提出了一種融合改進的EMGFT算法與生物力學分析的疲勞分析方法，并將其輸出結果與Borg量表進行映射。通過對15名受試者進行上肢負載屈肘運動試驗，驗證了本文方法的有效性。本文方法對綜合疲勞狀態的識別率達到了98.6%，明顯優于其他基于單一sEMG信號特征值疲勞評價算法。試驗證明，本文提出的基于sEMG信號與動作捕捉的上肢疲勞度評價模型能為上肢運動及訓練提供有效的疲勞監護，同時為上肢運動輔助設備的效能評估和測試提供理論技術支持。
但本文方法一方面僅通過單一肌肉與單一關節的特征值對上肢負載屈肘運動的疲勞狀態進行分析，且在融合判斷的過程中沒有引入權重參考。因此本文方法對于上肢運動過程中的疲勞判斷的準確率有待提高；另一方面，本研究試驗受試者均為（22.4±1.5）歲的健康男性，因此在接下來的相關研究中，應將性別以及年齡等差異所導致的疲勞狀態差異納入研究考慮范圍。
重要聲明
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
作者貢獻聲明：徐兆負責論文及試驗思路構思及分析撰寫，呂健負責試驗設計的核實和論文修改，潘偉杰負責論文修改及排版，何愷倫負責試驗實施和數據整理。
倫理聲明：本研究通過了貴州大學人體醫學實驗倫理分委員會的審批（批文編號：HMEE-GZU-2021-T005）。
