Comparison of Existing Methods for Characterizing Bi-Linear Natural Ankle Quasi-Stiffness

Natural ankle quasi-stiffness (NAS), also known as dynamic ankle joint stiffness, is a quantifiable mechanical property that is used to characterize the natural spring-like function of the human ankle [1–6] during the stance phase of gait. NAS is also used as design criteria for biomechanical assistive devices, such as ankle-foot orthoses (AFOs) [7–9]. NAS has historically been defined as the slope of a linear best-fit regression of sagittal ankle moment versus sagittal ankle angle [1] during the loading phase of stance (Fig. 1). The loading phase is characterized kinematically by the ankle dorsiflexing while the foot is flat on the ground, and kinetically by an increasing plantarflexion moment about the ankle that resists dorsiflexion. Many models of NAS often approach the loading phase as a single phase of stance with a single NAS value [1–6], though more complex models have been developed to better approximate the actual moment-angle behavior of the ankle that is inherently nonlinear [10–12]. These refined models provide better fits to measured data but have not been widely adopted. This lack of adoption of these more complex models is perhaps due to inconsistent methods and/or an absence of intuitive clarity or translatability of model outputs to physical devices at the time.

Bi-linear NAS (BL-NAS) is one such type of NAS model with increased complexity. BL-NAS attempts to capture the inherent nonlinearity of NAS by dividing the loading phase into two subphases—early loading (EL) and late-loading (LL) [11,12]. BL-NAS separates regions with potentially different stiffness properties and offers increased model complexity while retaining intuitively meaningful results (a steeper slope means higher stiffness). While two prior studies have developed and used BL-NAS models [11,12], no comparison between these methods' advantages and disadvantages has been conducted. Understanding and comparing these models can provide important insight into the joint-level mechanics of the human ankle, and how to adapt those mechanics to design better assistive ankle-foot devices. If these models provide better approximations of actual NAS data, then devices with these characteristics could provide increased benefit to users compared to current devices with single stiffness properties that may not reflect healthy NAS. This discrepancy between the natural stiffness behavior of the ankle and the stiffness behavior of a device might result in individual fighting against the device at worst, or otherwise experiencing nonoptimal stiffness assistance at best.

This study aimed to compare two BL-NAS models to each other, as well as to the standard single linear NAS method (SL-NAS). The two models assessed in this study were described by Shamaei et al. [11] (“Shamaei” model) and Crenna and Frigo [12] (“CF” model). Key outcome measures for the comparisons were EL and LL-NAS stiffness, root-mean-squared error (RMSE), and the sagittal ankle angles at which EL and LL begin, end, and transition from one to the other. EL and LL-NAS stiffness will identify if the two regions of loading have different stiffness properties and if a BL-NAS model is even appropriate. We hypothesize that EL and LL-NAS will be different within a given walking speed and that both will increase with speed. RMSE will identify which model more closely approximates measured data. We hypothesize that both BL-NAS models will have lower RMSE than standard single linear NAS (SL-NAS). Ankle angles at which EL and LL begin, end, and the transition will be used to assess timing and event differences between models. We hypothesize that the Shamaei and CF models will not differ in the ankle angles at which EL and LL begin, end, and transition.

Eleven young, healthy individuals (age = 23.6 (2.9 SD) years, height = 1.76 (0.10 SD) m, mass = 75.7 (16.3 SD) kg, 6 M/5F) with no history of neurological or musculoskeletal disease provided informed consent to participate in this Institutional Review Board approved the study. Subjects traversed a 20 m walkway while wearing a six degree-of-freedom marker set [13]. Kinetic data were captured at 1200 Hz by force plates (AMTI, Watertown, MA) embedded in the walkway. Kinematic data were captured at 120 Hz using a 12-camera passive marker-based motion capture system (Qualisys AB, Göteborg, Sweden). Subjects walked shod at three prescribed speeds: 0.6, 0.8, and 1.0 statures/second. For reference, the comfortable walking speed for healthy individuals is approximately 0.8 statures/s [14]. Statures/s was chosen to control speed because statures/s has been used in similar studies of locomotion to control for speed [10]. Speed was monitored using a real-time feedback program that displayed subjects' walking speed on a large monitor at the end of the walkway. Five good trials per limb per speed were captured, with a good trial defined as one in which the subject's entire foot was located on one and only one force plate while they walked within ±0.02 statures/s of the target speed. Footwear was not controlled to ensure that the subjects were not walking in any unfamiliar footwear that may influence their gait. However, all subjects wore conventional athletic/running shoes with typical midsole construction. No “minimal” or “barefoot” shoes were worn by any subjects, nor did any subjects wear overly cushioned or “maximal” shoes. No study to the authors' knowledge has systematically examined the effect of specific types of shoes on ankle joint stiffness during walking. However, footwear does seem to have minor effects on peak dorsiflexion angle but not peak ankle moment during running tasks [15], which may affect calculated ankle stiffness values.

After data collection, kinematic data with less than five consecutive frames missing were gap-filled using a cubic spline. Kinetic and kinematic data were filtered at 25 Hz and 6 Hz, respectively, using a fourth-order zero-lag Butterworth filter. Standard inverse dynamics, including sagittal ankle angle and sagittal ankle moment, were computed [13]. Ankle moments were normalized by subject mass as per International Society of Biomechanics recommendations [16], and all data were scaled in time to percent stance. For brevity, the terms “moment(s)” and “angle(s)” in this paper will always refer to sagittal ankle moments and sagittal ankle angles, respectively.

The Shamaei [11] and CF [12] models were then applied to our processed data. Though details of the models can be found in their respective papers, a brief outline of their methods follows. For the Shamaei model, EL (originally called “dorsi-flexion”) begins at the instant of initial ankle angle minimum at ∼5% stance and ends at the instant of the local minimum in the vertical ground reaction force at ∼50% stance. LL (originally called “dual-flexion”) begins at the end of EL and ends at the instant of peak dorsiflexion angle. The CF model differs in its definitions of the start and end of the EL and LL phases, originally called “early rising phase” and “late rising phase,” respectively. For the CF model, EL begins at the instant the plantarflexion moment reaches 5% of its peak moment and ends when the instantaneous slope (moment versus angle) between two data points exceeds 1.7 times the average slope up to that point. LL begins at the end of EL and ends when the dorsiflexion moment reaches 95% of the peak moment. The “transition” point, as we will refer to it, is the point in either model at which the EL phase ends and the LL phase subsequently begins.

One adjustment to the CF model was made for our study. If the 1.7× slope threshold was not met for the entire loading phase, or if it occurred within four frames of the end of the LL phase, the end of EL was forced to be located four frames from the end of LL, so that enough points in the LL phase existed to perform a linear regression. This was not mentioned in the original study, but it was necessary for our dataset. For both models, just like in their original respective studies, NAS values during each phase were calculated using a least-squares linear regression on each loading phase after the phases had been defined.

Outcome measures to compare the models for this study were EL and LL-NAS (N*m/deg/kg), EL and LL beginning, ending, and transition angles, and RMSE across the entire loading phase as defined by the respective model. RMSE was calculated by taking the squared error between the EL-NAS linear regression and measured data for the EL phase, and then the squared error between the LL-NAS linear regression and measured data during the LL phase—as defined by the respective model for that specific trial—then the square root of the combined average squared error across all data points in EL and LL was calculated. This resulted in a single value per trial that was directly compared to the SL-NAS RMSE. The SL-NAS model used for comparison defined the loading phase from the initial ankle angle minimum to peak dorsiflexion angle, based on Davis and De-Luca [1]. To note, SL-NAS beginning and ending definitions are the same as EL beginning and LL ending for the Shamaei model because it covers a wider range of data points per trial than the CF model. We also performed a type of sensitivity analysis for the transition point of both models by artificially changing it by ±1 deg (referred to as “Shamaei±” and “CF±”) to determine what effect, if any, it had on EL and LL-NAS values and RMSE for both models. All outcome measures were tested for significant differences between models and between speeds within each model using Kruskal–Wallis tests for non-normally distributed data with initial significance set to p < 0.05 before a posthoc Dunn's test to account for multiplicity.

There were 110 total trials per speed across 11 subjects (five trials per limb per subject). Both the Shamaei and CF models exhibited significantly and substantially lower RMSE compared to SL-NAS (Figs. 1 and 2, Table 1) at all three speeds, indicating a better model fit and confirming our first hypothesis (p < 0.001). The CF model had a slightly lower RMSE than the Shamaei model at all speeds (p < 0.021), but a higher standard deviation of RMSE than the Shamaei model (Table 1 and Fig. 2). Both models exhibited significant increases in LL-NAS with gait speed (p < 0.001), though no change in EL-NAS with speed (p > 0.223) (Table 1). The Shamaei model exhibited significantly lower EL and LL-NAS at 0.6 and 0.8 statures/s walking speeds than the CF model (p < 0.422), but neither stiffness value differed between models at 1.0 statures/s (p > 0.206) (Table 1). The beginning of EL occurred significantly later (higher dorsiflexion angle) at every speed for the CF model compared to the Shamaei model (p < 0.001). The transition angle was significantly higher for the CF model at the 0.6 statures/s speed only (p > 0.004). The angle at which LL ended was not significantly different between models (p > 0.380). These results provide partial support for our second hypothesis that the ankle angles at those events would not differ significantly between models.

Interestingly, the sensitivity analysis that was performed by manually shifting the transition point by ±1° resulted in slightly higher RMSE values for both models at all speeds, regardless of the direction of the shift. These differences were not significant at 0.6 statures/s (p > 0.245), but were at 0.8 and 1.0 statures/s (p < 0.013) (Table 1). The Shamaei± and CF± models generally did not result in any changes in EL or LL-NAS values at any gait speed (Table 1).

This study was the first to compare and contrast two existing models that define BL-NAS. Both BL-NAS models fit measured data better than the more commonly used SL-NAS model. Differences in how each BL-NAS model defined EL and LL phases resulted in some differences in stiffness values, phase beginning/ending, and the transition from EL to LL. While this study did not seek to make a definitive determination about which model is better or worse, we can examine some of the advantages and disadvantages of each and attempt to make distinctions between situations where one may be more appropriate than the other. These models could influence the design of ankle-foot orthotic and prosthetic devices to have stiffness properties that closely match a healthy and/or intact limb than currently available devices based on SL-NAS.

For any modeling study, the question often is: did one model better match actual measured data than another model? RMSE is a commonly used metric to assess model performance and was used in the present study. The CF model had slightly lower RMSE than the Shamaei model at all speeds, which were shown to be statistically significant. This would indicate that the CF model performed better than the Shamaei model. However, RMSE was calculated during the period from the beginning of EL to the end of LL for every stance phase as each model defined those events. Since each model defined those events differently, and each stance phase itself is different, the length of the period used to calculate RMSE was also different, as indicated by EL and LL angular excursion (Table 1). Most often, the Shamaei model encompassed a larger portion of the loading phase than the CF model, and this difference was statistically significant. This meant that some more complex, error-prone portions of the total NAS curve were simply not included in the CF model. This is best exemplified in Fig. 3(c), for which the CF model defined both EL and LL during the same linear part of the NAS curve and did not pick up on the bi-linear characteristics of the curve that occurred earlier. The Shamaei model, alternatively, covered the portion of that stance phase that did exhibit bi-linear behavior. For this specific curve (Fig. 3(c)), the RMSE for the Shamaei model was worse than the CF model only because the Shamaei model used more of the total curve. Mathematically the CF model appears to be better. However, when considering potential orthotic and prosthetic devices that are intended to better emulate natural ankle motion, including as much data from as many portions of stance will provide a more holistic approach to any future design.

Specifically related to how events are defined for the models, the main difference is the selection of the starting point for EL. The Shamaei model selects the initial minimum ankle angle, while the CF model selects the beginning of EL as the point at which the plantarflexion moment is 5% of its peak. Often, at the initial angle minimum, the ankle is generating a small net dorsiflexion moment (negative), and the CF model misses any period of EL for which this is the case. This often occurs when the ankle is in a slightly plantarflexed position (Fig. 1) and may explain the significant and substantial differences in average EL start angle, seen in Fig. 2.

On average, the transition point was generally not significantly different between models, except at 0.6 statures/s (Table 1 and Fig. 4). This selection of transition points, combined with the later EL start for the CF model, resulted in lower overall EL angular excursion for the CF model compared to the Shamaei model (Fig. 5). The CF model relies on the assumption that the transition point occurs when the instantaneous NAS slope exceeds 1.7 times the average slope up to that point. The justification given for this 1.7× threshold by the authors of that study was because it “well corresponded to the visual identification of a sharp increase of the slope” [12]. If a stance phase exhibits an obviously bi-linear NAS curve, this method works quite well (Fig. 3(a)), and both models choose transition points at approximately the same instant. However, for some stance phases, the CF transition point might not occur in a visually obvious location. Some stance phases exhibited loading phases with very smooth, steadily increasing plantarflexion moments and dorsiflexion angles, for which the instantaneous slope might not exceed the 1.7× threshold until the very end of the loading phase (Fig. 3(b)). In those cases, there was a long EL and short LL. Likewise, a sudden jump in ankle moment at the beginning of loading might result in the opposite situation—a short EL phase and a very long LL phase (Fig. 3(d)). Sometimes a stance phase exhibited a clearly bi-linear loading phase, but the visually obvious transition point occurred while the moment was still negative. In that case, the CF model instead placed the EL beginning in a nonobvious location (Fig. 3(c)). Most NAS curves resembled a mixture of Figs. 3(a) and 3(b). Figures 3(c) and 3(d) are outliers used to illustrate extremes of the CF model transition point as occurring later (c) or earlier (d) than average. This variability of the CF model in defining the EL and LL phases is contrasted by the relative consistency of the Shamaei model (Fig. 3). The Shamaei model automatically ensures that the transition point will occur nearly at 50% stance for healthy walking when the M-shaped ground reaction force (GRF) curve reaches a local minimum. This consistent transition point selection based on GRF may not be the case for some patient populations that exhibit a lot of heterogeneity, and who may not exhibit typical M-shaped GRFs during stance. Thus, the Shamaei model may not adapt as well as the CF model to nontypical ankle moment-angle profiles during stance and these results do not necessarily translate directly to patient populations. This represents a limitation of this study—that only young, healthy individuals were observed. These tradeoffs in consistency and adaptability of the models thus play a role in assessing which model to use, rather than just rote reliance on RMSE outputs (Fig. 6).

Regarding differences between EL and LL-NAS, the Shamaei model generally resulted in lower stiffness values compared to the CF model for both phases of loading. These differences were significant at 0.6 and 0.8 statures/s but were not at 1.0 statures/s. LL-NAS was shown to increase with gait speed, while EL-NAS did not (Table 1, Fig. 7). SL-NAS is known to increase with gait speed [5], and our results here support that finding (Table 1). When examined with a BL-NAS model, the increase in stiffness that is observed in an SL-NAS model can thus be attributed to late-loading changes, not necessarily behavior that occurs at the beginning of the loading phase.

Notably, this study did not directly compare higher-order polynomial regressions of the ankle moment-angle curve, as Nigro et al. have already done [10]. The reason why direct comparisons of BL-NAS models to polynomial regressions were not conducted is that this study only sought to identify differences between BL-NAS models and an SL-NAS model. BL and SL-NAS offer more intuitively meaningful and results than polynomial regressions since greater BL or SL-NAS regression coefficients immediately indicate a higher stiffness, and vice-versa. Higher-order polynomial regressions, in contrast, sometimes output negative coefficients, have squared or cubed terms, and other characteristics that obfuscate the overall concept of “stiffness.” Although BL-NAS RMSE may not be as low as, say, a cubic polynomial regression, the BL-NAS models make up for it in terms of intuitive clarity and translatability of results to physical devices.

The future directions of this study seek to translate the improved models of NAS that Crenna and Frigo [12] and Shamaei et al. [11] have created into physical devices. NAS is an easily translatable mechanical property to design passive spring-like devices that can assist impaired [17–19] or healthy [20] users, and translating these updated models to the current state-of-the-art of assistive devices is imperative to moving the field forward. An application of BL-NAS is to use these data to drive the design of novel ankle-foot assistive devices, such as AFOs. The concept of variable stiffness has been implemented into quasi-passive actuators for prosthetic ankle-feet [21–23] to closer match ankle-foot stiffness properties, but this concept has not yet been implemented into AFO design, to the authors' knowledge. Most passive-dynamic AFOs (PD-AFOs) consists of a foot plate or shell, a cuff, and a strut (posterior, anterior, or on the side). The strut is the main contributor to the stiffness properties of the PD-AFO. When the user dorsiflexes during the loading phase of stance, the strut bends and resists that motion. The stiffness of the strut dictates the amount of resistance to motion. This design can be effective for individuals with isolated plantarflexor weakness, often due to orthopedic or traumatic injury [17,24], or with neuromuscular impairments that result in lower-limb weakness, such as individuals poststroke [8]. However, tuning PD-AFO bending stiffness based on SL-NAS results in stiffness that is tuned to the entirety of the loading phase, which is inherently nonlinear [10]. Thus, SL-NAS may closely approximate true NAS during the middle of loading, but not at the beginning or end of loading. Arch & Reisman noticed that some individuals poststroke did not dorsiflex much when wearing a customized PD-AFO and posited that these individuals were unable to overcome the stiffness at the beginning of loading to initiate bending [8]. A solution to this problem would not simply be a matter of decreasing the PD-AFO bending stiffness overall, as such a stiffness would not provide enough support during the end of loading, and might result in excessive and/or uncontrolled dorsiflexion. AFOs with a single bending stiffness customized to the entire loading phase of an individual poststroke might be too stiff for the user to fully deform, and thus cannot attain proper peak dorsiflexion angle before push-off. If an AFO with low stiffness in early loading can be made, perhaps users may be able to initiate bending sooner, achieve a larger overall dorsiflexion excursion, and obtain limb kinematics that more closely approach those of a healthy population. Both BL-NAS models examined in this study indicate that AFOs with bi-linear stiffness should be less stiff during low dorsiflexion angles and more stiff during higher dorsiflexion angles to provide the right amount of support at the right time of stance. Comparing both BL-NAS models, the Shamaei model indicates that lower AFO bending stiffness during EL and LL might be appropriate compared to the CF model.

This study provides a thorough investigation of BL-NAS as a refined model for NAS when compared to SL-NAS. The two BL-NAS models examined both had substantially lower RMSE than the SL-NAS model. This improvement in modeling provides a justification for using the concept of BL-NAS to improve the designs of biomechanical assistive devices. Finally, and importantly, this analysis reveals the strengths and weaknesses of each model that are critical to consider when applying them to ankle-foot assistive device design.