无动力储能式截瘫助行外骨骼弹簧刚度优化
关鑫宇 , 季林红 , 王人成     
清华大学 机械工程系, 摩擦学国家重点实验室智能与生物机械分室, 北京 100084
摘要:无动力储能式截瘫助行外骨骼为不同身高、体重和损伤情况的截瘫患者提供精准化助行,储能弹簧作为外骨骼标准元件,其刚度为固定值。该文针对该外骨骼储能弹簧刚度参数选择问题,建立了截瘫患者穿戴无动力储能式截瘫助行外骨骼单侧下肢从站立中期到站立后期动力学模型,得到人体髋关节力矩与身高、体重、髋关节角度的关系,并提出以人体髋关节和弹簧力矩之和绝对值的积分为优化目标的弹簧刚度参数优化设计方法。研究发现:不同身高、体重组合情况下,外骨骼弹簧1刚度k1满足正态分布,弹簧2刚度k2不满足正态分布;k1k2的中位数3 180和1 279 N/m作为无动力储能式截瘫助行外骨骼储能弹簧刚度参数,适用于身高和体重在百分位数P1~P99间的患者。
关键词截瘫    行走    站立    外骨骼    储能    
Optimization of an unpowered energy-stored exoskeleton spring stiffness for spinal cord injuries
GUAN Xinyu, JI Linhong, WANG Rencheng     
Division of Intelligent and Biomechanical System, State Key Laboratory of Tribology, Department of Mechanical Engineering, Tsinghua University, 100084 Beijing, China
Abstract: The unpowered energy-stored exoskeleton can provide precise walking assistance for spinal cord injury patients for a specified body height, weight and injury level. Since the energy-stored springs are standard components of the exoskeleton, their stiffness should be constant. The aim of this study is to select the optimal stiffness of the exoskeleton springs. The hip joint moment was related to the body height, body weight and hip joint angle through an inverse dynamics model for spinal cord injury patients walking in the single support phase from mid stance to terminal stance. The optimization method minimized the integral of the absolute moment generated by the hip joint and springs. The results show that the stiffness of spring 1 (k1) has a normal distribution while the stiffness of spring 2 (k2) does not. The medians of the two spring stiffnesses (the median k1=3 180 N/m and the median k2=1 279 N/m) can apply to patients whose heights and weights are in the P1-P99 ranges.
Key words: spinal cord injury     walking     standing     exoskeleton     energy-stored    

截瘫指由胸腰段脊椎管内神经功能改变造成的损伤进而导致的下肢瘫痪,分为完全性和非完全性截瘫[1]。脊髓损伤后,患者由于长期卧床和缺乏运动,会相应地出现一些并发症如褥疮、骨质疏松、肌肉萎缩、血压低、深静脉血栓形成等等[2-3]。研究表明,站立和行走能够有效降低截瘫患者上述并发症的产生,并具有改善心血管系统、促进肠道规律性等益处。此外,站立和行走可以提高截瘫患者日常生活能力和社会参与程度[4]

无动力外骨骼目前在临床上被广泛用于辅助截瘫患者站立和行走,具有代表性的无动力式外骨骼包括膝-踝-足矫形器(KAFO)、往复式截瘫步行器(RGO)、改良式往复式截瘫步行器(ARGO)等[5]。这些外骨骼针对截瘫患者存在下肢肌肉无力及无法支配的状况,为关节提供了外部支撑,帮助患者保持站立,通过患者残存肌肉代偿作用,使患者实现迈步行走。由于损伤部位和损伤程度不同,患者保留的运动和感觉功能存在差异性,例如腰骶段脊髓损伤患者具有屈/伸髋能力,而胸段损伤的患者则不具备该能力[6]。学者们通过实验对无动力外骨骼适配性进行了研究,Sharma和Leimkuehler等发现KAFO适用于腰1(L1)或损伤平面更低的患者,而上述3种无动力外骨骼不能对胸段脊髓损伤患者之间的差异性提供精准化助行[7-8]

行走过程中,处于支撑期的一侧下肢从站立中期到站立后期髋关节做负功[9]。这些能量可储存于弹簧中,在摆动期辅助患者屈髋迈步。基于此原理,本文设计了无动力储能式截瘫助行外骨骼(见图 1),实验证明与KAFO和RGO相比,该外骨骼能够提高胸10(T10)截瘫患者行走的速度和增加其步长[10]。为满足不同身高、体重和损伤情况的患者精准化助行的需求,该外骨骼储能弹簧1和2两端均可沿相应的滑道移动,进行无级调节。弹簧1和2作为标准元件,其刚度和原长(弹簧1和2原长分别为50和75 mm)则设计为常量。

图 1 无动力储能截瘫助行外骨骼

针对外骨骼弹簧刚度参数选择问题,本文建立单侧下肢从站立中期到站立后期逆动力学模型,并通过优化的方法选择能够适用于身高和体重在百分位数P1~P99间的患者的外骨骼弹簧刚度参数。

1 下肢动力学建模 1.1 人体髋关节力矩

为研究人体髋关节力矩,本文在矢状面建立下肢动力学模型。模型由2个刚体和1个自由度构成:行走时,外骨骼膝关节锁定,无自由度,因此将大腿和小腿视为一个刚体;外骨骼小腿托为高分子材料,具有弹性,行走过程中脚踝产生背屈,因此将脚视为一个刚体。在一个步态周期中,单侧下肢从站立中期到站立后期的运动可简化为脚部固定不动,腿部绕脚踝旋转[11],下肢受力分析如图 2所示。

图 2 下肢动力学模型和受力分析示意图

对脚部进行受力分析,根据Newton-Euler公式得到:

$ {F_{{\rm{A}}x}} - {F_{{\rm{GRF}}x}} = 0, $ (1)
$ {F_{{\rm{A}}y}} - {F_{{\rm{GRF}}y}} - {m_f}g = 0, $ (2)
$ {M_{\rm{A}}} + {m_{\rm{f}}}g{L_{{\rm{fp}}}} + {F_{{\rm{GRF}}x}}{H_{\rm{f}}} - {F_{{\rm{GRF}}y}}{L_{{\rm{AG}}}} = 0. $ (3)

其中:FAxFAy分别为踝关节水平和竖直方向分力,FGRFxFGRFy分别为地面反力水平和竖直方向分力,mf为脚部重量,Hf为脚踝到地面高度,LfpLAG分别为脚部质心和足底压力中心到踝关节水平距离,MA为踝关节力矩。

同理,对腿部进行受力分析,根据Newton-Euler公式得到:

$ {F_{{\rm{H}}x}} - {F_{{\rm{A}}x}} = {m_1}{L_{1{\rm{d}}}}\ddot \theta \cos \theta - {m_1}{L_{1{\rm{d}}}}{{\dot \theta }^2}\sin \theta , $ (4)
$ {F_{{\rm{H}}y}} - {F_{{\rm{A}}y}} = - {m_1}{L_{1{\rm{d}}}}\ddot \theta \sin \theta - {m_1}{L_{1{\rm{d}}}}{{\dot \theta }^2}\cos \theta , $ (5)
$ \begin{array}{*{20}{c}} {{M_{\rm{H}}} - {M_{\rm{A}}} + {m_1}g{L_{1{\rm{d}}}}\sin \theta + }\\ {{F_{{\rm{H}}x}}{L_l}\cos \theta - {F_{{\rm{H}}y}}{L_1}\sin \theta = {I_1}\ddot \theta ,} \end{array} $ (6)
$ {I_1} = {m_1}R_1^2. $ (7)

其中:FHxFHy分别为髋关节水平和竖直方向分力;ml为腿部重量;Lld为腿部质心到脚踝距离;Ll为腿长;θ为髋关节伸展角度,单位为(°);MH为踝关节力矩;Il为腿部转动惯量;Rl为腿部回转半径。

根据式(1)—(7),得到髋关节力矩:

$ \begin{array}{*{20}{c}} {{M_{\rm{H}}} = {F_{{\rm{GRF}}y}}{L_{{\rm{AG}}}} - {m_{\rm{f}}}g{L_{{\rm{fp}}}} - {F_{{\rm{GRF}}x}}{H_{\rm{f}}} - }\\ {{m_1}g{L_{1{\rm{d}}}}\sin \theta - {F_{{\rm{GRF}}x}}{L_1}\cos \theta - {m_1}{L_{1{\rm{d}}}}{L_1}\ddot \theta + }\\ {{m_{\rm{f}}}g{L_1}\sin \theta - {F_{{\rm{GRF}}y}}{L_1}\sin \theta + {m_1}g{L_1}\sin \theta + {m_1}R_1^2.} \end{array} $ (8)

根据Johnson等[12]对截瘫患者使用无动力外骨骼行走机理研究,得到髋关节角度从站立中期到站立后期与步态周期百分比关系:

$ \theta = - 94.1 + 4.53t - 0.045{t^2}. $ (9)

其中t为步态周期百分比。

行走过程中足底压力中心(COP)到脚后跟距离与步态周期百分比近似成线性关系[13],进而得到COP到脚踝水平距离:

$ {L_{{\rm{AG}}}} = \left( {0.013t + 0.24} \right){L_{\rm{f}}} - {L_{{\rm{AH}}}}. $ (10)

根据Fineberg等[14]对截瘫患者使用外骨骼行走在一个步态周期内竖直地反力变化研究,得到髋关节从站立中到站立后期地面反力竖直方向分力与体重关系:

$ \begin{array}{*{20}{c}} {{F_{{\rm{GRF}}y}} = 0.01 \cdot \left( {215.33 - 14.32t + } \right.}\\ {\left. {0.38{t^2} - 0.003{t^3}} \right)Wg.} \end{array} $ (11)

其中:W为体重,单位为kg;g为重力加速度,本文取9.8 m/s2

图 3 一个步态周期髋关节角度、足底压力中心到脚后跟距离和地面反力竖直分力变化(阴影部分为站立中期到站立后期)

文[12, 15]发现截瘫患者行走过程中地面反力水平方向分力FGRFx非常小,因此忽略FGRFx作用;患者行走速度缓慢,因此忽略髋关节角加速度$\ddot \theta $影响。人体模型参数包括各部位长度、重量、质心位置和回转半径见表 1,其中,HLfLl分别为身高、脚长和腿长,单位均为m。将式(9)—(11)和人体模型参数代入式(8)得到人体髋关节力矩与人体身高、体重及髋关节角度关系:

$ \begin{array}{l} {M_{\rm{H}}} = \left( { - 2.99 \times {{10}^{ - 5}}{T^2} + 4.26 \times {{10}^{ - 4}}\sqrt {{T^3}} - } \right.\\ \;\;\;\;\;\;\;\;\;\;8.91 \times {10^{ - 4}}T - 1.16 \times {10^{ - 2}}\sqrt T + 6.34 \times \\ \;\;\;\;\;\;\;\;\;\;\left. {{{10}^{ - 2}}} \right)WgH + \left( { - 1.57 \times {{10}^{ - 3}}\sqrt {{T^3}} + } \right.\\ \;\;\;\;\;\;\;\;\;\;7.95 \times {10^{ - 3}}T + 2.63 \times {10^{ - 2}}\sqrt T - \\ \;\;\;\;\;\;\;\;\;\;\left. {3.67 \times {{10}^{ - 1}}} \right)WgH\sin \theta ; \end{array} $ (12)
$ T = 19.89 - \theta ,0 \le \theta \le 19.5. $ (13)
表 1 人体模型参数[15]
部位 长度 质量 回转半径(远心端)
Lf=0.152H
Hf=0.039H
Lfp=0.5Lf
LAH=0.328Lf
mf=0.0145W
Ll=0.491H
Lld=0.553Ll
ml=0.161W Rl=0.65Ll

1.2 弹簧力矩

以无动力储能式外骨骼髋关节转轴为原点,建立坐标系(见图 4)。根据文[16]提供数据得到中等身材人体髂前上棘A、髂后上棘P、大转子O在矢状面投影点坐标:A为(0.056, 0.05),P为(—0.075, 0.065),O为(0, 0)。OAOPOGOF为滑道1—4,APOGF为相应的滑道端点,G为(0.07, 0),F为(0, -0.1)。

图 4 无动力储能外骨骼储能元件示意图

B(B′)和D为弹簧1的端点,C(C′)和E为弹簧2的端点,令弹簧1端点B与滑道端点G重合(B′与G′重合),弹簧2端点C与滑道端点F重合(C′与F′重合),D(xDyD)和E(xEyE)分别沿APOP移动,设

$ {t_1} = \frac{{{l_{AD}}}}{{{l_{AP}}}}, $
$ {t_2} = \frac{{{l_{OE}}}}{{{l_{OP}}}}, $

其中lADlAPlOElOP分别为ADAPOEOP的长度。

$ {x_D} = 0.056 - 0.131{t_1}, $
$ {y_D} = 0.05 + 0.015{t_1}, $
$ {x_E} = - 0.075{t_3}, $
$ {y_E} = 0.065{t_3}. $

弹簧1和2相对于O点的力矩M1M2

$ {M_1} = {F_1}{L_1}, $
$ {M_2} = {F_2}{L_2}. $

其中:F1F2分别为弹簧1和2拉力,L1L2分别为弹簧1和2力臂,L1OB′D的距离;L2OC′E的距离。弹簧1和2原长分别为0.05和0.075 m,k1k2分别为弹簧1和2的刚度。F1F2L1L2

$ {F_1} = {k_1}\left( {\sqrt {{{\left( {{X_D} - {X_D}\cos \theta } \right)}^2} + {{\left( {{Y_D} - {Y_D}\sin \theta } \right)}^2}} - 0.05} \right), $ (14)
$ {F_2} = {k_2}\left( {\sqrt {{{\left( {{X_E} - {Y_C}\sin \theta } \right)}^2} + {{\left( {{Y_E} - {Y_C}\cos \theta } \right)}^2}} - 0.075} \right), $ (15)
$ {L_1} = \frac{{\left| {{X_B}{X_D}\sin \theta + {X_B}{y_D}\cos \theta } \right|}}{{\sqrt {{{\left( {{y_D} + {X_B}\sin \theta } \right)}^2} + {{\left( {{X_D} + {X_B}\cos \theta } \right)}^2}} }}, $ (16)
$ {L_2} = \frac{{\left| {{X_E}{y_C}\cos \theta - {y_E}{y_C}\sin \theta } \right|}}{{\sqrt {{{\left( {{y_E} - {y_C}\cos \theta } \right)}^2} + {{\left( {{X_E} - {y_C}\sin \theta } \right)}^2}} }}. $ (17)
2 弹簧刚度选择 2.1 优化目标

无动力截瘫助行外骨骼采用双弹簧目的在于:1)行走时,通过弹簧将单侧下肢从站立中期到站立后期过程中,髋关节所做的功转化为弹性势能并储存;2)站立时,弹簧1和2相对于O点的力矩相等或接近,增加患者站立稳定性。假设增加弹簧后,截瘫患者行走步态不发生改变。人体髋关节、弹簧1和2产生力矩的和为

$ M\left( {\theta ,{k_1},{k_2},{t_1},{t_2}} \right) = - {M_{\rm{H}}} - {M_1} + {M_2}. $ (18)

优化设计的目标是选择最优的k1k2t1t2参数组合,使髋关节从站立中期到站立后期所做功尽可能储存在弹簧中,定义优化目标函数为

$ {\rm{Min}}\int_0^{19.5} {\left| {M\left( {\theta ,{k_1},{k_2},{t_1},{t_2}} \right)} \right|{\rm{d}}\theta } . $
2.2 约束条件

无动力储能外骨骼弹簧刚度的设计约束有2类:站立稳定性约束和边界约束。截瘫患者站立状态下(θ=0)力线与正常人不同,其下肢及骨盆前倾,躯干后倾,患者髋关节力矩不为零[17]。为使患者站立时髋关节尽量不受弹簧力矩干扰,站立稳定性约束条件为

$ \left| {{M_1} - {M_2} - 0.017WH} \right| \le 1, $

边界约束条件为

$ 0 \le {t_1} \le 1, $
$ 0 \le {t_2} \le 1, $
$ 0 \le {k_1} \le 5500, $
$ 0 \le {k_2} \le 2500. $
2.3 优化方法和结果

根据《中国成年人人体尺寸(GB/T 10000-1988)》提供数据,中国成年人身高、体重P1和P99信息如表 2所示。将中国成年人身高、体重P1和P99对应的值作为边界值,涵盖了患者因脊髓损伤引起的变化后的体重阈值,设置身高步长0.025 m,体重步长1 kg,运用Matlab优化工具箱函数fmincon求解在不同身高、体重组合条件下,使$\int_0^{19.5} {\left| {M\left( {\theta, {k_1}, {k_2}, {t_1}, {t_2}} \right)} \right|} {\rm{d}}\theta $最小的k1k2t1t2参数组合。

表 2 中国成年人身高、体重参数
参数 P1 P99
身高/m 1.543 1.814
体重/kg 44.93 83.07

通过上述优化方法,得到不同身高、体重组合情况下弹簧刚度k1k2(见图 5)以及相对应的髋关节、弹簧力矩和绝对值的积分$\int_0^{19.5} {\left| {M\left( {\theta, {k_1}, {k_2}, {t_1}, {t_2}} \right)} \right|} {\rm{d}}\theta $ (见图 6)。

图 5 弹簧刚度分布

图 6 从站立中期到站立后期髋关节、弹簧力矩和的绝对值的积分

从上述优化结果中筛选出身体质量指数(BMI)在18.5~24.99范围内的身高、体重组合所对应的k1k2,并绘制直方图(见图 7)。本文对满足条件的k1k2进行S-W(Shapiro-Wilk)检验,判断其是否满足正态分布。结果显示k1呈正态分布(p=0.061),集中在2 500~3 500 N/m范围内,本文选择k1的中位数3 180 N/m作为无动力储能式截瘫助行外骨骼弹簧1刚度参数;k2不呈正态分布(p<0.05),中位数小于平均值,由于k2集中在700~1 600 N/m,选择k2的中位数1 279 N/m作为无动力储能式截瘫助行外骨骼弹簧2刚度参数。有关k1k2计算结果如表 3所示。

图 7 弹簧刚度分布直方图

表 3 k1k2相关计算结果
参数 中位数 平均数 最大值 最小值
k1/(N·m-1) 3 180 3 178 4 381 2 056
k2/(N·m-1) 1 279 1 339 640 2 554

3 结论

站立和行走能够有效降低脊髓损伤并发症的产生,对截瘫患者生理和心理层面均有积极作用。通过调节无动力储能式截瘫助行外骨骼储能弹簧端点位置,可为不同身高、体重和损伤情况的患者提供精准化助行。储能弹簧作为外骨骼标准化元件,其刚度为固定值,需通过优化方法进行选择。

本文针对该外骨骼储能弹簧刚度参数选择问题,首先建立了截瘫患者穿戴无动力储能式截瘫助行外骨骼单侧下肢从站立中期到站立后期动力学模型,得到人体髋关节力矩与身高、体重、髋关节角度关系,并提出以人体髋关节和弹簧力矩之和的绝对值的积分为优化目标的弹簧刚度优化设计方法。最后,选择弹簧1和2刚度的中位数3 180和1 279 N/m作为无动力储能式截瘫助行外骨骼储能弹簧刚度参数。

参考文献
[1] Kirshblum S C, Burns S P, Biering-Sorensen F, et al. International standards for neurological classification of spinal cord injury (revised 2011)[J]. The Journal of Spinal Cord Medicine, 2011, 34(6): 535–546. DOI:10.1179/204577211X13207446293695
[2] Maimoun L, Fattal C, Micallef J, et al. Bone loss in spinal cord-injured patients:From physiopathology to therapy[J]. Spinal Cord, 2006, 44(4): 203–210. DOI:10.1038/sj.sc.3101832
[3] Johnson W B. Walking mechanics of persons who use reciprocating gait orthoses[J]. Journal of Rehabilitation Research and Development, 2009, 46(3): 435.
[4] Colombo G, Wirz M, Dietz V. Driven gait orthosis for improvement of locomotor training in paraplegic patients[J]. Spinal Cord, 2001, 39(5): 252–255. DOI:10.1038/sj.sc.3101154
[5] Audu M L, To C S, Kobetic R, et al. Gait evaluation of a novel hip constraint orthosis with implication for walking in paraplegia[J]. IEEE Transactions on Neural Systems & Rehabilitation Engineering, 2011, 18(6): 610–618.
[6] Drake R, Vogl A W, Mitchell A W. Gray's Anatomy for Students[M]. New York: Bellevue Literary Press, 2014.
[7] Middleton J W, Sinclair P J, Smith R M, et al. Postural control during stance in paraplegia:Effects of medially linked versus unlinked knee-ankle-foot orthoses[J]. Archives of Physical Medicine & Rehabilitation, 1999, 80(12): 1558–1565.
[8] Abe K. Comparison of static balance, walking velocity, and energy consumption with knee-ankle-foot orthosis, walkabout orthosis, and reciprocating gait orthosis in thoracic-level paraplegic patients[J]. Journal of Prosthetics & Orthotics, 2006, 18(18): 87–91.
[9] Walsh C J, Paluska D, Pasch K, et al. Development of a lightweight, underactuated exoskeleton for load-carrying augmentation[C]//Proceedings 2006 IEEE International Conference on Robotics and Automation. Orlando, FL, USA:IEEE, 2006:3485-3491. http://www.academia.edu/3047024/Biomechanical_considerations_in_the_design_of_lower_limb_exoskeletons
[10] Guan X, Liu Y, Gao L, et al. Trunk muscle activity patterns of a person with spinal cord injury walking with different un-powered exoskeletons:A case study[J]. Journal of Rehabilitation Medicine, 2016, 48(4): 390–396. DOI:10.2340/16501977-2065
[11] Zajac F E, Neptune R R, Kautz S A. Biomechanics and muscle coordination of human walking:Part Ⅰ:Introduction to concepts, power transfer, dynamics and simulations[J]. Gait & Posture, 2002, 16(3): 215–232.
[12] Johnson W B, Fatone S, Gard S A. Walking mechanics of persons who use reciprocating gait orthoses[J]. Journal of Rehabilitation Research & Development, 2009, 46(3): 435–446.
[13] Chiu M C, Wu H C, Chang L Y. Gait speed and gender effects on center of pressure progression during normal walking[J]. Gait & Posture, 2013, 37(1): 43–8.
[14] Fineberg D B, Asselin P, Harel N Y, et al. Vertical ground reaction force-based analysis of powered exoskeleton-assisted walking in persons with motor-complete paraplegia[J]. The Journal of Spinal Cord Medicine, 2013, 36(4): 313–321. DOI:10.1179/2045772313Y.0000000126
[15] Winter D A. Biomechanics and Motor Control of Human Movement[M]. New York: John Wiley & Sons, 2009.
[16] Reynolds H M, Snow C C, Young J W. Spatial geometry of the human pelvis[J]. Spatial Geometry of the Human Pelvis, 1982.
[17] Lin K H, Lu T W, Hsu P P, et al. Postural responses during falling with rapid reach-and-grasp balance reaction in patients with motor complete paraplegia[J]. Spinal Cord, 2008, 46(3): 204–209. DOI:10.1038/sj.sc.3102100