利用肌肉協同研究上肢簡單和復雜運動的關聯性

黃永欣,文斌,毛菁菁,姚瑞,尹孟奇,徐進

(1. 西安交通大學生命科學與技術學院健康與康復科學研究所,710049,西安;2. 西安交通大學生物醫學信息工程教育部重點實驗室,710049,西安;3. 西安交通大學第一附屬醫院康復科,710061,西安)

肢體運動是人體重要的生理功能之一,中樞神經系統(central nervous system,CNS)通過控制相關肌肉運動的協同,以完成日常生活中的各種運動,由于運動模式多樣,使得運動控制問題變得極為復雜[1-2]。

在人體運動控制研究中,肌肉協同學說被廣泛用于研究神經肌肉控制問題,即CNS將相應部位的多塊肌肉以不同強度的協同作用結合在一起,完成運動控制的任務[3-4]。CNS通過肌肉協同機制實現對運動控制的簡化,協同被認為是神經運動模塊化控制的基本單元[5-6]。因此,研究運動過程中相關肌肉的協同作用對于神經肌肉控制模式研究具有十分重要的意義。

表面肌電信號(surface electromyography, sEMG)是經皮膚表面電極采集的生物電信號,其蘊含豐富的神經肌肉活動和運動意圖等信息,同時具有提取方便、準確和無創的優點,已經成為分析肢體運動過程中肌肉協同的有利手段[3,7]。

從多通道sEMG信號中獲取肌肉協同信息屬于盲源分離的過程,常見的矩陣分解算法有主成分分析(principal components analysis, PCA)、獨立成分分析(independent component analysis, ICA)和非負矩陣分解(non-negative matrix factorization, NMF)[6,8]等。因為NMF的結果具有非負約束和魯棒性強等優點,使得基于NMF的肌肉協同分析更具優勢[8]。

肌肉協同分析作為一種探索神經肌肉控制模式的有效方法,已經被廣泛用于臨床研究和機器人控制等領域。該方法可用于探索神經功能障礙患者的肢體運動控制模式,進而為康復治療提供指導。Steele等對比腦癱患者和健康人在步態活動中提取的肌肉協同模式,發現腦癱患者在步態活動中的肌肉變化用較少的肌肉協同即可表示[9]。Allen等通過比較帕金森患者在短期康復治療前后的步行與平衡的肌肉協同模式,發現了步行與平衡存在相似的協同模式與特有的協同模式,且康復治療后相似的協同模式數量增加[10]。

肌肉協同分析可以提取運動意圖,從而實現運動控制和識別。桂奇政等利用肌肉協同理論和支持向量回歸構建模型,實現了從sEMG信號到關節運動的連續估計[11]。鄭楠等利用肌肉協同和最小二乘法構建手勢識別算法,從神經協同的角度提取運動意圖,減少訓練數據的同時提高用戶無關手勢識別的正確率[12]。

研究特定肢體運動下的肌肉協同模式和不同肢體運動間肌肉協同模式的相似性以及差異,將能夠更好地了解、掌握肢體運動中的神經肌肉控制模式。Saito等研究健康人在不同速度和坡度的跑步機上跑步的肌肉協同差異,發現不同條件下的跑步均可由4個肌肉協同來解釋并存在相似的肌肉協同模式[13]。

與下肢運動相比,上肢運動模式多樣且復雜,運動控制問題更為復雜。上肢運動的肌肉協同分析開展較少,且已有研究大多是針對簡單運動。謝平等在研究健康人的腕屈伸運動的肌肉協同中發現了與腕部運動相關肌群的協同關系、共享的基本肌肉協同模式以及同一肌肉協同模式中協同性較高的肌肉間耦合關系較強[14]。Tang等通過對健康人的肘關節和肩關節的3個相似運動的肌肉協同分析,發現了3個相似運動中存在共享的基本肌肉協同模式[15]。王洪安等在伸手運動的肌間神經耦合關系研究中發現,卒中患者的肌肉協同數下降、β頻段內肌間耦合強度具有顯著差異以及γ頻段內大部分協同肌肉的耦合關系較強[16]。

此外,協同模式的關聯性研究較少。Cheung等在卒中患者的上肢健側運動與患側運動的肌肉協同模式關聯性研究中,通過構建協同模式的線性模型并利用非負最小二乘法求解模型系數,發現了肌肉協同模式的保存、合并和分離[17]。Hashiguchi等利用相同方法研究亞急性期卒中患者在康復過程中的步態活動,發現肌肉協同模式的合并可以作為亞急性期卒中患者運動協調異常的重要標志[18]。

迄今,復雜運動的肌肉協同模式研究較少,上肢簡單運動和復雜運動的肌肉協同模式之間的關聯性研究缺乏,給上肢運動的神經肌肉控制機制的全面掌握帶來困難。

因此,本文分別設計了上肢簡單和復雜運動范式,采集不同運動執行過程中的sEMG信號,通過肌肉協同分析探究簡單運動和復雜運動執行過程中肌肉協同模式之間的關聯性,從而為肢體運動的神經肌肉控制模式和肌-機接口技術研究提供參考。

1 實驗設計和信號采集

1.1 實驗設計

選擇日常生活中常見的上肢動作設計簡單和復雜運動范式。簡單運動范式包括:手伸展(HE)和手屈曲(HF)、腕屈曲(WF)和伸展(WE)、肘屈曲(EF)和伸展(EE)、肩前屈(SF)和外展(SA);復雜運動包括兩個,復雜運動Ⅰ(ComplexF)主要由肩前屈、肘屈伸和手屈曲組成,復雜運動Ⅱ(ComplexA)主要由肩外展、肘屈伸和手屈曲組成,如圖1、圖2所示。

(a)手伸展

(a)復雜運動Ⅰ

實驗利用Matlab軟件和Psychtoolbox工具箱呈現運動指示圖片以規范數據采集流程,如圖3所示。

首先在屏幕上出現一個“+”注視點,持續2 s;當運動指示圖片出現后,根據圖片提示做相應的動作并保持到圖片結束,簡單運動圖片持續10 s,復雜運動圖片持續7 s;間隔休息30 s。每個動作完成之后,必須恢復初始位置狀態再進行下一個動作。每種動作會重復3次,然后換下一種動作,每種動作結束后,間隔30 s用于休息和切換采集姿態,再繼續進行實驗。

另外,采集尺側腕屈肌(flexor carpi ulnaris,FCU)、指淺屈肌(superficial flexor digitorum,SFD)、橈側腕伸肌(extensor carpi radialis,ECR)、指伸肌(extensor digitorum,ED)、肱二頭肌(biceps brachii,BB)、肱三頭肌(triceps brachii,TB)、三角肌前束(anterior bundle of deltoid,ABD)和三角肌中束(middle bundle of deltoid,MBD)8個與運動相關肌肉的最大自主收縮(maximal voluntary contraction, MVC),用于后續肌電信號的MVC歸一化,MVC動作模式如圖4所示。

(a)尺側腕屈肌

1.2 表面肌電信號采集

使用美國Delsys公司的無線肌電信號采集系統進行8通道肌電信號采集,采樣頻率為1 259 Hz。傳感器順著肌肉纖維的方向,貼附于肌腹隆起處。同步采集實驗過程中右側上肢的8通道sEMG信號,肌電電極位置及相關肌肉如圖5所示。

圖5 上肢運動相關肌肉與肌電采集位置Fig.5 Related muscles of the movements of upper limbs and electrode position of sEMG signal

來自西安交通大學的20名健康青年參加了實驗,實驗參與者均無上肢肌肉和中樞神經等疾病,受試前24 h未進行劇烈運動,均為右利手,男女比例為11∶9,年齡為25.1±3.093(均值±標準差)。

2 基于非負矩陣分解的肌肉協同分析

2.1 信號預處理

為了更好地分析數據,首先,提取sEMG原始信號進行預處理并提取信號包絡[19-21],具體步驟如下:①通過截止頻率為20 Hz的4階巴特沃斯高通濾波器,用于去除信號的運動偽跡和基線漂移;②進行全波整流,使得信號滿足非負要求,符合生理上的肌肉激活度;③信號通過截止頻率為5 Hz的4階巴特沃斯低通濾波器去除高頻噪聲,完成信號包絡的提取。

其次,為了避免肌電幅值個體差異對肌肉協同分析的影響,對sEMG信號進行MVC歸一化處理[6,19]。MVC歸一化反映了該位置肌肉的激活程度,且保留了各通道之間的比例關系[22-23]。

MVC歸一化的具體步驟為:對被試者MVC動作時的sEMG信號做上述的提取包絡處理,然后計算滑動窗內最大平均值,窗長設置為250 ms。

最后,利用Teager-Kaiser能量算子[24]自動檢測運動起始點和終止點,提取sEMG信號片段,降采樣到200 Hz,并將處理后的片段進行串聯。

2.2 肌肉協同提取

1999年,Lee和Seung提出了NMF算法,該算法在非負的約束下可以實現非線性的降維[25]。因為NMF的分解結果具有非負約束性,其分解的結果更具實際意義[26]。因此,本文選擇NMF算法進行肌肉協同分析。

對于預處理后的sEMG信號X,肌肉協同提取算法需要找到r個肌肉協同及其對應的線性加權系數,使之能較好地重構信號X,即

Xm×n=Wm×rHr×n+Em×n

(1)

式中:W是肌肉協同矩陣;H是肌肉協同的加權系數矩陣;E是誤差矩陣;r為肌肉協同數;m是肌電導聯數;n為肌電信號長度。

NMF是一個最值優化問題,在肌肉協同研究中的求解算法主要有基于歐式距離的乘性迭代算法、基于交替非負最小二乘法的算法和基于梯度下降的加性算法[27]等?；跉W氏距離的乘性迭代算法具有較優的分解效果,但收斂性不穩定,可能陷入局部最優;基于梯度下降的加性算法雖然簡單但收斂速度慢,迭代過程易受迭代參數影響;基于交替非負最小二乘法的算法即固定一個矩陣對另一個矩陣進行優化,不斷交替進行,可收斂至穩定點[28]。

因此,為了避免陷入局部最優,本文使用了由兩步驟組成的NMF算法,具體如下。

(1)隨機選擇W和H,并重復使用20次基于歐式距離的乘性迭代算法運行NMF(假設噪聲符合高斯分布),選擇初步收斂最好的結果作為矩陣分解結果。

(2)使用初步收斂結果最好的W和H作為第二輪NMF分解的初始值,使用基于交替非負最小二乘法的算法運行NMF。

2.3 肌肉協同數確定

矩陣X經NMF后得到W和H矩陣,但NMF算法需要確定肌肉協同數r。通常計算原始矩陣X和重構矩陣X′=WH的變異性占比(RVAF)[9,16,29-30]來確定協同數r,計算公式為

(2)

RVAF越大,表示重構矩陣X′與原始矩陣X越接近。計算不同協同數r下的所有被試者平均的RVAF,具體結果如圖6所示,圖6中誤差線為標準差,Th為閾值。

圖6 不同動作下r與RVAF之間的關系Fig.6 The relationship between r and RVAF under different movements

以往的研究確定肌肉協同數時,RVAF閾值常設在0.80～0.95之間[21]。因此,為了更好地重構信號并保留原始信息,本文設定當滿足RVAF的均值大于等于95%,且隨著協同數r的增加RVAF的增量小于等于1%時[26],所分解的肌肉協同數r作為最佳分解數。滿足上述條件,得到簡單、復雜運動的最佳分解數和模式,見表1、表2,以手伸展分解出的協同模式HE-1為例進行解釋,其中HE表示手伸展,1為編號,其他類同。

表1 簡單運動的最佳分解數和協同模式

表2 復雜運動的最佳分解數和協同模式

2.4 協同模式排序

對同一動作下不同被試者提取到了類似的協同矩陣W,但提取的協同模式Wi,j以隨機順序出現,因此需要對協同模式根據相似度進行排序,具體步驟如下。

(1)選取被試者k的協同模式作為初始模板,其順序為Wk,1,Wk,2,…,Wk,r。

(2)計算每個被試者提取的協同模式與初始模板的余弦相似度[30],計算公式為

(3)

式中:Wi,j為第i個被試者提取的第j個協同模式;Wk,r為初始模板的第r個協同模式;cosθi,k,j,r為Wi,j與Wk,r的余弦相似度。

(4)從相似度矩陣Ri,k找到最大值maxj,r,即第i個被試者提取的協同模式Wi,j的順序位置為初始模板中的協同模式Wk,r的位置,并將相似度矩陣對應的第j行和第r列置為0,重復步驟(4),依次確定其他協同模式順序。

(5)計算所有被試者排序后的協同模式與該初始模板的平均相似度,選擇其他被試者作為初始模板并重復步驟(2)～(5)。

(6)找出平均相似度最高的初始模板作為最終排序模板。

2.5 關聯性分析

為了研究簡單、復雜運動肌肉協同模式的關聯性,首先構建線性組合模型如下式

(4)

其次,利用非負最小二乘法[17,31]求解系數mi,k,并對求解后的系數mi,k進行標準化,并選取0.2作為貢獻度閾值[17,31]。

最后,選取大于0.2閾值的非負線性系數及其對應的簡單運動協同模式構建重構運動模式,并通過余弦相似度量化重構運動模式與復雜運動協同模式之間的相似性,以探討簡單運動和復雜運動的關聯性。

3 結果及討論

3.1 簡單運動與復雜運動的協同模式

基于排序后的協同模式,將所有被試者在簡單運動和復雜運動協同模式下相關肌肉的權值系數進行平均,并對簡單運動和復雜運動的協同模式進行分析。

3.1.1 簡單運動分析

(1)手部動作。手部動作協同模式下相關肌肉的權值系數(均值±標準差),如圖7所示。手部動作分解的協同主要由FCU、SFD、ECR、ED 4塊肌肉構成。手伸展動作主要分解為一個協同,模式HE-1以ED為主,ED為前臂后群肌,其可以使得后四指掌指關節伸展,同時HE-1中顯示了FCU、SFD、ECR和ED 4塊肌肉的協同配合。手屈曲動作分解為兩個協同,HF-1以SFD為主,SFD為前臂前群肌,負責第2～5指的掌指關節和近側指骨間關節的屈曲;HF-2中FCU激活度最高,表明FCU對手屈曲也有一定作用。簡單手部動作的協同模式中均顯示了4塊肌肉的協同關系,表明主要屈肌和主要伸肌會相互配合。

(a)手伸展

(2)腕部動作。腕部動作協同模式下相關肌肉的權值系數(均值±標準差),如圖8所示。腕部動作分解的協同也主要涉及FCU、SFD、ECR、ED 4塊肌肉。腕伸展動作主要分解為兩個協同,模式WE-1以ED為主,體現該模式與手伸展運動有關,表明ED除了參與四指掌指關節的伸展,還參與腕關節的伸展運動;模式WE-2以ECR和ED為主,ECR為前臂后群肌,其主要功能就是伸展腕關節。腕屈曲動作分解為兩個協同,模式WF-1主要由SFD、FCU和ED貢獻,尤其是SFD,體現該模式與手屈曲運動有關,表明SFD雖然主要負責指關節屈曲,但是也會支持腕部屈曲;模式WF-2以FCU為主,FCU為前臂前群肌,是腕關節屈曲的主要肌肉之一。

(a)腕伸展

(3)肘部動作。肘部動作協同模式下相關肌肉的權值系數(均值±標準差),如圖9所示。

(a)肘伸展

肘伸展動作分解為兩個協同,模式EE-1以TB為主導,TB屬于上臂后側肌群,主要負責伸肘;模式EE-2中ED存在較大的肌肉激活度,表明ED對肘伸展也有一定作用。在肘屈曲動作中,模式EF-1中以ECR、FCU、ED和BB為主,表明4塊肌肉協同配合完成曲肘;模式EF-2以BB為主,BB位于上臂前側,當BB收縮時會使肘關節屈曲。

(4)肩部動作。肩部動作協同模式下相關肌肉的權值系數(均值±標準差),如圖10所示。在肩前屈動作的協同模式SF-2中,體現了ABD在肩前屈運動中起到的重要作用。與此類似,在模式SA-2中,表明MBD在肩外展動作起到了主要作用,同時TB和ABD協同配合進行肩部外展。此外,在簡單肩部動作的協同模式SF-1和SA-1中表現出MBD、ABD和TB等肌肉之間的協同關系。

(a)肩前屈

3.1.2 復雜運動分析

復雜運動協同模式下相關肌肉的權值系數(均值±標準差),如圖11所示。復雜運動Ⅰ分解出3個協同模式,模式ComplexF-1中FCU、SFD、ECR和ED有較大貢獻,體現了4塊肌肉之間的協同配合,反映了該范式主要與手部動作或腕部動作有關;模式ComplexF-2體現了FCU、ECR、ED和ABD的主要作用,表明與手部動作、腕部動作和肩部動作有關;模式ComplexF-3中ABD和MBD貢獻較大,尤其是ABD,說明該模式與簡單運動的肩前屈有關。

復雜運動Ⅱ分解出3個協同模式,模式ComplexA-1中FCU、SFD、ECR和ED有較大貢獻,尤其是FCU,表明與手屈曲、腕屈曲動作有關;模式ComplexA-2和模式ComplexA-3體現了MBD、ED、FCU和SFD等肌肉的主要作用,尤其是MBD,說明這兩個模式與手部動作、腕部動作和肩外展動作有關。

3.2 關聯性研究

通過線性擬合,并對系數進行標準化,線性組合的關系如圖12所示。以復雜運動模式ComplexF-1為例,其可以由非零系數(0.19、0.06、0.51、0.05、0.04、0.04和0.11)及對應的簡單運動模式(HF-2、WE-2、WF-1、EE-1、EF-2、SF-2和SA-2)構成。同理,可以得到其他復雜運動協同模式的構成。之后,選取大于0.2的系數及其對應的簡單運動協同模式,用于后續構建重構運動模式。

紅色方框中數據為大于0.2閾值的非負線性系數圖12 肌肉協同模式的線性組合關系Fig.12 Linear combination of muscle synergy patterns

在以往協同模式的相似性研究中,余弦相似度大于0.75定性為相似,大于0.9定性為非常相似[17,30]。

通過復雜運動協同模式與重構運動模式之間的相關性分析,發現了協同模式的保留和合并,如圖13所示。

(a)復雜運動Ⅰ(Fi,i=1,…,5,Fi為大于0.2閾值的非負線性系數)

在復雜運動Ⅰ中,模式ComplexF-1與對應的重構運動模式的余弦相似度為0.899,具有較高的相似度,同時,該重構運動模式僅由WF-1與非線性系數構成。因此,模式ComplexF-1可以由WF-1構成,體現了協同模式的保留;模式ComplexF-2與對應的重構運動模式非常相似,兩者的余弦相似度為0.984,因此模式ComplexF-2可以通過組合HF-2、WE-1和SF-2來重構,體現了復雜運動協同模式為簡單運動協同模式的合并;同樣,模式ComplexF-3與對應的重構運動模式的相似度為0.957,說明模式ComplexF-3可以由SF-2來重構,體現了復雜運動協同模式為簡單運動協同模式的保留。

重構運動模式與復雜運動Ⅱ的肌肉協同模式也具有很高的相似度(在0.909～0.995之間)。在復雜運動Ⅱ中,模式 ComplexA-1主要由HF-2和WF-1組成,表明與手部屈曲和腕部屈曲運動有關;模式ComplexA-2由HF-2、WE-1和SA-2組成,體現了協同模式的合并;模式 ComplexA-3主要由HF-2、WF-1和SA-2組成,表明與手部屈曲、腕部屈曲和肩部外展運動的組合有關,體現了復雜運動協同模式為簡單運動協同模式的線性組合。

4 結 論

本文基于表面肌電信號,通過采用非負張量分解和余弦相似度分析等方法,分別對簡單運動和復雜運動執行過程中的肌肉協同狀況進行比較分析。研究發現復雜運動中存在明顯的簡單運動肌肉協同模式的保留和合并過程。從肌肉協同的角度證明了人體復雜運動模式是從簡單運動提取的基本肌肉協同模式的基礎上線性組合而來,復雜運動與簡單運動之間還存在共享的基本肌肉協同模式。本文研究結果有助于進一步增強對神經肌肉運動控制機制的認識,為肌-機控制技術相關領域研究提供理論基礎。