Journal of Biomechanics 47 (2014) 596–601
Contents lists available at ScienceDirect
Journal of Biomechanics journal homepage: www.elsevier.com/locate/jbiomech www.JBiomech.com
Short communication
The effect of perturbing body segment parameters on calculated joint moments and muscle forces during gait Mariska Wesseling a,n, Friedl de Groote b, Ilse Jonkers a a b
KU Leuven, Department of Kinesiology, Human Movement Biomechanics, Tervuursevest 101, B-3001 Heverlee, Belgium KU Leuven, Department of Mechanical Engineering, Division PMA, Celestijnenlaan 300B, B-3001 Heverlee, Belgium
art ic l e i nf o
a b s t r a c t
Article history: Accepted 8 November 2013
This study examined the effect of body segment parameter (BSP) perturbations on joint moments calculated using an inverse dynamics procedure and muscle forces calculated using computed muscle control (CMC) during gait. BSP (i.e. segment mass, center of mass location (com) and inertia tensor) of the left thigh, shank and foot of a scaled musculoskeletal model were perturbed. These perturbations started from their nominal value and were adjusted to 7 40% in steps of 10%, for both individual as well as combined perturbations in BSP. For all perturbations, an inverse dynamics procedure calculated the ankle, knee and hip moments based on an identical inverse kinematics solution. Furthermore, the effect of applying a residual reduction algorithm (RRA) was investigated. Muscle excitations and resulting muscle forces were calculated using CMC. The results show only a limited effect of an individual parameter perturbation on the calculated moments, where the largest effect is found when perturbing the shank com (MScom,shank, the ratio of absolute difference in torque and relative parameter perturbation, is maximally 7.81 N m for hip flexion moment). The additional influence of perturbing two parameters simultaneously is small (MSmass þ com,thigh is maximally 15.2 N m for hip flexion moment). RRA made small changes to the model to increase the dynamic consistency of the simulation (after RRA MScom,shank is maximally 5.01 N m). CMC results show large differences in muscle forces when BSP are perturbed. These result from the underlying forward integration of the dynamic equations. & 2013 Elsevier Ltd. All rights reserved.
Keywords: Body segment parameters Inverse dynamics Dynamic simulations of motion Joint moments Muscle forces
1. Introduction Simulations of human motion calculate dynamics and muscle forces during different motions. These simulations rely on the use of musculoskeletal models. These models typically contain a description of segments, joints and muscle tendinous structures. For each body segment, body segment parameters (BSP) are defined: i.e. segment mass, center of mass location (com) and inertia tensor. In common practice, a generic parameter set derived from cadaveric work (e.g. Chandler et al., 1975) is used to scale BSP using segment length and body mass of the individual subject. To estimate subject-specific BSP, some more advanced techniques, like dual-energy X-ray absorptiometry, MRI and gammascanning, can be used (Ganley and Powers, 2004; Mungiole and Martin, 1990; Zatsiorsky et al., 1990). Differences found between methods to calculate BSP can be considerable (Ganley and Powers, 2004; Rao et al., 2006). However, the effect of erroneous estimation of BSP on calculated joint moments is still under debate. While some have argued that it is of high importance to accurately
n
Corresponding author. Tel.: þ 32 16 376 463; fax: þ 32 16 329 197. E-mail address:
[email protected] (M. Wesseling).
0021-9290/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.jbiomech.2013.11.002
estimate BSP (Andrews and Misht, 1996; Kingma et al., 1996; Rao et al., 2006), others report only small differences in calculated joint moments (Ganley and Powers, 2004; Pearsall and Costigan, 1999; Reinbolt et al., 2007; Silva and Ambrósio, 2004). These differences may relate to the different approaches. Whereas some authors (Pearsall and Costigan, 1999; Silva and Ambrósio, 2004) perturbed parameters individually, others examined the effect of simultaneous perturbation of the parameters (Ganley and Powers, 2004; Kingma et al., 1996; Rao et al., 2006). Also, the effects of individual and combined BSP perturbation on the dynamics of non-adjacent joints needs to be further explored. To minimize dynamic inconsistency resulting from modeling errors and marker kinematics inaccuracy, a residual reduction algorithm (RRA) was developed (Thelen and Anderson, 2006). This algorithm slightly changes the joint kinematics but specifically influences specific BSP of the model (i.e. segment mass and trunk com). As a result, the kinematics and ground reaction force (GRF) data will better satisfy the conservation of linear and angular momentum, i.e. be more dynamically consistent. The effect of erroneous estimation of BSP on calculated muscle forces during dynamic simulations of motion needs to be further explored. Some researchers have investigated the effect of perturbing the mass of a segment on calculated muscle forces (Dao
M. Wesseling et al. / Journal of Biomechanics 47 (2014) 596–601
et al., 2009; Van Den Bogert et al., 2012). These results suggest that an increase in mass increases muscle force, while the force profile is preserved. However, the effect of perturbing the com and inertial tensor on muscle forces calculated during a dynamic simulation has not been investigated. This study examines the influence of perturbing segment mass, com and inertial tensor on (1) joint moments calculated using an inverse dynamics procedure, before and after RRA and (2) lower limb muscle forces calculated using computed muscle control (CMC, Thelen et al., 2003) during gait.
2. Methods Data used to simulate a walking movement was taken from the example files installed with OpenSim 2.4.0 (Delp et al., 2007), called Gait2392. These data was collected on a subject (72.6 kg) walking on an instrumented treadmill (1.13 m/s) (John et al., 2012). Force data from the treadmill as well as 3d marker trajectories of 49 markers were used for the simulations. The musculoskeletal model consisted of 12 segments, i.e. pelvis, torso and left and right thigh, shank, talus, calcaneus and toes, 21 degrees of freedom and 92 musculotendon actuators (Delp et al., 2007). All analyses were performed in OpenSim 2.4.0. The musculoskeletal model was scaled using the marker trajectories of a static pose measurement. The scaled model was used for an inverse kinematics procedure based on measured 3d marker trajectories. The resulting kinematics was used in an inverse dynamics procedure to calculate the moments at the ankle, knee and hip joints. Next, RRA was applied to minimize the effects of errors in modeling and marker kinematics by changing the kinematics and adjusting the mass of the segments and torso com. The torso com is adjusted, since this segment is oversimplified and is expected to show the largest variation. This adjusted model was used as input for a computed muscle control (CMC) procedure that computed a set of muscle excitations, and therefore muscle forces, that drives the dynamic musculoskeletal model to track the original inverse kinematics solution (Fig. 1). In the sensitivity analysis, only muscles of the left leg were considered. In subsequent simulations the mass, com and inertia tensor of the left thigh, shank and foot were adjusted from 60% to 140% of the nominal value in steps of 10% (Table 1). The parameters were adjusted individually as well as for different combinations of BSP values, resulting in 728 adapted models. Thereafter, using the previously described inverse kinematics solution, the different steps of the workflow were repeated (Fig. 1). Sensitivity of inverse dynamics to a single parameter perturbation (MSp1 ) is expressed as MSp1 ¼
Mðp1 þ Δp1 Þ Mðp1 Δp1 Þ 2ðΔp1 =p1 Þ
where M(p1 þ Δp1) and M(p1 Δp1) are the absolute moments averaged over the gait cycle when perturbing parameter p1 by Δp1 and Δp1 respectively and Δp1/ p1 is the relative perturbation (De Groote et al., 2010). To determine the additional influence of simultaneously perturbing two parameters, a second order sensitivity (MSp1 þ p2 ) was used
MSp1 þ p2 ¼
597
As comparison of instantaneous differences in muscle forces between nominal and perturbed simulations could be misleading because of small temporal shifts in muscle force, CMC results are evaluated using the absolute difference in impulse, i.e. the time integral of the differences in muscle forces with the nominal situation.
3. Results For inverse dynamics, the largest influence of a single parameter perturbation was found for the hip flexion moment when perturbing the shank com (Fig. 2A). Changing the com in proximal direction had less influence on the calculated moments than changing it in the distal direction. For the mass and inertia tensor, a positive and negative perturbation had an almost equal effect (Appendix A). Further, a perturbation of the segmental mass had a smaller influence on inverse dynamics than a perturbation of the com, while the inertia tensor had the smallest influence (Fig. 2C). For combined perturbations, there was a small additional influence of combining a mass and com perturbation (Fig. 2B), indicating that the influence of com and mass perturbations in addition to single parameter perturbations was only small. For other combinations, the influence was even smaller (Fig. 2D). RRA only slightly adjusted the mass properties of the segments. The segment masses were all adjusted equally (with nominal values the adjustment was 0.012%, maximal was 5.30% when the mass and inertial tensor of the thigh were perturbed by 40% and com by þ 40%). The kinematics was changed maximally by 1.771, for the ankle angle, after a perturbation of all BSP of the foot by 40%. The torso com was adjusted in two directions, anteroposterior and mediolateral. In anteroposterior direction, it was adjusted most when all BSP of the thigh were perturbed by 40% (2.80 cm posterior). In mediolateral direction, a perturbation of the mass and inertia tensor of the thigh by þ 40% and the com by 40% led to the largest adjustment (4.34 cm left). Other large perturbations resulted in only small adjustments of the trunk com. In anteroposterior direction the smallest adjustments were made when all thigh BSP were perturbed by þ 40% (1.16 cm posterior, with nominal values the adjustment was 2.25 cm posterior). For the mediolateral direction, the adjustment is minimal when the mass and inertia tensor of the thigh were perturbed by 40% and the com by þ40% (2.50 cm left, with nominal values the adjustment was 3.40 cm). The sensitivity of the BSP after RRA showed that perturbation of a single parameter had a slightly smaller influence on the
Mðp1 þ Δp1 ; p2 þ Δp2 Þ Mðp1 þ Δp1 ; p2 Δp2 Þ Mðp1 Δp1 ; p2 þ Δp2 Þ þ Mðp1 Δp1 ; p2 Δp2 Þ 4ðΔp1 =p1 ÞðΔp2 =p2 Þ
where M(p1 þ Δp1, p2 þ Δp2), M(p1 þ Δp1, p2 Δp2), M(p1 Δp1, p2 þ Δp2) and M (p1 Δp1, p2 Δp2) are the absolute moments averaged over the gait cycle when perturbing parameters p1 and p2 by 7 Δp1 and 7 Δp2 and Δp1/p1 and Δp2/p2 are the relative perturbations. Sensitivity is considered to be high for MSp1 Z 10 Nm (De Groote et al., 2010) and MSp1 þ p2 Z 100 Nm. With these sensitivities, a parameter perturbation of 10% leads to a change in torque of 1 N m.
calculated moments after applying RRA (Fig. 2A and C). The additional influence for combined perturbations was again very small (Fig. 2B and D). CMC results showed large differences in individual muscle forces (Tables 2 and 3): large differences in impulse were observed, even
Fig. 1. Workflow used to calculate muscle forces. First, the model was scaled to the dimensions of the subject and an inverse kinematic procedure was done. Using this identical kinematic solution, models with different BSP values were used to calculate the inverse dynamics moments. RRA was applied before muscle forces were calculated using CMC.
598
M. Wesseling et al. / Journal of Biomechanics 47 (2014) 596–601
Table 1 Nominal values for the com, mass and inertial tensor of each segment. The com is expressed as the distance from the proximal joint center to the com. For the thigh and shank segments, the perturbation was only in proximal/distal direction. For the foot segment, the perturbation was made for the talus, calcaneus and toe segments at the same time in their own reference frame. The com perturbation for these segments was mainly in proximal/distal direction and slightly in medial/lateral and anterior/posterior direction. Com [m]
Mass [kg]
for small perturbations in BSP. Interestingly, the largest BSP perturbations did not correspond to the maximal differences in impulse. Large differences in impulse were mainly found for muscles acting around the ankle joint (Table 2 and Fig. 3).
4. Discussion
2
Inertia tensor [kg m ] Ixx
Iyy
Izz
Thigh Shank
0.20 0.18
8.98 3.58
0.13 0.049
0.034 0.0049
0.14 0.049
Foot Talus Calcaneus Toe
0.00 0.11 0.040
0.097 1.21 0.21
0.00097 0.0014 0.000097
0.00097 0.0038 0.00019
0.00097 0.0040 0.00097
This study analyzed the influence of perturbing the different BSP on inverse dynamic joint moments and muscle forces calculated using CMC. In accordance with previous research (Ganley and Powers, 2004; Pearsall and Costigan, 1999; Reinbolt et al., 2007; Silva and Ambrósio, 2004) this study found a limited effect of BSP perturbation on inverse dynamics. The largest effect was found by perturbing the shank com. For the com, a perturbation in distal direction had a larger influence on the calculated joint moments, indicating that the sensitivity for this parameter was
Fig. 2. Sensitivity of the absolute joint moments averaged over the gait cycle calculated using inverse dynamics (ID) and RRA for different BSP perturbations. Fig. 2A presents for each segment (thigh, shank and foot) the parameters with the highest sensitivity to a single parameter perturbation. Fig. 2C presents for each parameter (com, mass and inertia) the highest sensitivity to a single parameter perturbation irrespective of the segment. Likewise, Fig. 2B and D shows the highest second order sensitivity for each segment (Fig. 2B) and each parameter combination (Fig. 2D). The com is expressed in the local coordinate system of the segment. Also the minimal (Mminimal), maximal (Mmaximal) and nominal (Mnominal) absolute moments averaged over the gait cycle are represented.
M. Wesseling et al. / Journal of Biomechanics 47 (2014) 596–601
not linear (Appendix A). For the mass and inertia tensor, sensitivity was almost linear. Further, the additional influence of a combined perturbation of parameters was only very small. RRA slightly adapted the kinematics, segment masses and torso com, to improve dynamic consistency. Since the original model, kinematics, and GRF were not necessarily dynamically consistent due to modeling and measurement errors, some of the BSP perturbations potentially increased the dynamic consistency and therefore limited the adjustments made by RRA. In addition, RRA does not claim to correct the model but aims at increasing the dynamic consistency. This is important when CMC is subsequently used to calculate muscle forces based on forward integration to track the kinematics. Muscle forces calculated using CMC were largely affected by changes in BSP (Tables 2 and 3). Other authors also reported differences in muscle forces while perturbing the mass of segments (Dao et al., 2009; Van Den Bogert et al., 2012) however to a lesser extent. The larger differences found in this study could be inherent to the method used to calculate muscle forces: CMC combines a forward integration of the dynamic equations with a dedicated optimization method to calculate muscle excitations (Thelen et al., 2003). Perturbing BSP results in slightly different model parameters and hence a slightly different solution at the first time step. However, a small change in the initial conditions of a set of differential equations can produce a very large change in the final solution (Betts, 2001). This effect is augmented in every step of the integration. This is further aggravated by the variable step integrator used by CMC resulting in a completely different set of operations used to calculate muscle forces (Appendix B). When the simulation starts at a different time instant, the differences are only small at the beginning of the simulation, while further on large numerical errors are found (Fig. 4). When using static optimization (Appendix C), a technique that is not based on forward integration, muscle forces for the m. tibialis posterior are hardly affected by perturbations in BSP (Fig. 5). There are some limitations to the current study. A Monte Carlo analysis would have given more information about the distribution of the results and would be a valuable addition to the study. Further, motions with higher accelerations are expected to be more easily affected by perturbations in BSP and should be further
Table 2 Maximal differences in impulse with the unperturbed simulation after adjusting the BSP, with indication of the specific perturbations. Presented are the muscles that show the largest difference in impulse. Muscle
Nominal Segment Max impulse Perturbation impulse difference [N s] [N s] Mass Com [%] [%]
Inertia tensor [%]
Tibialis posterior 467
Thigh Shank Foot
323 213 438
40 20 40
0 40 40
10 20 30
Peroneus longus
136
Thigh Shank Foot
231 163 333
40 20 40
0 40 40
10 20 30
Tibialis anterior
326
Thigh Shank Foot
136 68.4 195
40 10 20
0 20 30
10 20 20
Table 3 Maximal differences in impulse with the unperturbed simulation after adjusting the BSP, with indication of the specific perturbations. Presented are the muscles that show the largest difference in impulse taken from the muscles acting around the hip, knee and ankle joint, with exception of the muscles presented in Table 2. Muscle
Nominal impulse [N s]
Segment Max impulse Perturbation difference [N s] Mass Com Inertia [%] [%] tensor [%]
Psoas major
408
Thigh Shank Foot
55.6 131 50.6
40 40 20
40 40 30
20 20 20
Soleus
379
Thigh Shank Foot
84.7 113 109
20 0 40
20 40 20
20 10 20
Gastrocnemius, medial head
677
Thigh Shank Foot
57.7 46.2 178
30 20 40
0 30 40
0 0 10
599
Fig. 3. Force calculated for the m. tibialis posterior using CMC for different perturbations of the thigh, shank and foot segments. The m. tibialis posterior is displayed, while the largest differences in impulse were found for this muscle. The continuous line represents the unperturbed situation (—). The other lines represent a perturbation of the different segments at which the largest difference in impulse was found for the thigh (- -), shank (- - -) and foot (——). The vertical lines represent respectively toe off (TO) and heel strike (HS).
600
M. Wesseling et al. / Journal of Biomechanics 47 (2014) 596–601
Fig. 4. Force calculated for the m. tibialis posterior using CMC for different perturbations of the thigh, shank and foot segments. The time instant of the start of the simulation is different from Fig. 3, namely at heel strike (the second vertical gray line in Fig. 3). The continuous line represents the unperturbed situation (—). The other lines represent the thigh (- -), shank (- - -) and foot (——) segments. The vertical line represents toe off (TO).
Fig. 5. Force calculated for the m. tibialis posterior using static optimization for different perturbations of the thigh, shank and foot segments. The simulation is started at the same time instant as for Fig. 3. The continuous line represents the unperturbed situation (—). The other lines represent the thigh (- -), shank (- - -) and foot (——) segments. The vertical lines represent respectively toe off (TO) and heel strike (HS).
investigated. Also, in different populations, e.g. children and obese subjects, results might be different. In conclusion, considering the influence of perturbing BSP on inverse dynamics, it is not important to accurately estimate these parameters for a walking movement. RRA made small changes to the model to increase the dynamic consistency of the simulation. When considering muscle forces calculated using CMC, large differences were found which result from the underlying forward integration. In addition, we showed that muscle forces calculated using CMC depend on the time instant the simulation is started.
Conflict of interest statement None.
Acknowledgment The support of the Agency for Innovation by Science and Technology (IWT-TBM) No. 100786) is gratefully acknowledged.
M. Wesseling et al. / Journal of Biomechanics 47 (2014) 596–601
Appendix A. Supplementary material Supplementary data associated with this article can be found in the online version at http://dx.doi.org/10.1016/j.jbiomech.2013.11.002.
References Andrews, J.G., Misht, S.P., 1996. Methods for investigating the sensitivity of joint resultants to body segment parameter variations. J. Biomech. 29, 651–654. Betts, J.T., 2001. Optimal control preliminaries. In: Burns, J.A. (Ed.), Practical Methods for Optimal Control Using Nonlinear Programming, second edition SIAM, Philadelpia, pp. 93–95. Chandler, R.F., Clauser, C.E., McConville, J.T., Reynolds, H.M., Young, J.W., 1975. Investigation of the inertial properties of the human body. Wright Patterson Air Force Base, OH. (AMRL Technical Report (TR-71-137). Dao, T.T., Marin, F., Ho Ba Tho, M.C., 2009. Sensitivity of the anthropometrical and geometrical parameters of the bones and muscles on a musculoskeletal model of the lower limbs. In: Proceeedings of the 31st Annual International Conference of the IEEE Engineering in Medicine and Biology Society. Minneapolis, Minnesota, USA. De Groote, F., Van Campen, A., Jonkers, I., De Schutter, J., 2010. Sensitivity of dynamic simulations of gait and dynamometer experiments to hill muscle model parameters of knee flexors and extensors. J. Biomech. 43, 1876–1883. Delp, S.L., Anderson, F.C., Arnold, A.S., Loan, P., Habib, A., John, C.T., Guendelman, E., Thelen, D.G., 2007. OpenSim: open-source software to create and analyze dynamic simulations of movement. IEEE Trans. Bio-Med. Eng. 54, 1940–1950. Ganley, K.J., Powers, C.M., 2004. Determination of lower extremity anthropometric parameters using dual energy X-ray absorptiometry: the influence on net joint moments during gait. Clin. Biomech. 19, 50–56.
601
John, C.T., Seth, A., Schwartz, M.H., Delp, S.L., 2012. Contributions of muscles to mediolateral ground reaction force over a range of walking speeds. J. Biomech. 45, 2438–2443. Kingma, I., Toussaint, H.M., Looze, M.P.D., Van Dieen, J.H., 1996. Segment inertial parameter evaluation in two anthropometric models by application of a dynamic linked segment model. J. Biomech. 29, 693–704. Mungiole, M., Martin, P.E., 1990. Estimating segment inertial properties: comparison of magnetic resonance imaging with existing methods. J. Biomech. 23, 1039–1046. Pearsall, D.J., Costigan, P.A., 1999. The effect of segment parameter error on gait analysis results. Gait Posture 9, 173–183. Rao, G., Amarantini, D., Berton, E., Favier, D., 2006. Influence of body segments' parameters estimation models on inverse dynamics solutions during gait. J. Biomech. 39, 1531–1536. Reinbolt, J.A., Haftka, R.T., Chmielewski, T.L., Fregly, B.J., 2007. Are patient-specific joint and inertial parameters necessary for accurate inverse dynamics analyses of gait? IEEE Trans. Bio-Med. Eng. 54, 782–793. Silva, M.P.T., Ambrósio, J.A.C., 2004. Sensitivity of the results produced by the inverse dynamic analysis of a human stride to perturbed input data. Gait Posture 19, 35–49. Thelen, D.G., Anderson, F.C., 2006. Using computed muscle control to generate forward dynamic simulations of human walking from experimental data. J. Biomech. 39, 1107–1115. Thelen, D.G., Anderson, F.C., Delp, S.L., 2003. Generating dynamic simulations of movement using computed muscle control. J. Biomech. 36, 321–328. Van Den Bogert, A.J., Hupperets, M., Schlarb, H., Krabbe, B., 2012. Predictive musculoskeletal simulation using optimal control: effects of added limb mass on energy cost and kinematics of walking and running. J. Sports Eng. Technol. 226, 123–133. Zatsiorsky, V., Seluyanov, V., Chugunova, L., 1990. In vivo body segment inertial parameters determination using a gamma scanner method. In: Berme, N., Cappozzo, A. (Eds.), Biomechanics of Human Movement: Application in Rehabilitation, Sports and Ergonomics. Bertec, Worthington, pp. 86–202.