下肢表面肌電信號的降維和映射分析*

章亞男,景銀平,沈林勇,宋 薇,錢晉武(上海大學機電工程與自動化學院,上海 200072)

在人體醫療、救災和軍工領域的機器人研究中,主要有康復機器人、外骨骼機器人等,通過采集人體信號,使機器人獲取和預測人類的運動意圖或行為,來完成醫生所要求的康復訓練和動作或是加強人體的力量減少人體體力消耗等。尤其是在醫療研究中,主動性的康復訓練對患者治療有更好的效果,基于肌電信號的外骨骼康復機器人,是利用人體表面肌電信號來控制機器人幫助偏癱患者行走并達到康復訓練的有效手段[1-2]。

表面肌電信號sEMG(surface Electromyography)是神經肌肉系統活動時從皮膚表面記錄下來的生理電信號,可以通過表面電極從被檢測肌肉的皮膚表面引導而獲得。由肌電信號產生的特性分析可知,肌電信號總是先于實際動作出現,因此通過分析肌電信號來預測判斷人體的動作或是人體動作的分類是可行的[3]。

肌電信號的差異性、相似性和重復性等特征決定了sEMG可以作為外骨骼機器人等的控制信號[4]。目前,國內外對于sEMG的特征提取和模式識別的研究分析方法多樣,例如Canal等[5]采用短時傅里葉變換和小波變換對sEMG進行特征提取,得到了較為平穩的肌電信號,但是其特征提取性能受到所選信號長度的限制。Michieletto等[6-7]基于高斯混合模型,通過EMG訓練來估計單個關節角度,在計算精度和時間上都有較高的表現。Veer等[8-9]利用后向傳播神經網絡算法分析人體sEMG,實驗分類結果表明準確率達90%以上。陳彥釗等[10-11]提出了多步支持向量機分類策略,實現了較高的上肢動作分類準確率,使得整個上肢分類過程得以簡化。

總體而言,國內外對肌電信號的研究大都集中于特征提取和動作模式分類識別方面,對于把肌電信號作為康復機器人的輸入控制信號,實時預測人體行走運動信息的研究較少,本文主要研究人體行走時的髖關節、膝關節和踝關節角度和運動相關的肌肉sEMG的映射關系,對采集得到的sEMG進行信號處理,并用神經網絡預測關節角度,實現sEMG和關節角度的實時映射,主要應用于如圖1所示的外骨骼機器人中[12],利用肌電信號作為輸入控制信號來映射人體各關節運動,實現外骨骼機器人的主動運動。

圖1 外骨骼康復機器人

1 肌電信號預處理和特征提取

1.1 實驗方案

通過對人體健康四肢各關節運動范圍[13]的研究和對人體行走、上/下樓梯等常見動作的運動數據采集實驗分析,發現下肢三關節的運動大部分發生在矢狀面內,而在冠狀面和橫截面運動幅度較少。除此之外,本研究中的研究目標是將sEMG作為外骨骼機器人的輸入控制信號,為簡化結構,外骨骼機器人的設計和建模中均只考慮矢狀面的運動,故本研究中簡化為僅對矢狀面內的運動進行研究分析。髖關節、膝關節和踝關節的主要動作和運動相關的肌肉功能如表1所示。

表1 人體肌肉運動功能表

選定比目魚肌、腓腸外肌、腓腸內肌、脛骨前肌、股二頭肌、股內側肌、股外側肌、股直肌為行走過程的sEMG的采集源。選擇6名健康的受試者作為實驗對象,各受試者分別以0.8 km/h、1.2 km/h、1.6 km/h 3種速度在平地上直線行走120 s,肌電儀和運動采集系統同步采集行走時的各相關數據。肌電信號采集頻率為1 500 Hz,運動角度采集頻率為100 Hz。在神經網絡訓練階段,要求作為輸入的sEMG和輸出的三關節角度一一對應,采樣點需滿足式(1):

(1)

圖2為sEMG肌肉標記點。

(1)股直肌;(2)股外側肌;(3)股內側肌;(4)脛骨前肌;(5)股二頭肌;(6)腓腸內肌;(7)腓腸外肌和(8)比目魚肌圖2 sEMG采集標記點

1.2 sEMG預處理

由表面電極測量得到的原始肌電信號有很多噪聲的存在,需要進行濾波處理。圖3為脛骨前肌原始肌電信號。

圖3 脛骨前肌原始肌電信號

肌電信號的濾波降噪方法主要有移動中值法、移動均值法、小波濾波、經驗模態分解等[14]。小波濾波是對信號進行多尺度分解,采用一個小波基函數去逼近不同尺度上的信號,可以較好的刻畫信號的非平穩特征,但是相對的計算量過大,且降噪效果易受到信噪比的影響。經驗模態分解是一種數據驅動的自適應分解方法,降噪效果受噪聲幅度影響較大,如果噪聲幅度太大,會出現噪聲能量泄露到有用信息模態中,出現模態混疊問題,而且還存在邊界誤差的問題,去噪后信號細節丟失較多,且上述方法運算復雜,不適用對實時性要求高的信號處理。

移動中值濾波基于排序統計理論,基本思想首先確定一個以某個信號為中心點的鄰域,然后將鄰域中各信號值排序,取其中間值作為中心信號的新值。

移動均值濾波法基本原理是先確定一定的間隔,將間隔內的所有數據點加起來,然后除以該間隔內數據點的個數,取平均值。

(2)

式中:Es(i)是第i個數據點的平均值,T是窗口大小。

表2為不同數據窗口下的移動中值濾波和移動均值濾波信噪比值對比,隨著數據窗口的增加,信噪比的值有所不同,均呈逐漸減小趨勢,在數據窗口相同情況下,中值濾波的信噪比略大于均值濾波信噪比。

表2 不同數據窗口中值濾波和均值濾波信噪比值

圖4 移動中值濾波與移動均值濾波的脛骨前肌sEMG

圖4為移動中值濾波與移動均值濾波的sEMG效果圖,移動中值濾波較好的去掉脈沖噪聲。而移動均值濾波整體上較好的保留了原始信號,大大減小了隨機噪聲,尤其是尖銳的階躍響應信號。

考慮到本研究中主要利用sEMG信號預測人體下肢的運動角度實現外骨骼機器人的實時控制,因此對控制的實時性和快速響應要求較高。綜合考慮我們選用計算簡單、運算速度最快的移動均值濾波算法,考慮到運動信號和肌電信號的采集頻率關系以及高效性等,取T=5。

1.3 sEMG特征值提取

本研究中所分析的肌電信號是作為康復機器人的輸入控制信號,因此我們采用時域特征進行分析和特征值提取。主要特征值有肌電積分值、均方根值RMS(Root Mean Square)、標準差等。其中肌電積分值主要反映肌電信號隨時間的強弱變化。均方根值反映在時間維度上肌電信號振幅的變化特征,直接反映sEMG的瞬時功率,相比肌電積分值具有更加直接的意義。

(3)

式中:E(i)是經整流濾波處理后肌電信號的i第樣本點的電壓值,NR=15,表示對sEMG每15個采樣點提取一個特征值。圖5為脛骨前肌RMS肌電圖,提取RMS特征值后的sEMG所表達的運動特征更加明顯,為近似規律性的周期變化。

圖5 脛骨前肌的RMS肌電圖

2 降維處理

將整流濾波和特征值提取后的8通道sEMG,作為神經網絡輸入來預測三關節角度,對高維sEMG降維處理是必要的。除了降低數據量,加快后續計算速度,還體現在原始變量之間相互關聯,會導致預測結果不連貫;過多的變量會妨礙查找規律的建立;僅在變量層面上分析可能會忽略變量之間的潛在聯系。目前的降維方法主要有主成分分析法PCA(Principal Component Analysis)、核主成分分析KPCA(Kernel PCA)和線性判別分析LDA(Linear Discriminant Analysis)等[15]。LDA降維法主要是為了使得降維后的數據點盡可能容易地區分開。PCA降維法主要是將高維的數據映射到低維空間,并保留住盡可能多的原數據特性,是降維常用的方法之一。KPCA是PCA方法的推廣,基于核映射實現,相比于PCA,能夠挖掘到數據集中蘊含的非線性信息,對于非線性數據具有較好的處理能力,但KPCA的計算更加復雜,運算需花費大量時間,降維的性能受核函數的影響大。本研究中采用PCA降維方法。

2.1 標準化

特征值提取后的肌電信號需要標準化,否則會影響PCA的誤差計算及導致信息的缺失,不利于PCA降維算法的收斂。標準化中基準動作或是基準數據的選擇是關鍵,如果選擇恰當不僅可以減小sEMG個體差異給實驗數據分析帶來的誤差,而且可以簡化實驗和分析,提高實驗一致性和可重復性。

本研究中引入肌肉最大自主等長收縮MVC(Maximum Voluntary Contraction)作為標準化的基準。MVC的測量需要特定的動作完成,主要是提取對應肌肉束的最大肌電信號值。設定每個動作的持續時間是6 s,間隔15 s,重復5次完成,圖6為標準化MVC的基準動作示意圖,圖中黑色箭頭為運動方向,白色箭頭為阻力方向。

圖6 標準化MVC基準動作示意圖

標準化后的表面肌電信號為

(4)

式中:%MVC是標準化后的sEMG,Eavg是sEMG特征值,Emax是最大自主等長收縮的最大值。

圖7為脛骨前肌MVC肌電信號。對獲得的MVC進行與原始肌電信號相同的預處理操作,圖8為完成標準化的肌電信號。

圖7 脛骨前肌MVC肌電信號

圖8 標準化后的肌電信號

2.2 PCA數據降維

PCA降維方法的主要思想是:找出k個綜合變量來代替原來的m個變量,使得這些綜合變量盡可能地代表原來變量的信息量,而且彼此之間互不相關[16]。即將標準化后的sEMG從數據量m×N降為k×N,(m>k,m為sEMG原數據維數或通道數,m=8,k為降維后的sEMG維數,N是sEMG數據采樣點數,N=6 000)。運算步驟主要有:

①特征去中心化。即計算的不同維度之間的協方差矩陣Cx。

(5)

即

(6)

(7)

式中:X為原始sEMG數據集,mX表示X的均值向量,Cij(i,j=1,2,…,m)是協方差值。

②計算協方差矩陣的特征向量和特征值

由式(6)計算得知,協方差矩陣Cx是方陣,故可計算出它的特征向量和特征值。且為實對稱矩陣,解特征方程

λE-Cx=0

(8)

求出特征值(λ1,λ2,…,λm),并使其按大小順序排列,即λ1≥λ2≥…≥λm≥0,然后分別求出對應于特征值λi的特征向量A。

③計算主成分貢獻率及累計貢獻率

貢獻率的大小決定了k值的大小,即降維后矩陣維數的大小。k越大,即使用的原始sEMG中的特征向量越多,更多的保留了原來的數據特性,導致降維誤差越小。

主成分貢獻率ηi為:

(9)

累計貢獻率η為:

(10)

累計貢獻率η≥85%,則說明是肌電信號的主要信息,可替代原始sEMG作為神經網絡的輸入信號,在本研究中,η≥90%。

④新的主分量變換

Xk=Ak(X-mX)

(11)

式中:Xk是新的主分量數據集,Ak為前k個特征值對應的特征向量經轉置得到的矩陣Ak×N。表3為PCA降維處理得到的主成分累計貢獻率,由表可知,取k=6。

表3 PCA處理得到的主成分累計貢獻率

圖9為PCA降維后的獲得的6個主成分肌電信號,按貢獻率大到小從圖9(a)～圖9(f)依次為第一、第二、第三、第四、第五和第六主成分,第一主成分中,肌電信號的特征變化最為明顯,信號的分布更有周期性和規律性,突變信號減少。第六主成分的貢獻率最小,信號的變化較為平緩。

圖9 PCA降維后的主成分分量肌電信號

3 實驗結果與分析

3.1 神經網絡訓練及誤差分析

BP神經網絡是目前應用最為廣泛和成熟的神經網絡訓練和預測模型之一[17],本研究建立BP雙隱層神經網絡來訓練肌電信號,以各實驗對象在1.2 km/h的步速為例,分別對未降維處理的sEMG,經KPCA降維的sEMG和PCA降維的sEMG進行神經網絡訓練,通過對比訓練處理時間和預測關節角度均方根誤差RMSE和相關系數 γ 進行分析。

(12)

(13)

表4 未降維處理訓練誤差分析

表5 KPCA處理訓練誤差分析

表6 PCA處理訓練誤差分析

表7 3種數據處理方式比較

表8為采用BP神經網絡和徑向基神經網絡訓練PCA降維后肌電信號,預測獲得的髖關節誤差和相關系數等性能比較。BP神經網絡的訓練時間僅為徑向基神經網絡的1/3,其相關系數和均方根誤差均優于徑向基神經網絡。

表8 BP神經網絡和徑向基神經網絡性能比較

圖10為未經過降維處理的神經網絡訓練輸出髖關節角度誤差。在受試者正常行走過程中,髖關節角度是周期性變化的,預測值在波峰和波谷的時候誤差較大,絕對誤差最大達15°左右。圖11為KPCA處理后的訓練數據,髖關節相對誤差與未處理時相比有所下降。圖12為PCA降維處理后的訓練數據,髖關節在波峰和波谷的誤差明顯減小,絕對誤差最大為5°左右,預測的角度信號的平滑性和穩定性均有改善。綜上所述,本研究中采用PCA降維方法對sEMG進行處理和分析。

圖12 PCA降維的髖關節角度誤差分析

3.2 肌電信號映射實驗

神經網絡預測輸出的關節角度作為康復機器人的實時控制輸入信號,要求經過訓練后的網絡具有可重復性,因此需要對同一實驗對象的多次實驗數據進行可重復性預測,以驗證實驗的可復現性,而且預測值需和實際測量值對比分析其誤差,以驗證其可行性。

保存3.1節中已訓練好的6組神經網絡,然后重新采集6個實驗對象的5組肌電信號分別調用已保存的網絡進行各關節角度預測,分析其平均均方誤差、平均最大誤差角度和相關性,研究實驗可重復性效果。實驗采集方法和肌肉位置同2所示,圖13為肌電信號重復性采集和操作實驗。表9為經PCA降維處理后的sEMG可重復性預測誤差分析。

表9 PCA可重復性實驗預測誤差分析單位:(°)

圖13 肌電信號和運動角度采集實驗

圖14 PCA降維的膝關節預測角度分布

圖15 PCA降維的髖關節預測角度分布

4 結論

本文通過對人體行走步態中下肢相關肌肉的肌電信號進行預處理和PCA降維處理分析,提出基于MVC的標準化方法。通過對比分析和實驗,經過PCA降維處理的sEMG相比KPCA處理和未降維的sEMG,BP神經網絡可以更加快速和準確的預測關節角度,最后通過可重復性和映射實驗驗證了肌電信號作為機器人的控制信號,可以實現人體肌電信號和下肢運動的實時映射控制。本文主要就適用于實時控制下肢機器人對肌電信號進行了分析處理和研究,以便將預測得到的關節角度輸入到下肢康復機器人,實現如圖1中機器人的主動運動。目前正進行后續的研究,將建立更加完善神經網絡分析模型對處理后的sEMG進行訓練和髖關節、膝關節和踝關節角度的預測,降低各關節的預測誤差和提高精確度,開發肌電設備到機器人控制系統的交互接口,實現在線實時的預測和下肢機器人的主動運動控制。