An Optimization Method Tracking EMG, Ground Reactions Forces, and Marker Trajectories for Musculo-Tendon Forces Estimation in Equinus Gait

Abstract

In the context of neuro-orthopedic pathologies affecting walking and thus patients' quality of life, understanding the mechanisms of gait deviations and identifying the causal motor impairments is of primary importance. Beside other approaches, neuromusculoskeletal simulations may be used to provide insight into this matter. To the best of our knowledge, no computational framework exists in the literature that allows for predictive simulations featuring muscle co-contractions, and the introduction of various types of perturbations during both healthy and pathological gait types. The aim of this preliminary study was to adapt a recently proposed EMG-marker tracking optimization process to a lower limb musculoskeletal model during equinus gait, a multiphase problem with contact forces. The resulting optimization method tracking EMG, ground reactions forces, and marker trajectories allowed an accurate reproduction of joint kinematics (average error of 5.4 ± 3.3 mm for pelvis translations, and 1.9 ± 1.3° for pelvis rotation and joint angles) and ensured good temporal agreement in muscle activity (the concordance between estimated and measured excitations was 76.8 ± 5.3 %) in a relatively fast process (3.88 ± 1.04 h). We have also highlighted that the tracking of ground reaction forces was possible and accurate (average error of 17.3 ± 5.5 N), even without the use of a complex foot-ground contact model.

Introduction

Walking is often considered to be the most important activity in daily living (Chiou et al., 1985). The ability to move without pain, fatigue, or major gait deviation is closely related to quality of life (Cuomo et al., 2007; van Schie, 2008). Many neuro-orthopedic pathologies (e.g., cerebral palsy, stroke) induce impairments (i.e., paresis, muscle overactivity, soft tissue contractures, and bone deformities) that compromise normal movement. Consequently, the goal of many therapeutic interventions is to minimize gait deviations in patients. In order to improve these interventions, understanding the possible mechanisms of these gait deviations, and being able to identify the causal motor impairments is of primary importance (Davids et al., 2004; Gough and Shortland, 2008; Wren et al., 2011). Currently, the relationship between motor impairments and gait deviations is unclear (Bonnefoy-Mazure et al., 2016), and there is a lack scientific evidence for these relationships due to the inherent complexity of the human neuromusculoskeletal system during dynamic tasks such as walking (Armand et al., 2017). Compared to existing approaches (e.g., pathologic models, experimental procedures with human subjects, robots with human-like gait), in silico neuromusculoskeletal simulations of normal and pathological gait could provide additional insight into gait deviations (Armand et al., 2017). The advantage of neuromuscular simulations as a method is that large numbers of simulations may be performed relatively quickly, and without the ethical issues involved with performing invasive and lengthy experiments in vulnerable patient groups such as those with neuromuscular deficits.

To date, simulations reported in the literature are often limited to the analysis of the consequences of isolated impairments on gait such as muscle weakness (van der Krogt et al., 2009; Thompson et al., 2013) or muscle spasticity (Jansen et al., 2014). Most of these studies only report possible muscular compensations (adaptations) that occur due to muscular redundancy–results are achieved by tracking the normal gait kinematics and then applying a perturbation to the model. Very few studies have been based on a numerical framework allowing kinematic adaptations in response to more varied perturbations. Within the TLEMsafe project (https://www.tlemsafe.eu), Fluit et al. (2014b) combined an optimized inverse model and a ground reaction force predictive model to simulate lower limb kinematics after the removal of the rectus femoris and the vastus lateralis from the model. Santos et al. (2017) also proposed a numerical framework based on direct collocation and an optimal control package to simulate lower limb kinematics after the introduction of a weakening of the triceps surae and the tibialis anterior, or after increasing the ankle joint stiffness. These two approaches represent a first step toward the simulation of pathological gait. However, neither were able to reproduce muscle co-contractions. While this capacity was not necessarily needed in these studies, this feature is essential to establish a pathological gait simulator that would be able to reproduce physiological gait adaptations biofidelically.

From a methodological point of view, inverse dynamics-based approaches (such as static optimization) are commonly used due to their computational efficiency (Erdemir et al., 2007), but are not appropriate for predictive simulations. Moreover, static optimization method underestimates or neglects antagonist co-contractions unless hybrid approaches are used (Brookham et al., 2011; Son et al., 2012). On the other hand, forward dynamics-based approaches are often criticized for being time-consuming–several studies report convergence times in the hundreds of hours range (Anderson and Pandy, 2001). Despite this disadvantage, these methods have the potential to predict new movements, such as an adaptation in response to a perturbation. For example, state-of-the-art algorithms used in conjunction with existing musculoskeletal models–like direct collocation (Santos et al., 2017) and direct multiple shooting (Bélaise et al., 2018a,b)—can be used to solve forward dynamics problems in a timely manner. Recently, Bélaise et al. (2018a,b) introduced an EMG-marker tracking optimization method to predict musculo-tendon forces in a co-contraction case. Based on simulated datasets of upper limb movements, the authors showed the importance of tracking both marker trajectories and EMG, in particular to reproduce muscle co-contractions. To the best of our knowledge, such an approach has never been applied on experimental gait records with muscle co-contractions.

The objective of our project is to establish a computational framework appropriate for predictive simulations of healthy and pathological gait, that is able to reproduce muscle co-contractions, and that allows for the introduction of various kind of perturbations on the model (e.g., therapy-related, surgery-related, pathology-related perturbations). This preliminary study represents a first step toward this project by adapting the computational framework proposed by Bélaise et al. (2018a,b) to a lower limb musculoskeletal model during gait. This framework has been tested for this purpose on a type of pathological gait known as equinus gait.

Methods

Lower Limb Musculoskeletal Model

A generic three-dimensional musculoskeletal model of the lower limb [Lower Extremity Model, OpenSim (Delp et al., 1990)] was adapted for our study (Figure 1). This model consists of five rigid segments: the pelvis, right thigh, patella, shank, and foot. Twenty-six markers were associated with these segments by virtual palpation to reproduce the experimental marker locations (Table 1; see section Dataset). To simplify the dynamic optimizations to a two-dimensional motion in this preliminary study, the original degrees of freedom (DoF) were reduced to three DoFs for the pelvis-ground joint (vertical translation, translation in the direction of walking, pelvis tilt) and one DoF (flexion-extension modeled as a hinge joint) at the hip, knee, and ankle joints. Joints were actuated by the muscle torques resulting from 17 muscle lines of action (Table 2), and the pelvis DoFs were actuated by three generalized forces applied on the pelvis. The path, optimal length, maximal isometric force, tendon slack length, and pennation angle of each muscle lines of action were derived from the original model (Delp et al., 1990).

Figure 1

Segment lengths were scaled to the dataset used in this study (see section Dataset) using OpenSim 3.3 (Delp et al., 2007) by minimizing the distance between experimental and model markers placed on bony landmarks (Table 1). All components of the model that depend on bone lengths (e.g., muscle attachment points, optimal fiber length), segment masses, and inertial parameters were also scaled. The resulting scaled model was transferred to the bioRBD musculoskeletal modeling package (https://github.com/pyomeca/biorbd) based on the Rigid Body Dynamic Library (Felis, 2017). This model was defined by 29 states (six generalized joint positions and their six related velocities, 17 muscle activations) and 20 controls corresponding to the 17 muscle neural excitations plus three generalized forces driving the three pelvis DoFs.

Equations of Motion and Activation Dynamics

The generalized accelerations of the rigid multibody system were computed using a forward dynamics approach for given generalized joint positions q, joint velocities , and generalized forces τ:

where M is the inertia matrix, C is the external contact Jacobian matrix, R is the Lagrange multipliers vector corresponding to the ground reaction forces (GRF), N is the non-linear effects (Coriolis and centrifugal forces) vector, and G is the gravity effects vector. It was assumed that contact points have a null acceleration and velocity throughout the entire contact phase. In line with equinus gait, one fixed contact point was defined on the forefoot for the entire contact phase (see section Dynamic Optimizations). Generalized forces were divided into driving the three pelvis DoFs, and corresponding to the net joint torques due to the musculo-tendon forces F_mt, where is the generalized muscular lever arms matrix and L_mt the vector of muscle line of action lengths. F_mt were computed from muscle activations a using a Hill-type muscle model with a generic force-length-velocity relation f (Zajac, 1989):

where is the maximal isometric forces vector, l_m is the muscle fiber lengths vector, and v_m is the muscle fiber velocities vector. Again, musculo-tendon forces were divided into F_mt1 (with related activations a₁ and excitations e₁), corresponding to the muscles for which electromyographic (EMG) records were available, and F_mt2 (with related activations a₂ and excitations e₂) where EMG measurements were unavailable. Muscle activation dynamics was implemented as a set of first-order differential equations (Buchanan et al., 2004):

where e(t) is the muscle neural excitations at time t. Time constants t_act and t_deact (for activation and deactivation) were set at 10 and 40 ms, respectively (Thelen et al., 2003).

Dynamic Optimizations

As proposed in Bélaise et al. (2018a), controls and state variables were simultaneously optimized using an EMG-marker tracking optimization process. Because the model has a reduced muscular redundancy and because a generic model was used as opposed to a subject-specific model, the optimal maximal isometric forces were also identified during this process.

The optimization consisted of the minimization of the differences between predicted M_p and measured M_m marker trajectories in the sagittal plane, and between predicted e_1p and measured e_1m (EMG envelop) muscle neural excitations (corresponding to the muscles for which EMG records are available). This tracking optimization was extended to the minimization of the differences between predicted R_p and measured R_m GRF in the sagittal plane. This tracking was necessary to impose physiologic generalized forces to the pelvis (τ₁), i.e., generalized forces that compensate for the missing upper part of the body and contralateral lower limb.

To predict the activity of the muscles for which tracking was not possible (i.e., muscles for which EMG records were not available), the objective function J was written to find the least squared muscle activations a₂ that produced the prescribed marker trajectories, muscle neural excitations, and GRF (during stance phase only):

where w_M, w_e, w_R, and w_L are weighting factors adjusted to the relative importance of each term, T_i is the duration of the current stage (see section Simulations) and N_i is the related number of time frames.

This objective function was minimized under three sets of constraints. First, boundary conditions were applied on the state and the control variables. In this study, the range of motion of each DoF and related velocities were set to physiologic values (Table 3), while activations and excitations were bounded between 0 and 1. Second, the velocity of the contact point was constrained to be null at the first frame and its acceleration to be null at each time frame (see section Equations of Motion and Activation Dynamics). Third, periodicity was ensured by constraining the first and last time point of the cycle to have similar values in terms of hip, knee, ankle joint angles, and velocities, pelvis velocities, muscle excitations, and GRF.

Simulations

Each dynamic optimization was solved using a direct multiple shooting algorithm with MUSCOD-II (Leineweber et al., 2003). Three phases were defined in the gait step: (1) the stance phase (with an external contact between foot and ground), (2) the swing phase (no external contact), and (3) the first frame of the next stance phase following the impact between the foot and the ground. These stages were divided into 25, 25, and 1 multiple shooting intervals, respectively. For the sake of simplicity, the first stage started just after the collision impact between the foot and the ground. The duration of each stage was fixed to the measured value.

The initial guess was set to the measured values for the joint positions and velocities, 1% for all activations and excitations, and 0 for the controls corresponding to the generalized forces related to the pelvis DoFs. Weighting factors were set to w_M = 30 (except for the foot markers for which w_M = 50 to ensure the correct position of the contact point), w_e = 1, w_R = 0.05, and w_L = 1. These weighting factors were adjusted empirically to set values around 1, in order to ensure optimization convergence and produce simulation results close to the experimentally measured data.

Dataset

The previously defined method was evaluated on a dataset of emulated equinus gait. All data were recorded on a healthy volunteer (male, 35 years old, 165 cm, 66 kg) without any neuro-orthopedic conditions. This participant gave written informed consent prior to his inclusion and the protocol was conformed to the Declaration of Helsinki and approved by the National Research Ethics Committee of Luxembourg (201805/01).

The 3D trajectories of 26 reflective cutaneous markers (bilateral iliac anterior and posterior spines, right leg great trochanter, medial and lateral femoral epicondyles, peroneal head, tibial tuberosity, medial and lateral malleoli, 1st, 2nd, and 5th proximal and distal metatarsal heads, calcaneum, completed by a four-marker cluster on the thigh and on the shank) (Figure 1) were recorded using a 10-camera optoelectronic system (OQUS-4, Qualisys AB, Sweden) sampled at 200 Hz. Markers were placed by anatomical palpation (Table 1) following the recommendation of van Sint Jan (2007) by an experienced user. GRF and moments were recorded using two side-by-side force plates (OR6-5, AMTI, USA) sampled at 2,000 Hz. The EMG activity of nine right leg muscles (tibialis anterior, soleus, gastrocnemius medialis, vastus medialis, rectus femoris, semitendinosus, biceps femoris long head, gluteus medius, gluteus maximus) was collected with a wireless electromyographic system (DTS clinic, Noraxon, USA) sampled at 2,000 Hz. The EMG surface electrodes were placed following the recommended standard of the Surface EMG for a Non-Invasive Assessment of Muscles (SENIAM) project (Hermens et al., 2000).

All data were imported under Matlab (R2018a, The MathWorks, USA) using the ezc3d package (https://github.com/pyomeca/ezc3d). Marker trajectories were interpolated when necessary using a cubic spline and smoothed by a 4th order low-pass Butterworth filter with a cutoff frequency of 6 Hz. Generalized kinematics () were computed using an extended Kalman filter (Fohanno et al., 2014) following the segmental coordinate systems defined in the original generic three-dimensional lower limb musculoskeletal model [Lower Extremity Model, OpenSim (Delp et al., 1990)]. GRF were smoothed by a 4th order low-pass Butterworth filter with a cut-off frequency of 15 Hz. Raw EMG signals were band pass filtered (4th order) between 30 and 300 Hz, rectified, and EMG envelops were obtained by a 4th order low-pass Butterworth filter with a cut-off frequency of 25 Hz. EMG envelops were then normalized to their respective maximal voluntary activation (Gaudet et al., 2018).

The participant was asked to mimic an equinus gait by producing voluntarily controlled co-contractions of the muscles crossing the ankle joint to restrain ankle dorsiflexion. Eight trials were recorded and the related right steps were analyzed in this study.

Analysis

In order to evaluate the capacity of the model to reproduce the measured gait pattern and muscle excitations under the mechanical properties and constraints imposed to the model, a set of goodness-of-fit parameters were employed. Root mean square error (RMSE) and coefficient of determination (R²) were computed to assess the differences in intensity and shape, respectively, between measured and estimated excitations, joint angles and GRF.

Only estimated muscle excitations corresponding to the measured EMG envelops (Table 2) are presented in this analysis. The coefficient of determination (CC) (Giroux et al., 2013) was computed for the muscles for which EMG data was recorded. This method uses active/inactive state concordance between the estimated muscle excitations and normalized EMG envelopes to compute a coefficient of concordance defined as the percentage of concordance elements.

Results

The convergence time of the eight optimizations using MUSCOD-II was 3.88 ± 1.04 h on an Intel® Core™ i5-3570 CPU @3.4 GHz. Estimated and measured muscle excitations, musculo-tendon forces, joint angles, and GRF are reported in Figures 2–5, respectively. Goodness-of-fit parameters (RMSE, R² and CC) are reported in Table 4.

Figure 2

Figure 3

Figure 4

Figure 5

Considering the tracked muscle excitations (i.e., muscles for which EMG records are available, see Table 2), the temporal muscle activity of the model was in good overall agreement with the experimental measurements, with an average CC of 76.8 ± 5.3%. RMSE values were generally low with an average value of 0.2 ± 0.1 (values were adimensioned between 0 and 1). However, RMSE was found to be higher for the gastrocnemius medialis (0.3 ± 0.1) and the tibialis anterior (0.4 ± 0.1). For all muscles, the correlation remained low with an average R² at 0.02 ± 0.52. Regarding all other model muscles, for which EMG was not tracked, the optimized muscle excitations were higher than those estimated for muscles for which EMG was tracked. This is the case for R_GLUT_MAX2 and R_GLUT_MAX3 compared to R_GLUT_MAX1, R_GLUT_MED1, and R_GLUT_MED3 compared to R_GLUT_MED2, R_SEMIMEM compared to R_SEMITEN, R_VAS_INT, and R_VAS_LAT compared to R_VAS_MED, and for R_GAS_LAT compared to R_GAS_MED. The same results are observed on the estimated musculo-tendon forces. These forces are ranged between 0 and 2,000 N, with the highest peak forces obtained for R_SEMITEN, R_VAST_INT, R_VAST_LAT, and R_SOLEUS.

With regard to the pelvis position/orientation and joint angles, the model estimations were generally in agreement with the experimental measurements. Average RMSE were 5.4 ± 3.3 mm for the pelvic translations, and 1.9 ± 1.3° for pelvic rotations and joint angles. However, RMSE was found to be higher for the ankle joint (4.0 ± 0.9°). Considering all degrees of freedom, the correlation remained very high with an average R² at 0.94 ± 0.09.

For the GRF, the model estimations were found to be in agreement with the experimental measurements (the average RMSE is 17.3 ± 5.5 N). For these forces, the correlation remained high with an average R² at 0.97 ± 0.03.

Discussion

The main objective of this study was to adapt an EMG marker tracking optimization process to solve a forward dynamics problem on a 3D musculoskeletal model of the lower limb during equinus gait. To the best of our knowledge, the use of a direct multiple shooting algorithm on a musculoskeletal model with the tracking of measured EMG, marker trajectories, and GRF has never been performed to date. As already shown by Bélaise et al. (2018a,b), this approach allows for an accurate reproduction of joint kinematics and ensures temporal fidelity in muscle activity with improved computational time compared to traditional forward dynamic approaches. We have also highlighted that the tracking of GRF could be performed accurately, even without the use of a complex foot/ground contact model.

Limitations

A primary limitation of this preliminary study is that it was based on a small number trials for a single task, performed by a single participant. As such, limited conclusions can be drawn from this paper.

A second limitation is that only one contact point was defined at the forefoot and it was only constrained to null velocity and acceleration during the whole contact phase. While this approach was in line with an equinus gait and was able to accurately reproduce the tracked GRF, this definition cannot be applied during normal gait trials, during which several contact points should be defined (Fluit et al., 2014a). Elastic contact elements (Peng et al., 2018), artificial muscle-like actuators (Fluit et al., 2014a) or distance and velocity-dependent force models (Jung et al., 2016) should be adapted to the present model to extend its use to normal gait.

The proposed musculoskeletal model also only consisted of the pelvis and the right lower limb. A forward dynamic approach was made possible by replacing the forces and moments produced by the opposite lower limb and the upper part of the body by a set of generalized forces acting on the pelvis. With such an approach, the individual contribution of the opposite lower limb and the upper part of the body to the muscle activity and joint contact forces of the right lower limb cannot be evaluated individually. It would thus be an important next step to complete the missing body segments of the present musculoskeletal model in order to obtain a full body musculoskeletal model. Several full body musculoskeletal models (Rajagopal et al., 2016; Bassani et al., 2017) have been proposed in the literature and could be transferred to the bioRBD musculoskeletal modeling package. These segments could be actuated by joint torques instead of muscles.

In addition to this, while we used a 3D musculoskeletal model, DoFs were reduced to only allow a two-dimensional motion in the sagittal plane. It is however, established that walking is a locomotion task that is performed in sagittal, coronal, and transversal planes (Perry and Burnfield, 2010). In particular, patients often develop compensatory movements in the coronal plane when pathological impairments result in reduced foot clearance capacity in the sagittal plane (Chantraine et al., 2016). Despite this simplification, the accuracy of kinematic tracking observed in the present results suggest that there is potential for the EMG marker tracking optimization process to perform 3D gait motion simulations.

Finally, most of the parameters of the Hill-type muscle model were kept with generic property definitions in this study. Only optimal maximal isometric forces were identified during dynamic optimizations, while muscle optimal lengths, tendon slack lengths, and maximal isometric forces are usually identified in similar studies (Sartori et al., 2014; Pizzolato et al., 2015). This may explain the high excitations observed in muscles contained in a group, for which tracking of excitation was applied to other muscles. The introduction of further muscle parameters could be implemented in future studies, and would only be expected to impact the convergence time of the simulations.

Muscle Activity

The present study support the results of Bélaise et al. (2018a,b), which demonstrated that EMG tracking could be an efficient way to reproduce measured muscle excitations in simulations. While some amplitude differences appeared for certain tracked muscles, the temporal activity was generally reproduced with good fidelity. This outcome is crucial to ensure the ability of the model to produce muscle co-contractions. Similar approaches have already been proposed in the literature–EMG-driven musculoskeletal models have also been used to accurately reproduce muscle excitation patterns observed on EMG records (Shao et al., 2009; Sartori et al., 2012). However, these models are constrained to have as many muscle lines of action as EMG available in the dataset. To overcome this limitation, some authors have proposed hybrid approaches that combine EMG-driven and static optimization methods (Lloyd and Besier, 2003; Moissenet et al., 2014; Sartori et al., 2014). The drawback with this strategy is that by minimizing the difference between the motor joint moments computed by EMG-driven and inverse dynamics, kinematics may not be accurately reproduced. In that sense, the EMG-marker tracking algorithm proposed by Bélaise et al. (2018a,b) is a novelty. As this method tracks joint kinematics based on marker trajectories rather than joint moments as in a forward dynamic approach, the error diffusion is minimized and the simulation outputs reproduce the experimental kinematics more faithfully.

In our trials of emulated equinus gait, co-contractions of the ankle dorsiflexors, and plantarflexors can be observed during early stance to stabilize the joint in this specific posture. It is interesting to observe that, while the gastrocnemius medialis and the soleus (muscles for which EMG records were tracked) were contracted during this phase, the gastrocnemius lateralis (a muscle for which EMG records were not measured) was not in our simulations. Although use of EMG is somewhat limited to available hardware, muscle locations and the signal quality, measurement should be prioritized toward muscles presumed to be active during the task being investigated. In our case, focusing on a greater number of muscles crossing the ankle joint would have brought more relevant information to the model. A similar recommendation has already been proposed by Sartori et al. (2014); these authors suggested to prioritize EMG use on muscles “that reflect the patient's non-physiological muscular behavior.”

Although the EMG-GRF-marker tracking algorithm was able to reproduce physiological muscle activity, two points must be considered. First, we observed that when EMG was not tracked, the optimized muscle excitations and musculo-tendon forces were higher than the ones estimated when muscle EMG was tracked. As pointed out in the limitations of the study, this over-estimation is perhaps related to the use of generic muscle model parameters. For example, if the parameters applied to a muscle group would tend to limit its capacity to produce a motor joint moment, a higher muscle excitation would be required to reproduce the experimental measurements. This effect is further exacerbated if the excitation of a muscle in this group is constrained to a low level in accordance with the experimental EMG tracking, as the excitation of the other muscles of the group will have to compensate for this reduced excitation in the other muscle. Second, due to equinus (results in increased plantarflexion), the capacity of the triceps surae to produce a plantarflexion moment is reduced (Delp et al., 1990), in particular during the push-off phase. Thus, hip flexor recruitment may be increased in this gait pattern to pull the leg forward (Romkes and Schweizer, 2015). However, in our model the primary hip flexors, i.e., the iliopsoas muscles, were not included. van der Krogt et al. (2012) showed that an increased activation of the rectus femoris may be developed to compensate for a weakness of the primary hip flexors. In our case, the absence of the iliopsoas muscles (equivalent to a complete reduction of strength in these muscles) may have induced the increased rectus femoris excitations observed during the simulations compared to the experimental measurements. Because the rectus femoris is a bi-articular muscle (i.e., hip flexor, knee extensor), an increased excitation of this muscle used to assist in hip flexion would simultaneously act to reduce knee flexion, which would have require compensation from knee flexors in order to maintain experimental kinematics (van der Krogt et al., 2012). This could explain the non-physiological activity of the triceps surae (ankle plantarflexors) during pre-swing and early swing observed in our simulations, and the increased tibialis anterior activity (ankle dorsiflexor) used to balance ankle flexion due to the increased activity of ankle plantarflexors. All these observations support the need for a more comprehensive, full body musculoskeletal model, as already discussed in section Limitations.

Kinematics and Ground Reaction Forces

Unlike inverse dynamics-based optimization approaches (i.e., static optimization), where measured kinematics and calculated joint torques are the input constraints to the optimization problem, in a forward dynamics approach it is essential to assess the accuracy of reproduced kinematics. This is generally assessed by tracking experimental kinematics (Erdemir et al., 2007; Chèze et al., 2015). In the present study, marker trajectories were tracked rather than joint kinematics, as proposed by Bélaise et al. (2018a,b). This tracking was able to produce accurate marker trajectories as highlighted in Bélaise et al. (2018a) with a tracking residual of 0.31 ± 0.32 cm during elbow flexion. This approach, used in the present study gave accurate kinematics with a maximum RMSE obtained at the ankle joint < 5°, a threshold recognized as critical for clinical interpretation (McGinley et al., 2009). The main errors appeared at the end of the stance phase in most of the DoFs. This issue may be associated with the non-physiological activity of the triceps surae and the rectus femoris during this phase, as discussed in section Muscle Activity. Our simulations tended to reduce the plantarflexion induced by the emulated equinus gait. By increasing ankle dorsiflexion, the triceps surae moment arm was increased (Delp et al., 1990) and a minimal ankle plantarflexion moment was kept.

Interestingly, GRF were estimated accurately without an advanced foot-ground contact model definition. The use of simple generalized forces applied on the pelvis, designed to compensate for the absence of the opposite lower limb and the upper part of the body in our model, acted as efficient reserve actuators (Modenese et al., 2013) to provide the forces required to track experimental GRF. This approach may thus present a valid means by which to manage external forces and moments in the dynamic equation when the interactions between the upper limb and/or the contralateral lower limb with the ipsilateral lower limb are not known. Otherwise, as already discussed in section Limitations, a full body musculoskeletal model would be recommended.

Clinical Perspectives

The simulation of equinus gait represents an important clinical issue in the context of toe walking, a common gait deviation observed in many pathologies such as cerebral palsy, myopathy, and neuropathy (Armand et al., 2007). Numerical simulations may present a useful tool in this context of identifying potential biomechanical causes of this deviation, such as: pre-tibial muscle weakness, inadequate ankle dorsiflexors activity, ankle plantarflexors contracture, and/or spasticity, excessive voluntary ankle plantarflexion in compensation for quadriceps weakness, knee flexor contracture caused by overactivity of the hamstring, combined spasticity of the hamstring and ankle plantarflexors, and leg length discrepancy (Armand et al., 2007; Perry and Burnfield, 2010). In more general terms, there is a need for a numerical framework allowing for the introduction of pathology (Santos et al., 2017), treatment, or surgical intervention-related (Fox et al., 2009) perturbations in a model, and the analysis of their impact on the structures of the musculoskeletal system during daily living activities. However, before clinical applications, models have to be evaluated and validated (Hicks et al., 2015). It will thus be necessary that we assess the capacity of our approach to produce physiological musculo-tendon forces and joint contact forces. Validation datasets, such as the ones made available by Bergmann et al. (2016) and Fregly et al. (2012), should thus be tested on our numerical framework in the future.

Conclusion

In conclusion, we have improved the recent EMG-marker tracking optimization method to a multiphase cyclic movement with GRF. This numerical framework was successfully tested on a dataset of equinus gait for which our approach was able to estimate lower-limb kinematics, GRF and muscle activity with reasonable accuracy.

References

• 1

  AndersonF. C.PandyM. G. (2001). Dynamic optimization of human walking. J. Biomech. Eng.123, 381–390. 10.1115/1.1392310

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 2

  ArmandS.MoissenetF.de CoulonG.Bonnefoy-MazureA. (2017). Identifying and understanding gait deviations: critical review and perspectives. Mov. Sport Sci. Sci. Mot.2017, 77–88. 10.1051/sm/2017016

  • CrossRef
  • Google Scholar

• 3

  ArmandS.WatelainE.RouxE.MercierM.LepoutreF.-X. (2007). Linking clinical measurements and kinematic gait patterns of toe-walking using fuzzy decision trees. Gait Posture25, 475–484. 10.1016/j.gaitpost.2006.05.014

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 4

  BassaniT.StucovitzE.QianZ.BriguglioM.GalbuseraF. (2017). Validation of the AnyBody full body musculoskeletal model in computing lumbar spine loads at L4L5 level. J. Biomech.58, 89–96. 10.1016/j.jbiomech.2017.04.025

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 5

  BélaiseC.Dal MasoF.MichaudB.MombaurK.BegonM. (2018a). An EMG-marker tracking optimisation method for estimating muscle forces. Multibody Syst. Dyn.42, 119–143. 10.1007/s11044-017-9587-2

  • CrossRef
  • Google Scholar

• 6

  BélaiseC.MichaudB.Dal MasoF.MombaurK.BegonM. (2018b). Which data should be tracked in forward-dynamic optimisation to best predict muscle forces in a pathological co-contraction case?J. Biomech.68, 99–106. 10.1016/j.jbiomech.2017.12.028

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 7

  BergmannG.BenderA.DymkeJ.DudaG.DammP. (2016). Standardized loads acting in hip implants. PLoS ONE11:e0155612. 10.1371/journal.pone.0155612

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 8

  Bonnefoy-MazureA.SagawaY.PomeroV.LascombesP.De CoulonG.ArmandS. (2016). Are clinical parameters sufficient to model gait patterns in patients with cerebral palsy using a multilinear approach?Comput. Methods Biomech. Biomed. Engin.19, 800–806. 10.1080/10255842.2015.1064112

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 9

  BrookhamR. L.MiddlebrookE. E.GrewalT.DickersonC. R. (2011). The utility of an empirically derived co-activation ratio for muscle force prediction through optimization. J. Biomech.44, 1582–1587. 10.1016/j.jbiomech.2011.02.077

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 10

  BuchananT. S.LloydD. G.ManalK.BesierT. F. (2004). Neuromusculoskeletal modeling: estimation of muscle forces and joint moments and movements from measurements of neural command. J. Appl. Biomech.20, 367–395. 10.1123/jab.20.4.367

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 11

  ChantraineF.FilipettiP.SchreiberC.RemacleA.KolanowskiE.MoissenetF. (2016). Proposition of a classification of adult patients with hemiparesis in chronic phase. PLoS ONE11:e0156726. 10.1371/journal.pone.0156726

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 12

  ChèzeL.MoissenetF.DumasR. (2015). State of the art and current limits of musculo-skeletal models for clinical applications. Mov. Sport Sci. Sci. Mot.90, 7–17. 10.1051/sm/2012026

  • CrossRef
  • Google Scholar

• 13

  ChiouI. I.BurnettC. N.ChiouI. I.BurnettC. N. (1985). Values of activities of daily living. A survey of stroke patients and their home therapists. Phys. Ther.65, 901–906. 10.1093/ptj/65.6.901

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 14

  CuomoA. V.GamradtS. C.KimC. O.PirpirisM.GatesP. E.McCarthyJ. J.et al. (2007). Health-related quality of life outcomes improve after multilevel surgery in ambulatory children with cerebral palsy. J. Pediatr. Orthop.27, 653–657. 10.1097/BPO.0b013e3180dca147

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 15

  DavidsJ. R.OunpuuS.DeLucaP., aDavisR. B. (2004). Optimization of walking ability of children with cerebral palsy. Instr. Course Lect.53, 511–522. 10.2106/00004623-200311000-00028

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 16

  DelpS. L.AndersonF. C.ArnoldA. S.LoanP.HabibA.JohnC. T.et al. (2007). OpenSim: open-source software to create and analyze dynamic simulations of movement. IEEE Trans. Biomed. Eng.54, 1940–1950. 10.1109/TBME.2007.901024

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 17

  DelpS. L.LoanJ. P.HoyM. G.ZajacF. E.ToppE. L.RosenJ. M. (1990). An interactive graphics-based model of the lower extremity to study orthopedic surgical procedures. IEEE Trans. Biomed. Eng.37, 757–767. 10.1109/10.102791

  • CrossRef
  • Google Scholar

• 18

  ErdemirA.McLeanS.HerzogW.BogertA. J. (2007). Model-based estimation of muscle forces exerted during movements. Clin. Biomech. Bristol. Avon.22, 131–154. 10.1016/j.clinbiomech.2006.09.005

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 19

  FelisM. L. (2017). RBDL: an efficient rigid-body dynamics library using recursive algorithms. Auton. Robots41, 495–511. 10.1007/s10514-016-9574-0

  • CrossRef
  • Google Scholar

• 20

  FluitR.AndersenM. S.KolkS.VerdonschotN.KoopmanH. F. J. M. (2014a). Prediction of ground reaction forces and moments during various activities of daily living. J. Biomech.47, 2321–2329. 10.1016/j.jbiomech.2014.04.030

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 21

  FluitR.AndersenM. S.VerdonschotN.KoopmanH. F. J. M. (2014b). Optimal inverse dynamic simulation of human gait. Gait Posture39S, S1–S141. 10.1016/j.gaitpost.2014.04.060

  • CrossRef
  • Google Scholar

• 22

  FohannoV.BegonM.LacoutureP.ColloudF. (2014). Estimating joint kinematics of a whole body chain model with closed-loop constraints. Multibody Syst. Dyn.31, 433–449. 10.1007/s11044-013-9366-7

  • CrossRef
  • Google Scholar

• 23

  FoxM. D.ReinboltJ. A.ÕunpuuS.DelpS. L. (2009). Mechanisms of improved knee flexion after rectus femoris transfer surgery. J. Biomech.42, 614–619. 10.1016/j.jbiomech.2008.12.007

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 24

  FreglyB. J.BesierT. F.LloydD. G.DelpS. L.BanksS. A.PandyM. G.et al. (2012). Grand challenge competition to predict in vivo knee loads. J. Orthop. Res. Off. Publ. Orthop. Res. Soc.30, 503–513. 10.1002/jor.22023

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 25

  GaudetS.TremblayJ.Dal MasoF. (2018). Evolution of muscular fatigue in periscapular and rotator cuff muscles during isokinetic shoulder rotations. J. Sports Sci.36, 2121–2128. 10.1080/02640414.2018.1440513

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 26

  GirouxM.MoissenetF.DumasR. (2013). EMG-based validation of musculo-skeletal models for gait analysis. Comput. Methods Biomech. Biomed. Engin.16, 152–154. 10.1080/10255842.2013.815878

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 27

  GoughM.ShortlandA. P. (2008). Can clinical gait analysis guide the management of ambulant children with bilateral spastic cerebral palsy?J. Pediatr. Orthop.28, 879–883. 10.1097/BPO.0b013e31818e197c

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 28

  HermensH. J.FreriksB.Disselhorst-KlugC.RauG. (2000). Development of recommendations for SEMG sensors and sensor placement procedures. J. Electromyogr. Kinesiol.10, 361–374. 10.1016/S1050-6411(00)00027-4

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 29

  HicksJ. L.UchidaT. K.SethA.RajagopalA.DelpS. L. (2015). Is my model good enough? Best practices for verification and validation of musculoskeletal models and simulations of movement. J. Biomech. Eng.137:020905. 10.1115/1.4029304

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 30

  JansenK.De GrooteF.AertsW.De SchutterJ.DuysensJ.JonkersI. (2014). Altering length and velocity feedback during a neuro-musculoskeletal simulation of normal gait contributes to hemiparetic gait characteristics. J. Neuroeng. Rehabil.11:78. 10.1186/1743-0003-11-78

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 31

  JungY.JungM.RyuJ.YoonS.ParkS.-K.KooS. (2016). Dynamically adjustable foot-ground contact model to estimate ground reaction force during walking and running. Gait Posture45, 62–68. 10.1016/j.gaitpost.2016.01.005

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 32

  LeineweberD. B.BauerI.BockH. G.SchlöderJ. P. (2003). An efficient multiple shooting based reduced SQP strategy for large-scale dynamic process optimization. Part 1: theoretical aspects. Comput. Chem. Eng.27, 157–166. 10.1016/S0098-1354(02)00158-8

  • CrossRef
  • Google Scholar

• 33

  LloydD. G.BesierT. F. (2003). An EMG-driven musculoskeletal model to estimate muscle forces and knee joint moments in vivo. J. Biomech.36, 765–776. 10.1016/S0021-9290(03)00010-1

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 34

  McGinleyJ. L.BakerR.WolfeR.MorrisM. E. (2009). The reliability of three-dimensional kinematic gait measurements: A systematic review. Gait Posture29, 360–369. 10.1016/j.gaitpost.2008.09.003

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 35

  ModeneseL.GopalakrishnanA.PhillipsA. T. M. (2013). Application of a falsification strategy to a musculoskeletal model of the lower limb and accuracy of the predicted hip contact force vector. J. Biomech.46, 1193–1200. 10.1016/j.jbiomech.2012.11.045

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 36

  MoissenetF.ChèzeL.DumasR. (2014). Introduction of a set of EMG-based muscular activations in a multi-objective optimisation when solving the muscular redundancy problem during gait. Comput. Methods Biomech. Biomed. Engin.17, 132–133. 10.1080/10255842.2014.931542

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 37

  PengY.ZhangZ.GaoY.ChenZ.XinH.ZhangQ.et al. (2018). Concurrent prediction of ground reaction forces and moments and tibiofemoral contact forces during walking using musculoskeletal modeling. Med. Eng. Phys.52, 31–40. 10.1016/j.medengphy.2017.11.008

  • CrossRef
  • Google Scholar

• 38

  PerryD. J.BurnfieldD. J. (2010). Gait Analysis: Normal and Pathological Function. Thorofare, NJ: SLACK Inc.

  • Google Scholar

• 39

  PizzolatoC.LloydD. G.SartoriM.CeseracciuE.BesierT. F.FreglyB. J.et al. (2015). CEINMS: A toolbox to investigate the influence of different neural control solutions on the prediction of muscle excitation and joint moments during dynamic motor tasks. J. Biomech.48, 3929–3936. 10.1016/j.jbiomech.2015.09.021

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 40

  RajagopalA.DembiaC. L.DeMersM. S.DelpD. D.HicksJ. L.DelpS. L. (2016). Full-body musculoskeletal model for muscle-driven simulation of human gait. IEEE Trans. Biomed. Eng.63, 2068–2079. 10.1109/TBME.2016.2586891

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 41

  RomkesJ.SchweizerK. (2015). Immediate effects of unilateral restricted ankle motion on gait kinematics in healthy subjects. Gait Posture41, 835–840. 10.1016/j.gaitpost.2015.02.015

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 42

  SantosG. F.GomesA. A.SaccoI. C. N.AckermannM. (2017). Predictive simulation of diabetic gait: Individual contribution of ankle stiffness and muscle weakening. Gait Posture58, 208–213. 10.1016/j.gaitpost.2017.07.124

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 43

  SartoriM.FarinaD.LloydD. G. (2014). Hybrid neuromusculoskeletal modeling to best track joint moments using a balance between muscle excitations derived from electromyograms and optimization. J. Biomech.47, 3613–3621. 10.1016/j.jbiomech.2014.10.009

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 44

  SartoriM.ReggianiM.van den BogertA. J.LloydD. G. (2012). Estimation of musculotendon kinematics in large musculoskeletal models using multidimensional B-splines. J. Biomech.45, 595–601. 10.1016/j.jbiomech.2011.10.040

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 45

  ShaoQ.BassettD. N.ManalK.BuchananT. S. (2009). An EMG-driven model to estimate muscle forces and joint moments in stroke patients. Comput. Biol. Med.39, 1083–1088. 10.1016/j.compbiomed.2009.09.002

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 46

  SonJ.HwangS.KimY. (2012). A hybrid static optimisation method to estimate muscle forces during muscle co-activation. Comput. Methods Biomech. Biomed. Engin.15, 249–254. 10.1080/10255842.2010.522187

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 47

  ThelenD. G.AndersonF. C.DelpS. L. (2003). Generating dynamic simulations of movement using computed muscle control. J. Biomech.36, 321–328. 10.1016/S0021-9290(02)00432-3

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 48

  ThompsonJ. A.ChaudhariA. M. W.SchmittL. C.BestT. M.SistonR. A. (2013). Gluteus maximus and soleus compensate for simulated quadriceps atrophy and activation failure during walking. J. Biomech.46, 2165–2172. 10.1016/j.jbiomech.2013.06.033

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 49

  van der KrogtM. M.DelpS. L.SchwartzM. H. (2012). How robust is human gait to muscle weakness?Gait Posture36, 113–119. 10.1016/j.gaitpost.2012.01.017

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 50

  van der KrogtM. M.DoorenboschC. A. M.BecherJ. G.HarlaarJ. (2009). Walking speed modifies spasticity effects in gastrocnemius and soleus in cerebral palsy gait. Clin. Biomech.24, 422–428. 10.1016/j.clinbiomech.2009.02.006

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 51

  van SchieC. H. M. (2008). Neuropathy: mobility and quality of life. Diabetes Metab. Res. Rev.24(Suppl 1), S45–51. 10.1002/dmrr.856

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 52

  van Sint JanS. (2007). Color Atlas of Skeletal Landmark Definitions. Edinburgh; New York, NY: Churchill Livingstone/Elsevier. 1623.

  • Google Scholar

• 53

  WrenT. A. L.GortonG. E.ÕunpuuS.TuckerC. A. (2011). Efficacy of clinical gait analysis: a systematic review. Gait Posture34, 149–153. 10.1016/j.gaitpost.2011.03.027

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 54

  ZajacF. E. (1989). Muscle and tendon: properties, models, scaling, and application to biomechanics and motor control. Crit. Rev. Biomed. Eng.17, 359–411.

  • Pubmed Abstract
  • Google Scholar