Introduction

Researchers have used experimental data collected in a gait laboratory to study a wide range of movement impairments, including those caused by stroke [1], pelvic sarcoma [2], knee osteoarthritis [3], cerebral palsy [4], Parkinson’s disease [5], lower limb amputation [6], and incomplete spinal cord injury [7]. Commonly collected data include motion capture and ground reaction data, which have enabled the calculation of joint motions and joint loads derived from those quantities. However, there is growing interest in taking experimental data collection outside the laboratory without the need for force plates. Furthermore, eliminating the need for force plates is desirable inside the laboratory as well. The experimental hardware typically used in a gait laboratory is often expensive [3], people perform activities of daily living differently in a laboratory setting compared to in their daily lives [8, 9], and patients could receive objective feedback on therapies performed outside a clinical setting [10, 11].

Recent research advances are making it possible to take human movement measurement outside the laboratory, though significant limitations remain. Significant progress in hardware and software development (e.g., IMU systems [12] and markerless video motion capture [3, 13]) has made joint angle data available “in the wild.” Unfortunately, ground reaction data are more difficult to collect outside the laboratory. While instrumented shoe insoles exist [14,15,16,17,18], they do not capture all components of ground reaction forces and moments, plus most remain confined to a laboratory setting. Instrumented shoes have also been developed that can measure ground reaction forces and center of pressure, but the design influences subjects’ walking motions [19]. Without reliable ground reaction data, joint loads cannot be calculated via inverse dynamics for motions involving contact with the ground.

Researchers are seeking to overcome these limitations by using computational methods to predict ground reactions and joint loads accurately and reliably for walking from motion data alone. Artificial neural networks have been used to predict joint moments directly from more easily measured quantities like motion and electromyography (EMG) data [20]. Similarly, machine learning methods have been used to predict ground reaction forces from kinematic measurements [21, 22]. Physics-based full-body tracking optimizations employing deformable foot-ground contact (FGC) models have also been used to predict ground reactions and joint moments from measured motions [3]. These simulations typically use Hertzian elastic contact models with Hunt–Crossley nonlinear damping [23,24,25] to estimate ground reaction forces and moments between the feet and the ground. Physics-based simulations are advantageous for their interpretability, utilization of validated physical laws, and ability to extrapolate to physical situations where no data are available for training machine learning models.

While the contact force equations used by deformable FGC models are well established, the parameter values used in these equations are almost always generic [26] and not validated to reproduce experimentally measured ground reactions from specific subjects. Furthermore, simulation studies that use generic FGC models almost never report ground reaction moment information in any form (i.e., three time-varying moments calculated about a fixed point, a time-varying vertical moment and a time-varying center of pressure location, or an array of three time-varying forces distributed across the bottom of the foot). However, ground reaction moment information is critical for calculating joint moments during walking, and no researcher would think of throwing away ground reaction moment information when performing an inverse dynamics analysis for walking. Consequently, inverse dynamics joint loads of clinical interest estimated using generic FGC models may be unrealistic [3].

A likely explanation for unrealistic joint load results is the lack of foot-ground contact models calibrated to subject data. Historically, calibrating such models to actual data has been challenging because of the large number of parameter values that must be tuned (e.g., stiffness, damping, and friction) [27]. However, the Neuromusculoskeletal Modeling (NMSM) Pipeline [28] now makes it possible to calibrate personalized foot-ground contact models quickly and easily using experimental ground reaction and foot motion data. This automated calibration process removes the need for manual tuning of contact parameters [27]. When used in tracking optimizations of walking, FGC models calibrated using this software have closely reproduced experimental measurements of all six ground reaction components [1, 28]. Personalized FGC models have never been more accessible, but how they compare to widely used generic FGC models [26] for predicting ground reactions and inverse dynamics joint moments from kinematic data alone has not yet been evaluated.

This study compared the predictive capabilities of generic versus personalized deformable foot-ground contact models using two approaches. The first approach performed calibration optimizations that tracked experimentally measured ground reaction and foot kinematic data to test the models’ ability to reproduce both types of data simultaneously. The second approach performed tracking optimizations that estimated ground reaction forces and moments and lower-extremity inverse dynamics joint moments during walking from kinematic data alone while solving for dynamic consistency. For a robust comparison, the study was performed using experimental video motion and ground reaction walking data collected from three subjects in different states of health. Personalized FGC models have potential research applications inside and outside the laboratory and may help clinicians design better treatments for various movement impairments and assess treatment outcomes.

Materials and Methods

Given the complexity of the computational methods performed in this study, a flowchart detailing the interaction of the various computational tasks is provided in Supplementary Figure 1.

Experimental Data Collection

This modeling study used previously collected walking data obtained from three subjects. Subject 1 was a healthy male (age 47 years, height 1.7 m, mass 66.7 kg) [29]. Subject 2 was a male stroke survivor (age 79 years, LE Fugl-Meyer Motor Assessment 32/34 pts, right-sided hemiparesis, height 1.7 m, mass 80.5 kg) [1]. Subject 3 was a male with a right-side pelvic sarcoma who had not yet received surgery (age 46 years, height 1.73 m, mass 82.5 kg) [2]. All three subjects provided written informed consent, and data collection was approved by the institutional review boards of the University of Florida, Rice University, the University of Texas Health Science Center, and MD Anderson Cancer Center.

Although walking data were collected from the three subjects at different times and in different labs, the data collection protocols were nearly identical, supporting the use of these three data sets for the present study. All datasets included marker-based motion capture and ground reaction data. Motion capture data were collected using either a Vicon system (Vicon Corp., Oxford, UK) for Subjects 1 and 2 or a Qualisys system (Qualisys AB, Gothenburg, Sweden) for Subject 3, where marker placement utilized a modified Cleveland Clinic full-body marker set with markers added to the feet [30]. Ground reaction data were collected using a split-belt instrumented treadmill (Bertec Corp., Columbus, OH, USA) with belts tied for all three subjects walking at a self-selected speed (1.4 m/s for Subject 1, 0.5 m/s for Subject 2, and 1.0 m/s for Subject 3). One representative gait cycle with high periodicity starting with right heel strike was used for each subject. Subjects also performed a standing static pose to support skeletal model scaling.

Skeletal Model Development

For each subject, a skeletal model with personalized lower-body joints was created starting from a generic full-body OpenSim [31] skeletal model [28] with 31 DOFs. The generic model possessed a 6 DOF ground-to-pelvis joint, two 3 DOF hip joints, two 1 DOF knee joints, two 1 DOF ankle joints, two 1 DOF subtalar joints, two 1 DOF toe joints linking the hindfoot and toes, a 3 DOF back joint, two 3 DOF shoulder joints, and two 1 DOF elbow joints. The model also included forearm pronation and knee adduction coordinates, though these coordinates were locked to values consistent with the standing static pose. The models were scaled based on the static pose motion capture data using the OpenSim Scale Model tool.

The Joint Model Personalization (JMP) Tool in the NMSM Pipeline [28] was used to personalize select pelvis and lower-body joint functional axes, body scale factors, and marker locations using each subject’s walking data, reducing inverse kinematics (IK) marker distance errors for gait motion data. Knee, ankle, and subtalar joint axis parameters, pelvis, femur, and tibia scaling, and femur and tibia marker positions were optimized. Joint angles were calculated by the OpenSim Inverse Kinematics tool for each subject’s gait motion using experimental marker trajectories and skeletal models updated by JMP. The resulting motion data were filtered with a fourth-order zero phase lag Butterworth filter with a cutoff frequency of 7/tf Hz, where tf (seconds) equals the duration of the gait cycle for each subject [32].

Calibration Optimizations

Calibration optimizations were performed for both the generic and automatically personalized FGC models using the Ground Contact Model Personalization (GCP) Tool within the NMSM Pipeline [28]. These optimizations were used to evaluate how well both models could reproduce experimental ground reaction and foot motion data while tracking both types of input data. Prior to use in calibration optimizations, ground reaction data were filtered with a second-order zero phase lag Butterworth filter with a cutoff frequency of 7/tf Hz. For computational speed, GCP automatically extracts the two feet in each full-body OpenSim skeletal model and converts them into a feet-only OpenSim skeletal model, where hindfoot kinematics for each foot are redefined with a 6 DOF ground-to-hindfoot joint. This step allows GCP to make slight modifications to foot kinematics directly, avoiding the need to model the lower limb during GCP optimizations. For each type of contact model, the GCP cost function penalized errors in tracking experimental ground reaction forces and moments as well as foot kinematics.

For personalized FGC models, the GCP tool evaluated how well measured ground reaction and foot kinematic quantities could be tracked simultaneously while automatically creating personalized contact geometry. A grid of viscoelastic contact elements with linear springs and nonlinear damping [24] was placed under each hindfoot and toes body using 37 locations within the bounds of each foot. The GCP tool determined the values of physical contact model parameters during these runs, including stiffness coefficients unique to each contact element and damping, dynamic friction, and resting spring length values common to all contact elements. As the personalized models had varying stiffness coefficients, the cost function additionally penalized Gaussian-weighted deviations in stiffness coefficients to ensure that nearby contact elements had similar stiffnesses without physically unrealistic discontinuities. Although stiffness coefficients were element-specific, they were shared between sides for each subject so that the contact models were symmetric. An additional cost function term with low weight was used to enforce near-periodic motions. All GCP cost function terms are listed in Supplementary Table 1. The optimization also adjusted the planar position of each force plate in the lab coordinate system to promote consistency between kinematic and ground reaction data [33, 34], addressing a common source of error with split-belt instrumented treadmills [35]. Ground reaction data with updated force plate positions were used to calculate joint moments with the OpenSim Inverse Dynamics (ID) tool. These joint moments were filtered with the same process as IK data, and they were used for evaluating the quality of predicted joint moments from simulations in future steps.

For generic FGC models, the GCP tool evaluated how well measured ground reaction and foot kinematic quantities could be tracked simultaneously while using generic contact geometry. Instead of automatically adding a grid of springs, models were manually created based on a generic contact model (visualized in Supplementary Figure 2) with six sparsely, nonuniformly placed spherical viscoelastic contact elements with nonlinear springs and nonlinear damping [24] commonly used in OpenSim [26]. Physical parameters of the generic model were not optimized; only minimized kinematic adjustments were allowed. These tests used the same ground reaction and kinematic tracking terms as the personalization process. The ground reaction data included the updated force plate positions from the personalization process for consistency. Results were quantified using root mean square error (RMSE) calculated between GCP results and experimental data, with smaller errors representing more accurate reproduction of measured quantities.

Tracking Optimizations

Tracking optimizations were performed for both the generic and automatically personalized FGC models using the Tracking Optimization (TO) Tool within the NMSM Pipeline [28]. For each subject, these optimizations were used to evaluate the ability of each type of contact model to predict ground reaction forces and moments while tracking only kinematic data and solving for a dynamically consistent walking motion. Coordinate positions, velocities, and accelerations were included as design variables, and torque controls were used for hip, knee, ankle, and subtalar coordinates. The cost function tracked only coordinate positions and velocities for all joints except the ground-to-pelvis joint, leaving ground reactions and joint moments to be predicted during the optimization. The constraints prevented significant changes to initial coordinate positions and velocities. It also imposed dynamic consistency by keeping torque controls within 0.1 Nm of their respective ID joint loads during the optimization and by constraining the residual ID loads of the 6 DOF ground-to-pelvis joint to be less than 1 N or 0.1 Nm in magnitude, consistent with past work [1, 28]. All TO cost function and constraint terms are listed in Supplementary Tables 2 and 3. As with the calibration optimizations, results were quantified using root mean square error (RMSE) calculated between TO results and experimental data.

Results

Calibration Optimizations

For all three subjects, GCP tracking errors were lowest when using personalized instead of generic FGC models (Table 1). Overall, kinematic and ground reaction quantities were roughly 2–3 times more accurate with personalized models. RMS errors in ground reaction forces and moments and in foot translations and rotations were generally 40 to 60% lower when using personalized rather than generic FGC models (Supplementary Table 4). As a visual example, tracked quantities and stiffness coefficients for the right foot of Subject 1 are plotted in Figures 1 and 2. For generic FGC models, GCP runs took 1.83 ± 0.47 minutes (mean ± standard deviation), while for personalized FGC models, GCP runs took 9.87 ± 0.59 minutes.

Table 1 Root mean square error (RMSE) (mean ± standard deviation) in GCP tracked quantities on the right and left sides of all subjects for the calibration optimizations
Fig. 1
Fig. 1
Full size image

Right foot reference ground reactions and GCP calibration optimization results with generic and personalized FGC models for Subject 1. The distribution of stiffness coefficients on the personalized model is shown on the right

Fig. 2
Fig. 2
Full size image

Right foot reference kinematics and GCP calibration optimization results with generic and personalized FGC models for Subject 1

Tracking Optimizations

For all three subjects, TO tracking and prediction errors were lowest when using personalized instead of generic FGC models (Table 2). Motion tracking was roughly 4–5 times more accurate with personalized models, while ground reaction predictions were roughly 3–4 times more accurate and inverse dynamics load predictions were roughly 2–3 times more accurate. RMS errors in ground reaction forces and moments and in lower-body joint angles and joint moments were generally 60 to 80% lower when using personalized rather than generic FGC models (Supplementary Table 5). As a visual example, tracked motions, predicted ground reactions, and predicted joint loads for Subject 1 are shown in Figures 3, 4, and 5, respectively. Figure 6 shows the positions of skeletal geometry over the course of the gait cycle for Subject 1 for the nominal motion alongside motions predicted using each contact model.

Table 2 Root mean square error (RMSE) (mean ± standard deviation) in tracked and predicted quantities on the right and left sides of all subjects for the Tracking Optimizations
Fig. 3
Fig. 3
Full size image

Right foot reference ground reactions and Tracking Optimization results with generic and personalized FGC models for Subject 1

Fig. 4
Fig. 4
Full size image

Right leg reference kinematics and Tracking Optimization results with generic and personalized FGC models for Subject 1

Fig. 5
Fig. 5
Full size image

Right leg reference joint loads and Tracking Optimization results with generic and personalized FGC models for Subject 1. No reference load is plotted for the toes because the reference ID results did not have ground reactions distributed between the hindfoot and toes bodies, so there is no meaningful reference for this joint load

Fig. 6
Fig. 6
Full size image

Animation strip of reference and resulting motions from Tracking Optimization for Subject 1. Note the differences in foot positions at 25 and 75% and in trunk orientation at 100%

When comparing TO results from the three subjects, joint angle and ground reaction errors for the generic and personalized FGC models exhibited similar trends across all three subjects. Overall, Subject 2 had the best TO results and Subject 3 the worst. Hip flexion angle, especially for Subject 3, exhibited the largest reduction in kinematic tracking errors when using personalized FGC models. Using generic FGC models, TO runs took 37.9 ± 9.44 minutes, while using personalized FGC models, they took 32.1 ± 16.2 minutes.

Discussion

This study compared the predictive abilities of generic and personalized deformable FGC models for estimating ground reactions and joint loads and modeling dynamically consistent motions. The three subjects included provided a wide variety of test conditions, showing FGC model performance for a healthy subject, a stroke survivor, and a subject with a pelvic sarcoma. For all three subjects, the calibration and tracking optimization runs successfully converged for both generic and personalized FGC models. In calibration results, the personalized FGC models had lower RMS errors in all tracked quantities than did the generic models except in hindfoot Z translation, which differed on average by less than a millimeter. The improved performance of personalized models in calibration optimizations is an intuitive result since the physical parameter values were adjusted in the personalized models to reduce errors. However, it is more notable that the personalized models performed significantly better in tracking optimizations, greatly reducing errors in ground reaction and ID joint load predictions and also reducing errors in tracked joint angles when solving for dynamically consistent walking motions. These results suggest that although both generic and personalized deformable FGC models can produce reasonable ground reactions from walking motions and can support dynamically consistent walking simulations, personalized FGC models produce better tracking optimization results and have the potential to produce significantly more accurate predictive simulations for individual subjects.

While performing this study, we made additional observations that may provide further insight into the use of deformable FGC models for tracking optimizations and predictive simulations. First, FGC model personalization can reduce tracking optimization computation time. For all subjects, TO simulations using personalized FGC models converged in fewer iterations than simulations using generic FGC models. For two of the subjects, tracking optimizations with personalized FGC models also converged in less time, saving more time than was required to calibrate their contact models. For the third subject, the tracking optimization with a personalized FGC model converged in a comparable amount of time. A speed increase in tracking optimizations would be especially beneficial for users intending to run multiple predictive simulations for the same subject. Second, a TO problem formulation more consistent with past work [3] produced unrealistic toe motions. This formulation included tracking joint accelerations and removing tracking from the toes coordinate to let it change freely. Acceleration tracking yielded slight reductions in ground reaction tracking errors, though joint angles were slightly less accurate. Allowing the toes coordinate to move freely also slightly (though less consistently) reduced ground reaction tracking errors overall. However, the toes moved in a highly unrealistic fashion, attempting to flex and extend by more than 20 degrees from their nominal motions.

The performance of our generic FGC model in tracking optimizations is consistent with past work, though it is important to note that the weighting of cost function terms is key to how errors in simulation results are distributed. For simulations predicting ground reactions, it is possible to have higher accuracy in ground reaction and joint load predictions at the cost of less accurate motion tracking [3, 36]. As seen in the present study, irregular shapes may appear when plotting ground reactions over time from results generated using generic FGC models [37], though these plots are not always available in published work. Similar errors are observed in other motions, such as sit-to-stand simulations [3, 38]. It would additionally be useful to compare the results of this study’s ground reaction moment predictions to past work. Use of personalized FGC models yielded dramatic improvements in predicted ground reaction moments, but they are often not reported, even for studies using three-dimensional models.

Because the generic and personalized FGC models used in this study possessed different numbers and distributions of contact elements, it is unclear whether personalized FGC models performed better due to personalization or number/distribution of contact elements. To explore this issue, we repeated GCP and TO tool runs using (1) a dense generic model and (2) a sparse personalized model. The dense generic model used contact element locations from the personalized model, with generic stiffness values scaled down to match the increased number of contact elements, while the sparse personalized model used the generic model contact element locations. Some GCP results improved while others got worse when using the dense element grid with the generic model, while all GCP results got worse when using the sparse element grid with the personalized model (Supplementary Table 6). Notably, the generic model still did not approach the accuracy of the personalized model simply by using the dense element grid. In contrast, TO would not converge when using the dense element grid with the generic model, while TO converged to a less accurate solution when using the sparse element grid with the personalized model (Supplementary Table 7). These findings indicate that simply increasing the contact element density in the generic model will not yield a more successful model, while the personalized model can scale to different grid densities. However, a denser grid of contact elements is still needed to predict ground reactions and joint moments reliably. Thus, returning to the original issue, within the context of predictive simulations of walking, (1) the number of contact elements is important—too few and ground reactions cannot be predicted accurately; and (2) personalization of contact elements is also important—a larger number of unpersonalized contact elements still does not predict ground reactions reliably.

It is difficult to assess the clinical significance of the accuracy improvements produced by personalized FGC models compared to generic FGC models, particularly for tracking optimizations. These optimizations also provide an indication of how each FGC model would perform for predictive simulations of a patient’s post-treatment walking function. Clearly, sixty to 80% reductions in ground reaction forces and moments and in lower-body joint angles and joint moments can only be beneficial when seeking to predict post-treatment walking function for a particular patient. For the generic FGC model, inverse dynamics joint moment errors produced by TO were comparable to or larger than the errors one might expect from using a lower-extremity skeletal model with non-personalized joint centers and functional axes [39]. Furthermore, depending on the clinical application, even small improvements in tracking errors may improve the clinical utility of predictive simulations. For example, predicted load magnitudes that differ by as little as a few percent can greatly change the number of cycles that tissue can endure before failing [40].

The need to eliminate residual ID loads in the 6 DOF ground-to-pelvis joint can be explained by the approach we used to generate Tracking Optimization solutions. Within a direct collocation optimal control framework, predictive simulations of human movement can be generated using either forward dynamic or inverse dynamic simulation methods [25]. At the level of skeletal dynamics, forward dynamic methods take force and torque controls as inputs and calculate joint motions as outputs, whereas inverse dynamic methods take joint motions as inputs and calculate force and torque controls as outputs. Either way, motions and controls are predicted. Studies have demonstrated that for human movement problems, inverse dynamic methods tend to converge more reliably in fewer iterations and with less computation time than do forward dynamic methods [41, 42]. For this reason, the NMSM Pipeline utilizes an inverse dynamic method to perform predictive simulations [28]. However, for some human movement situations, an inverse dynamic method requires elimination of “residual loads” not present with a forward dynamic method. Residual loads are non-physical external loads acting on the model segment connected to ground (termed the “root segment”) and account for dynamic inconsistencies between the motion of the model and the external loads applied to the model. When an inverse dynamic analysis was performed on one of our OpenSim walking models, OpenSim calculated inverse dynamic loads for all joints in the model, including the 6 DOF ground-to-pelvis joint. Since no external loads were applied to the pelvis (the root segment) during our human subject experiments, the 6 non-zero ID loads (3 forces and 3 torques) calculated for the ground-to-pelvis joint represented residual loads required to balance the dynamic equations. Thus, for our Tracking Optimization solutions to be dynamically consistent, they had to drive these 6 residual loads to zero to within some specified tolerance. This requirement was achieved through the root segment residual load constraint used in our Tracking Optimization problem formulation.

Though relatively large differences in calibration and tracking optimization errors between generic and personalized FGC models make the conclusions of this study clear, this work has limitations that are important to consider when applying these findings. First, the data used in this study represent a best-case scenario for both calibration and prediction. For all three subjects, the motion capture data came from a multi-camera marker-based system, which is currently the “gold standard” in the biomechanics field and the most accurate motion capture method without using ionizing radiation [43, 44]. For reasons such as budget and the ability to measure motion outside limited laboratory settings, researchers are becoming increasingly interested in other motion capture methods such as IMU and markerless video motion capture [3, 12, 13]. Both the GCP and TO tools are compatible with joint angle motion data derived from any motion capture method, so this study and its procedures can be applied to these motion capture methods as well. However, lower-quality motion capture data will likely lead to a lower-quality model calibration, so it is unknown how well an FGC model personalized using these data would perform. The TO tool would still need to resolve any dynamic inconsistency in the captured motion, so it is likely that using either a generic or a personalized FGC model with lower-quality motion data would produce higher predictive and tracking errors than those found in this study. Second, all simulations in this study were torque-driven, so muscles and their limitations on possible motions and control were not accounted for. The inclusion of muscles can be especially important for predicting new motions [1], so it will be valuable to confirm that the trends in the results of this study remain valid when using more complete musculoskeletal models. Third, FGC models were used in this study to predict dynamically consistent versions of the same motions used in calibration. It would be useful to evaluate the predictive capabilities of FGC models for new situations, such as walking at a faster speed or running. Fourth, FGC models were calibrated and tested using a single gait cycle for each subject. Initial testing during NMSM Pipeline development suggested that while use of multiple gait cycles to calibrate FGC models is possible, it greatly increased the required computation time while yielding FGC model parameter values that were similar to those found when calibrating with only a single gait cycle. Furthermore, FGC models calibrated using only one gait cycle also worked well for other gait cycles from the same subject. We do not expect that using more gait cycles would significantly impact our results, but it would be valuable for a future study to quantify the influence of additional calibration data. Fifth, although our personalized FGC models have parameters with physical meaning, we have not evaluated our calibrated parameter values against corresponding experimental measurements of stiffness, damping, and friction. Some model parameter values, such as spring stiffness and resting length, may be non-unique, with multiple potential solutions possessing low cost function values. Future work should compare calibrated FGC model parameter values with experimental measurements of these quantities. Sixth, the simulation results presented in this study may represent local but not global optima. We used physically realistic initial guesses for our optimizations to keep results close to where the global optimum is expected to be, but more thorough testing with different initial guesses could reveal better, previously unknown solutions. Seventh, only three subjects were included in this study. This limitation is the result of the limited availability of subject data and the requirements of the computational procedures. Additional subjects would enable statistical testing of our results. Eighth, although our three subjects were tested using similar protocols, inconsistencies between motion capture systems from different manufacturers in different laboratories may have led to individual subjects having different errors in their data.

In conclusion, we found that personalized FGC models predict dynamically consistent walking motions more accurately than do generic FGC models for subjects in a variety of states of health. The study demonstrates that the automatic personalization process for FGC models is robust, fast, and easy to perform for different subjects. Given the significant improvement in joint motion, ground reaction, and joint moment accuracy during tracking optimizations in particular, foot-ground contact model personalization appears to be worth the effort. Personalized FGC models may enable higher-quality predictive simulations in biomechanical walking studies of movement impairments in an array of settings.