一種智能膝關(guān)節(jié)假肢模型的阻尼參數(shù)優(yōu)化方法

張意彬，呂 杰，喻洪流

根據(jù)2006年第二次全國殘疾人抽樣調(diào)查顯示，中國各類的殘疾人群總數(shù)約8 296 萬，其中肢體殘疾人2 412 萬人，占比29.07%，截肢者約186 萬人，下肢截肢者占比達(dá)85%。 到目前為止，僅有44 萬人能夠安裝下肢假肢[1]。由于經(jīng)濟(jì)原因，絕大部分患者僅能穿戴傳統(tǒng)機(jī)械型假肢（passive prosthesis），只有少數(shù)患者能夠穿戴智能型假肢（active prosthesis）。 無論是機(jī)械型假肢，還是智能型假肢都是由阻尼系統(tǒng)來控制行走支撐相與擺動相的膝關(guān)節(jié)角度。健康人在行走過程中， 可以利用肌肉的控制完成支撐相和擺動相股骨、脛骨平臺的相對轉(zhuǎn)動[2]。與此類似，接受腔所連接的膝關(guān)節(jié)假肢可以在大腿殘肢的帶動下完成膝關(guān)節(jié)的屈曲，在臀部伸肌與阻尼系統(tǒng)的相互配合下，完成健康人的行走步態(tài)。臀部伸肌可以部分代償缺失的膝關(guān)節(jié)伸肌作用，阻尼系統(tǒng)則可以起到膝關(guān)節(jié)伸肌和屈肌的調(diào)節(jié)作用[3，4]。

目前多種膝關(guān)節(jié)阻尼模型已廣泛應(yīng)用到步態(tài)分析中，包括步態(tài)模型的建立及不同步速下的膝關(guān)節(jié)角度預(yù)測等，其中大部分研究的重點(diǎn)是利用阻尼系統(tǒng)控制膝關(guān)節(jié)的運(yùn)動[5~10]。 阻尼控制系統(tǒng)可分為機(jī)械性阻尼控制系統(tǒng)和微處理器阻尼控制系統(tǒng)。機(jī)械性阻尼控制系統(tǒng)是利用阻尼調(diào)整旋鈕的轉(zhuǎn)動來調(diào)節(jié)阻尼參數(shù)，而微處理器阻尼控制系統(tǒng)則是利用內(nèi)置的軟件對阻尼參數(shù)進(jìn)行調(diào)節(jié)，通過調(diào)整阻尼參數(shù)，控制阻尼器的能量消耗以維持膝關(guān)節(jié)在整個步態(tài)周期中的屈曲。在步態(tài)的支撐相，通過增大阻尼可以提高膝關(guān)節(jié)的穩(wěn)定性；在擺動期，通過減小阻尼使膝關(guān)節(jié)的屈曲變得更加容易[11]。 由此可見，通過適當(dāng)?shù)恼{(diào)節(jié)阻尼參數(shù)可以提高膝關(guān)節(jié)假肢的功能性，讓穿戴者得到近乎于健康人的步態(tài)， 使其在行走過程中具有更好步態(tài)對稱性。以往利用人工方法尋找合適的阻尼控制參數(shù)不僅效率低，而且準(zhǔn)確率也不高，易發(fā)生錯誤。 因此，建構(gòu)一種有效的膝關(guān)節(jié)假肢模型來自動尋找合適的阻尼控制參數(shù)就顯得更為重要。

1 材料與方法

1.1 實(shí)驗(yàn)材料

1.1.1 智能膝關(guān)節(jié)假肢樣機(jī)

智能膝關(guān)節(jié)假肢（AiKneeOne）是由上海理工大學(xué)康復(fù)工程與技術(shù)研究所的課題組自行研制的智能型膝關(guān)節(jié)假肢樣機(jī)（圖1），該樣機(jī)主要由膝關(guān)節(jié)假肢與阻尼控制單元兩部分組成，通過多組慣性傳感器及足底壓力傳感器獲取假肢的運(yùn)動狀態(tài)，可對多種行走模式（平路行走、上下坡、上下樓梯等）進(jìn)行智能識別，并利用微粒群優(yōu)化的反向傳輸算法實(shí)現(xiàn)了膝關(guān)節(jié)假肢在平路行走中的擺動相步速跟隨控制，詳細(xì)可參考近期已發(fā)表的相關(guān)文獻(xiàn)[12]。

圖1 智能膝關(guān)節(jié)假肢樣機(jī)（AiKneeOne）Fig. 1 Diagram of intelligent above-knee prosthesis prototype(AiKneeOne)

1.1.2 臨床材料

選擇測試對象5 例（表1），其中男性3 例，女性2例；年齡29~55 歲，平均年齡44.20 歲（標(biāo)準(zhǔn)差10.45歲）；身高1.52 ~ 1.78 m，平均身高1.650 m（標(biāo)準(zhǔn)差0.098 m）；體質(zhì)量55~75 kg，平均體質(zhì)量63.60 kg（標(biāo)準(zhǔn)差8.59 kg）；大腿長0.37~0.44 m，平均大腿長0.404 m（標(biāo)準(zhǔn)差0.027 m）；小腿長0.43~0.51 m，平均小腿長0.470 m（標(biāo)準(zhǔn)差0.029 m）；步幅0.56~0.66 m，平均步幅0.610 m（標(biāo)準(zhǔn)差0.038 m）。且所有測試對象在實(shí)驗(yàn)開始前均了解了實(shí)驗(yàn)的具體內(nèi)容，并同意參加實(shí)驗(yàn)。

1.2 方法

1.2.1 智能膝關(guān)節(jié)假肢（AiKneeOne）的模型建構(gòu)

人體下肢的肌肉、骨骼及韌帶等組織在大腦的協(xié)調(diào)控制下，完成了基本的運(yùn)動功能，包括站立、坐下、行走、跑步和上下樓梯等基本動作。 骨骼作為中心支架，通過關(guān)節(jié)連接在一起，肌肉附著于骨頭表面，作為動力源，三者共同組成了一個動力系統(tǒng)。 對于下肢而言，髖關(guān)節(jié)、膝關(guān)節(jié)及踝關(guān)節(jié)共同組成了一個多桿剛性系統(tǒng)，下肢可以簡化成如圖2（a）所示的多連桿的剛性系統(tǒng)（假肢模型model），阻尼單元（阻尼器）既可以作為動力源，在運(yùn)動過程中提供能量輸出，也可以作為能量緩沖器，提供阻尼力，減少和吸收振動，穩(wěn)定膝關(guān)節(jié)，其作用相當(dāng)于肌肉和韌帶的功能。 連桿相當(dāng)于骨骼，連桿之間以鉸鏈形式進(jìn)行連接，相當(dāng)于關(guān)節(jié)。該剛性連桿系統(tǒng)屬于平面機(jī)構(gòu)，而人體的骨骼屬于空間連桿系統(tǒng)，這里采用簡化模型，因?yàn)樾凶卟綉B(tài)分析更多的是關(guān)注其在矢狀面的運(yùn)動，忽略冠狀面和橫斷面的運(yùn)動。

圖2 膝關(guān)節(jié)假肢的耦合建模Fig.2 Diagram of above-knee prosthesis coupled modelling

圖2（b）為下肢步態(tài)采集裝置，大腿段、小腿段所安裝的角速度傳感器可分別獲取大腿、小腿的擺動角速度及角度信息，安裝于膝關(guān)節(jié)的角度傳感器可以獲取膝關(guān)節(jié)屈曲、伸展角度信息，足底薄膜壓力傳感器可以獲取足尖處及足跟處的足底壓力信息，目前該裝置已取得國家發(fā)明專利授權(quán)[13]。

根據(jù)人體行走的步態(tài)分析可知，整個行走周期可劃分為支撐相與擺動相。 在支撐相，膝關(guān)節(jié)主要提供扭矩，以保持膝關(guān)節(jié)在承重時的穩(wěn)定性，而在擺動相，膝關(guān)節(jié)更易于屈曲，利用小腿的擺出，為下一次邁步做好準(zhǔn)備。 由于兩個步態(tài)相位的受力狀態(tài)不同，因此需要對支撐相與擺動相分別建模。 對于支撐相，采用倒立擺模型進(jìn)行分析（圖3a）；而對于擺動相，則直接建立正擺模型（圖3b）。 在支撐相，體質(zhì)量（作用力向下）和足底反力（作用力向上）作用于假肢的阻尼器。在擺動相， 假肢的小腿部分逆重力方向反向轉(zhuǎn)動，沒有足底反力的作用，即f=0。

圖3 支撐相與擺動相模型Fig.3 Diagrams of stance and swing phase model

圖3 中虛線部分是連桿（L2）轉(zhuǎn)動θi角后阻尼器的變化位置。 阻尼器的安裝端采用的是曲柄滑塊結(jié)構(gòu)，采用該結(jié)構(gòu)方式，可以將阻尼器活塞的直線移動轉(zhuǎn)變成膝關(guān)節(jié)的轉(zhuǎn)動。 在整個步態(tài)周期內(nèi)，曲柄滑塊機(jī)構(gòu)可以根據(jù)行走狀態(tài)提供相應(yīng)的膝關(guān)節(jié)轉(zhuǎn)動角度。當(dāng)連接桿以θi角度發(fā)生變化時，則有：

其中：φi代表假肢模型的膝關(guān)節(jié)角度變化； 其他相關(guān)變量及參數(shù)詳見表2。 表中L1、L2的值取自智能膝關(guān)節(jié)假肢樣機(jī)（AiKneeOne）。

表2 假肢模型的相關(guān)變量與參數(shù)Tab.2 Related variables and parameters of prosthesis model

1.2.2 阻尼參數(shù)預(yù)估的優(yōu)化方法

為了使假肢側(cè)與健肢側(cè)的步態(tài)保持一致（步行過程中，雙側(cè)下肢的膝關(guān)節(jié)角度變化相同），在假肢模型（model）中采用Grid Search Optimization 優(yōu)化法[14]，以假肢模型與真實(shí)人體膝關(guān)節(jié)角度的最小均方根誤差為目標(biāo)，對阻尼參數(shù)（ωd，ζ）進(jìn)行迭代優(yōu)化：

其中：MMSE 為最小均方誤差 （minimum mean square error）；φ 代表假肢模型的膝關(guān)節(jié)角度；ρ 為真實(shí)人體膝關(guān)節(jié)（knee）的角度。 圖4 和圖2（b）給出了模型的阻尼參數(shù)優(yōu)化的具體方法。

圖4 阻尼參數(shù)優(yōu)化方法Fig.4 Diagram of optimization method of damping parameters

圖5 中虛線代表足跟處壓力傳感器獲取的壓力值， 點(diǎn)線代表足尖處壓力傳感器測得的壓力值，因此疊加兩者的測量值，即為整個足底的受力情況（實(shí)線）。 利用Matlab 函數(shù)庫中峰值檢測（Peakdet.M），Threshold=50，可以來查找波峰和波谷，從而確定P1與P2點(diǎn)。足底反力的起始點(diǎn)用P0表示，即足跟觸地。第一、二次出現(xiàn)的足底反力波峰點(diǎn)分別用P1和P2表示，即足跟最大受力點(diǎn)和足尖最大受力點(diǎn)。 Pend代表擺動相結(jié)束。 步態(tài)采集頻率為100 Hz，即Δt=0.01。在P0-P2（支撐相）之間，定義等時間間隔Δtj（tj-tj-1)內(nèi)的足底反力增量為Δfj， 增量為正值表示足底壓力增加，增量負(fù)值表示足底壓力減小。 例如，Δf1、Δf2均為正值， 表示t0到t1時刻、t1到t2時刻的足底壓力均是增加的；而Δf8為負(fù)，表示t7到t8時刻，足底壓力減小。

圖5 足底反力GRF 的變化曲線Fig.5 Curve of GRF(plantar reaction force)

根據(jù)S.S.Rao Mechanical Vibration 的二階振動理論[15]，即傅里葉級數(shù)展開型的阻尼運(yùn)動，其支撐相（P0- P2），L1與L2之間的角度θ 的變化由如下step函數(shù)表示：

其中：i 為采樣點(diǎn)；j 為總采樣間隔，即i = 1，2，…，j - 1，Δfi=Δfi-Δfi-1；k′彈簧的凈彈性系數(shù)， 其關(guān)系式為：

其中：L 代表膝關(guān)節(jié)軸至作用于大腿的人體重心距離[16]；k 為彈簧的彈性系數(shù)，可表示為：

在擺動相（P2至Pend），L1與L2之間的角度θ 的變化由如下step 函數(shù)表示：

其中：A 為P2（擺動初期）時的足底壓力最大值；M 為測試者的體質(zhì)量；ε 為步態(tài)相位差別常數(shù)，這里ε=0（支撐相與擺動相分別仿真）。

1.2.3 實(shí)驗(yàn)測試

測試對象被要求分別以0.8 m/s、1.2 m/s、1.6 m/s的步速平路行走5 min。 利用圖1（b）所示的下肢步態(tài)采集系統(tǒng)對行走步態(tài)的膝關(guān)節(jié)角度（ρ） 及足底壓力（f） 信息進(jìn)行實(shí)時采集， 其中足底壓力f 代入公式（4）、（7），可以分別求出支撐相、擺動相的θ 變化，再利用公式（1）、（2）進(jìn)一步求得φ。根據(jù)公式（3）提出的優(yōu)化函數(shù)，結(jié)合圖4 給出的優(yōu)化方法，對阻尼參數(shù)進(jìn)行迭代優(yōu)化。

2 結(jié)果

2.1 假肢模型的膝關(guān)節(jié)角度仿真結(jié)果

在P0~ P2區(qū)間，需將Δfi、ζ、ωd代入公式（4）對支撐相的膝關(guān)節(jié)角度（θstance）進(jìn)行仿真。 而在P0~Pend區(qū)間，則僅需將ωd代入公式（7）即可對擺動相的膝關(guān)節(jié)角度（θswing）進(jìn)行仿真。 圖6 給出了經(jīng)阻尼參數(shù)（ζ，ωd）優(yōu)化仿真后，不同步速的假肢模型（model）與真實(shí)人體膝關(guān)節(jié)的角度變化曲線。

圖6 膝關(guān)節(jié)假肢模型與人體下肢膝關(guān)節(jié)不同步速的膝關(guān)節(jié)角度變化比較Fig. 6 Comparison of knee joint angle change curves between knee prosthesis model and human lower limb knee joint at different walking speeds

ζ 和ωd的取值范圍分別是：在支撐相，ζR[0，1] =0，0.1，0.2，0.3，0.4，0.5，0.6，0.7，0.8，0.9，1.0；ωdR[1，20] =1，2，3，5，6，7， …，19，20； 在擺動相，ζR[0，0.2]=0，0.001，0.005，0.010，…，0.200；ζR的取值范圍取決于考慮阻尼環(huán)境下（ζ < 1）的情況，通過減小阻尼器的阻尼力來協(xié)助完成擺動相的屈曲態(tài)。 鑒于此，在仿真過程中，將阻尼比的范圍控制在0~0.2。 ωdR的范圍由ωd決定，一個步態(tài)周期內(nèi)，阻尼活塞在支撐期往復(fù)1 次，而擺動相則為半次，因此支撐相，擺動相， 其中：T 代表健肢完成一個完整的步態(tài)周期所需時間。圖6 給出了經(jīng)阻尼參數(shù)（ζ，ωd）優(yōu)化后的膝關(guān)節(jié)假肢模型與真實(shí)人體下肢膝關(guān)節(jié)的軌跡曲線對比。

利用均方根誤差（root mean squared error，RMSE）及相關(guān)系數(shù)可以描述真實(shí)人體膝關(guān)節(jié)角度和假肢模型角度之間的相似性，如公式（8）、（9）所示。

其中：Dataknee,i為真實(shí)膝關(guān)節(jié)的角度數(shù)據(jù)；Datamodel,i經(jīng)阻尼參數(shù)迭代優(yōu)化的假肢模型角度；n 為總采樣點(diǎn)。結(jié)果表明：假肢模型與真實(shí)膝關(guān)節(jié)的角度在支撐相的均方根誤差為3.8°±0.4°；擺動相相對小些，為1.7°±0.2°。 兩者在支撐相的相關(guān)性為0.86±0.06，擺動相為0.990±0.005，表現(xiàn)出較好的步態(tài)相似性（步態(tài)對稱性）。

2.2 假肢模型的阻尼參數(shù)優(yōu)化結(jié)果

將5 例測試對象的阻尼參數(shù)的最終優(yōu)化值導(dǎo)入Graph pad prism8.4.3，進(jìn)行數(shù)值統(tǒng)計分析（圖7、8）。圖中可以看出：隨著步速的變化，預(yù)估參數(shù)值的中值也相應(yīng)發(fā)生變化。 支撐相的阻尼頻率（ωd）的范圍位于5~20，阻尼比ζ 在0.001~0.990；而擺動相的ωd范圍則位于3.7~8.6，阻尼比ζ 在0.000 1~0.110 0。 在支撐相，ζ 隨著步速的增加而減小， 這是由于阻尼器對增強(qiáng)力的阻力減小造成的， 它與屈曲運(yùn)動的方向相反。這樣有利于維持由于步態(tài)增快而導(dǎo)致的支撐相踝關(guān)節(jié)角度增加。 在擺動相，ζ 值相近，這有利于控制膝關(guān)節(jié)的脛骨向前運(yùn)動。 無論是支撐相， 還是擺動相，ωd的值都會隨著步速的增加而增加。合適的阻尼參數(shù)不僅可以提高支撐相的膝關(guān)節(jié)穩(wěn)定，而且還有助于擺動相的假肢小腿前向擺動。

圖7 支撐相的阻尼參數(shù)（ωd，ζ）優(yōu)化結(jié)果Fig.7 Diagrams of optimized results of damping parameter(ωd，ζ)in stance phase

圖8 擺動相的阻尼參數(shù)（ωd，ζ）優(yōu)化結(jié)果Fig.8 Diagrsms of optimized results of damping parameter(ωd，ζ)in swing phase

3 討論

筆者研究的目的是尋找一種阻尼參數(shù)優(yōu)化的自動調(diào)節(jié)方法，提高智能膝關(guān)節(jié)假肢與健肢側(cè)下肢的步態(tài)對稱性。利用該方法可以計算出不同步速下的最優(yōu)阻尼參數(shù)預(yù)估值。優(yōu)化的阻尼參數(shù)值可以為智能膝關(guān)節(jié)假肢的穿戴調(diào)試提供參考依據(jù)，減少大腿截肢者穿戴智能膝關(guān)節(jié)假肢的適應(yīng)性訓(xùn)練時間。 另外，沒有安裝經(jīng)驗(yàn)的假肢裝配師也可以利用此方法快速確定阻尼參數(shù)調(diào)整參考值。筆者研究中所采用的智能膝關(guān)節(jié)假肢為研發(fā)樣機(jī)（AiKneeOne），但該假肢模型的應(yīng)用范圍并不局限于此，對于部分結(jié)構(gòu)類似的智能膝關(guān)節(jié)假肢也可以使用， 只要調(diào)整L1、L2、L 等參數(shù)的值即可。目前，實(shí)驗(yàn)所采集的數(shù)據(jù)均取自健康人，后續(xù)會考慮采集大腿截肢患者的數(shù)據(jù)。