A general relationship links gait mechanics and running ground reaction forces

Running ground reaction forces are of fundamental physical and biological importance. Acutely, they determine the body's state of motion, limb-loading rates and tissue stresses. Over time, they influence the health and functional status of the tissues, limb and runner. However, the mechanical basis of the vertical force–time waveform patterns described most extensively for human runners (Cavanagh, 1987; Munro et al., 1987) continues to be a matter of significant disagreement (Chi and Schmitt, 2005; Clark et al., 2014; Denoth, 1986; Derrick, 2004; Lieberman et al., 2010; Nigg, 2010; Shorten and Mientjes, 2011). The current discordance is at least partially attributable to the mechanical complexities of the limbs and bodies responsible for the forces present at the limb–ground interface. Human and other vertebrate runners are mechanically complex in their body and limb-segment morphology, tissue properties, neural control of muscle forces, and joint and limb stiffnesses. These features allow body mass components to accelerate differently with respect to one another and the ground. Because the total waveform corresponds to the acceleration of the body's entire mass, the summed accelerations of different mass components must somehow determine the instantaneous forces and waveform patterns observed.

At present, the most common approach to explaining the force–time waveform patterns of human runners are lumped element, spring–mass systems that include four or more hypothetical mass components with multiple springs and dashpots (Fig. 1A). Most current versions include 14 or more input variables derived via forward dynamics simulations (Chi and Schmitt, 2005; Liu and Nigg, 2000; Ly et al., 2010; Nigg and Liu, 1999; Nikooyan and Zadpoor, 2011; Zadpoor and Nikooyan, 2010). Per their intended purpose, these models are able to provide close, post facto fits to the rising edge of the force–time waveforms that result from rear-foot strike mechanics at jogging speeds under a variety of surface, footwear and other conditions (Ly et al., 2010; Zadpoor and Nikooyan, 2010). However, these models do not attempt to predict the falling edge of the waveform, they do not explain the differently shaped waveforms that typically result from fore-foot strike mechanics, and their ability to fit waveforms from intermediate and fast running speeds is completely unknown.

A scientific justification for including numerous mass components to account for vertical ground reaction force–time waveforms was importantly provided by a direct motion-to-force experiment conducted by Bobbert et al. (1991) a quarter century ago. These investigators demonstrated that the total waveform can indeed be reasonably predicted from the summed accelerations of seven body mass components at modest running speeds. Recently, we have theorized that an alternative approach may allow the number of masses needed for full waveform prediction to be reduced from seven to only two.

Our approach (Clark et al., 2014) divides the body's mass into two, anatomically based, invariant, mass components: the first corresponding to the mass of the lower limb (m₁, 8% total body mass) and the second corresponding to the remainder of the body's mass (m₂, 92% total body mass; Fig. 1B). The model theoretically allows the full vertical force–time waveform to be predicted from the force contributions corresponding to the two body mass components. Impulse 1 results from the vertical collision of m₁ with the running surface and impulse 2 results from the vertical accelerations of m₂ throughout the ground contact period (Clark et al., 2014). Our introductory effort (Clark et al., 2014) was able to account for essentially all of the variation present in four vertical force–time waveforms (R² range: 0.95–0.98, mean=0.97±0.01) selected for their amplitude, duration and shape heterogeneity. However, the close fits achieved resulted from post facto selection of the input parameters to maximize the goodness of each fit.

Here, we undertook a direct, experimental test of the hypothesis that the motion of two discrete body mass components is sufficient to predict running vertical ground reaction force–time waveforms in full. The two-mass model requires only body mass and three stride-specific measures as inputs: contact time, aerial time and lower limb acceleration. Because these three inputs can be readily acquired from a variety of video or other inexpensive motion-sensing technologies, our model potentially offers economical options for generating running ground reaction force waveforms without a force platform. Additionally, the concise running force–motion linkage provided could be applied to footwear, prosthetic and orthotic design, or used to inform gait interventions designed to reduce injuries or enhance running performance.

where m_b is body mass and g is gravitational acceleration (g=9.8 m s⁻²), and the quantity t_c+t_a equals the step time, t_step.

where m₁ is 8.0% of the body's mass (Plagenhoef et al., 1983; Winter, 1990), Δt₁ is the time interval between touchdown and the vertical velocity of m₁ slowing to zero, Δv₁ is the change in vertical velocity of m₁ during Δt₁, and F_1,avg is the average force during the total time interval (2Δt₁) of impulse J₁, here defined as the impact interval. A single ankle marker is used to measure Δv₁ and Δt₁, which determine the vertical acceleration of the lower limb mass m₁ (Fig. 2). The lower limb attains a relatively constant positive velocity after the impact interval (Fig. 2C), resulting in a near-zero acceleration of m₁ (Fig. 2D) and negligible force (Fig. 2E). Thus, J₁ can be represented by a finite impulse during the impact interval. Kinematic data for the ankle marker position during the impact interval for representative rear-foot strike (RFS) and fore-foot strike (FFS) subjects appear in Fig. 3. For highest accuracy of the Δt₁ measurement, a velocity threshold and projection method is used to eliminate minor fluctuations in the ankle marker velocity profile near zero (see  Appendix). To assess the accuracy of this projection method, the measured ankle marker velocity at the projected time to zero was quantified for all footfalls.

where F_2,avg is the average force of impulse J₂ during the contact time t_c.

The raised cosine bell (RCB) curve function was used to generate the F(t) waveforms of J₁ and J₂ for both foundational and empirical reasons. The RCB function is unique among all bell curve functions in that it can be derived from the first two terms of the Fourier series (see  Appendix). Analyses of vertical ground reaction force–time waveforms from jumping in place indicate that using a function consisting of two bell curves provides a superior representation and requires a lower number of Fourier terms than a half-sine curve (Racic and Pavic, 2009). Direct waveform comparisons indicate that the RCB function provides superior descriptions versus both the half-cosine and Gaussian functions (see cos² data in table 2 of Sim et al., 2008) when vertical ground reaction force impulses are generated at frequencies ≥2 Hz.

where A=2F_avg is the peak amplitude, B is the center time of the peak and C is the half-width time interval. For higher accuracy, the waveform function for impulse J₂ includes an offset to account for the force plate threshold and an asymmetry adjustment to account for the center of mass height difference at takeoff and touchdown. A J₂ force peak corresponding to the minimum height of the center of mass was set at 0.47t_c based on the observations of Cavagna et al. (1977, see their table 4) for human running (see  Appendix for additional details).

Model-predicted versus measured vertical force waveforms were tested for both individual footfalls and trial-averaged force data. For the individual footfall predictions, the input parameters of body mass, ground contact time, subsequent aerial time, m₁ vertical velocity at touchdown and impact interval time Δt₁ were used to generate the model-predicted waveform, which was then compared with the measured waveform from that footfall. For the trial-averaged waveform comparisons, the input parameters from each right footfall were averaged for that trial. A minimum of three and a maximum of six right footfalls were included in the trial average, depending on the number of steps completed during the trial. The measured input parameters from each trial average were used to generate a model-predicted waveform, which was then compared with the measured trial-averaged waveform. An example of this data treatment appears in Fig. 4, including a series of original measured waveforms (Fig. 4A), the trial-average model-predicted waveform (Fig. 4B), and the model-predicted versus trial-average measured waveform (Fig. 4C). Both single-footfall and trial-averaged predictions were assessed to evaluate whether the number of footfalls included influences the predictive accuracy of the model.

To validate the anatomical mass fractions in our model, waveforms were alternatively generated with literature-suggested impact-mass minimum and maximum values (Derrick, 2004). A smaller m₁ quantity of 1.5% was used (representing the approximate mass of the foot; Plagenhoef et al., 1983; Winter, 1990; Hamill and Knutzen, 2009) with a corresponding m₂ quantity of 98.5%. A larger m₁ quantity of 16% was used (representing the approximate mass of the entire stance limb; Plagenhoef et al., 1983; Winter, 1990; Hamill and Knutzen, 2009) with a corresponding m₂ quantity of 84.0%. We predicted that across both speed and foot-strike mechanics, values along the rising edge of the waveform would be under-predicted with an m₁ quantity of 1.5% and over-predicted with an m₁ quantity of 16%.

The predictive accuracy of the model was assessed on the 500 individual footfalls acquired using both the R² statistic and the root mean square error (RMSE) statistic, global values that quantify goodness of fit in relative and absolute terms, respectively. These statistical assessments are broadly used for a variety of purposes, including quantifying the degree of overlap present in time-series data per prior practice (Cohen, 2013; Clark et al., 2014; Clark and Weyand, 2014; Morin et al., 2005). The footfall sample size was sufficiently large to detect very small differences in model performance across foot-strike type and speed using the R² and RMSE statistics.

A two-way ANOVA was used to evaluate whether goodness-of-fit and RMSE values varied as a function of speed (slow, medium, fast) and foot-strike (rear- versus fore-foot categories) classifications (P≤0.05, Table 1) on all 500 footfall waveforms acquired.

A total of 42 subjects, 23 men and 19 women, volunteered and provided written informed consent. The consent process and all experimental procedures were approved by, and conformed to, the approval terms granted by the Institutional Review Board of Southern Methodist University. All subjects were between 18 and 37 years of age and engaged in regular physical activity at the time of testing. The mean age and size characteristics of the men and women were as follows: men: age=23.3±5.0 years, range=18–37 years; height=1.79±0.07 m, range=1.69–1.95 m; mass=81.1±8.5 kg, range=71.0–101.5 kg; and women: age=22.5±1.7 years, range=18–36 years; height=1.68±0.06 m, range=1.55–1.78 m; mass=63.3±9.4 kg, range=43.4–82.0 kg. Subjects included recreationally trained individuals, intercollegiate team-sport athletes, and professional track and field athletes, four of whom were Olympic medalists in sprint or hurdle events.

Data were collected across a range of speeds (3.0 m s⁻¹ to top speed) on a three-axis, custom-built high-speed force treadmill (AMTI, Watertown, MA, USA) capable of speeds of over 20 m s⁻¹. To ensure that the model was being evaluated at speeds that required a normal running gait, only trials above a Froude speed of 1.0 [v/√(gL₀)>1.0] were included in the statistical analysis. For both submaximal and maximal tests, subjects followed testing procedures similar to those described in Weyand et al. (2000, 2010). Thirty-nine of the 42 subjects completed trials up to top speed, and these subjects were habituated to high-speed treadmill running during one or more familiarization sessions before undergoing top speed testing. All subjects wore standardized black compression shirts and shorts, and the same model of running shoes (Nike Waffle Racer version 7.0, Beaverton, OR, USA).

Ground reaction force data were acquired at 1000 Hz and were post-filtered using a low-pass, fourth-order, zero-phase-shift Butterworth filter with a cutoff frequency of 25 Hz (Winter, 1990). For each footfall, the vertical ground reaction force applied during the contact period was determined from the time during which the vertical force signal exceeded a threshold of 40 N. Additionally, trial-averaged vertical ground reaction force waveforms were generated for individual subjects at different trial speeds by averaging the force from each millisecond of the contact period for the right-foot waveforms that had corresponding kinematic data.

For each trial, 3.46 s of video data were collected using a high-speed video system consisting of three Fastec Imaging HiSpec 2G cameras (Mikrotron GmbH, Unterschleissheim, Germany) operating at 1000 frames s⁻¹. Subjects wore reflective markers on the heel and ball of the foot on the lateral aspect of the right running shoe, as well as on the lateral aspect of the joint axes of rotation of the right ankle, knee and hip to capture these respective positions during the trials. A single-frame video file of each subject was recorded prior to testing to serve as a reference for the kinematic analyses. To assess the predictive accuracy of the model across different foot-strike types, footfalls were classified as RFS if the first surface contact occurred on the posterior half of the foot, and FFS if the first surface contact occurred on the anterior half.

The procedures used in extracting three-dimensional marker coordinates from the high-quality multiple camera videos consisted of correcting image distortions, calibrating the three-dimensional space and digitizing the markers. The calibration and digitization routines extensively used functions from the MATLAB Image Processing Toolbox (MathWorks, Natick, MA, USA). The calibration and digitization MATLAB routines were developed by the Hedrick Lab at the University of North Carolina (Hedrick, 2008). A resolution of 0.7 mm was measured under dynamic conditions using the high-speed video system with the custom MATLAB image correction, calibration and digitization routines. Data acquisition timing for the AMTI DigiAmp amplifier and Fastec Imaging cameras was synchronized through hardware interfaces. The digitized marker data were filtered at 25 Hz using the same filter as for the force data.

The goodness-of-fit agreement (R²) between the 500 ground reaction force waveforms we measured and those predicted by our two-mass model approached unity as hypothesized (Table 1). The corresponding error of prediction expressed in force units standardized to the body's weight was slightly greater than 0.2W_b and was equal to 11.5% of the mean stance-averaged vertical force from all 500 trials (1.82±0.23W_b). The R² and RMSE statistics for the 108 trial-averaged waveform values were nearly identical to those for the 500 individual waveforms. Specifically, the goodness of fit between model-predicted and measured trial-averaged ensemble waveforms was only 0.01 greater (R²=0.96±0.03) than the individual footfall value, while the RMSE of prediction was only 0.02W_b less (0.19±0.07W_b) than the corresponding individual footfall values.

The goodness-of-fit between the model-predicted and measured waveforms (R²) was nearly identical (ΔR²=0.01) for the 177 RFS versus 323 FFS waveforms (Table 1). Because of the large number of footfalls tested and minimal variability in model predictive accuracy, the 0.01 greater R² value for the FFS versus RFS waveforms was significantly different statistically. Corresponding RMSE differences for the RFS versus FFS waveform predictions were small (ΔRMSE=+0.04W_b) and not statistically different.

The goodness-of-fit between the model-generated and actual waveforms was slightly, but significantly, better for waveforms acquired from slow and medium speeds versus those from fast speeds. The R² values, when averaged for the waveforms from the three speed ranges (R² total, Table 1) regardless of foot-strike type, were slightly, but significantly lower (ΔR²=−0.03) for the faster versus medium and slower speed trials. Similarly, the RMSE of the model fits to the slow and medium speed waveforms were not different from one another; however, both were significantly smaller than the error of prediction for the waveforms from the fastest speeds. Similar across-speed patterns were present for both the R² and RMSE statistics when the waveforms were analyzed within either the RFS or FFS mechanics categories (Table 1).

The waveforms generated with the two-mass model accurately predicted the more rapid rising edges of the RFS versus FFS waveforms, including transient rising-edge peaks when present, regardless of the running speed, and total waveform amplitude and duration (Figs 5A,C,E and 6A,C,E). Waveform shape variability across foot-strike types was accurately predicted from the shorter deceleration periods of the m₁ mass segment for RFS versus FFS (Δt₁, Table 2) with little difference in the overall J₁ impulse values (Fig. 5B,D,F versus Fig. 6B,D,F; Δv₁, Table 2) at similar speeds. From slower to faster speeds, J₁ impulses predicted by the model became greater for both RFS and FFS waveforms as a result of the greater pre-impact ankle velocities (Table 2). The close fits to the middle and later portions of all the waveforms resulted largely or exclusively from the J₂ impulse predicted from the model because the m₁ impact deceleration event concluded during the first quarter to half of the total contact period.

The measured ankle marker velocity at the time our projection technique identified a zero value was, on average, −0.05±0.04 m s⁻¹ for the 500 individual trials. All but two of the footfalls had a measured velocity within ±0.20 m s⁻¹ of a literal zero value (i.e. 0.00 m s⁻¹).

Poorer agreement between model-predicted and measured waveforms resulted from using m₁ segment values that were either smaller or larger than the originally assigned anatomical model value of 8.0% of m_b. As hypothesized, waveforms generated using lesser m₁ values (m₁=1.5% with m₂=98.5% of m_b) consistently under-predicted the force values measured along the rising edge of the waveforms (Figs 7A,D,G and 8A,D,G). Conversely, waveforms generated using greater m₁ values (m₁=16.0% with m₂=84.0% of m_b) consistently over-predicted measured rising-edge force values (Figs 7C,F,I and 8C,F,I). Differences were more pronounced at the faster trial speeds because of the greater m₁ segment decelerations and correspondingly larger J₁ impulses at faster speeds.

The R² goodness of fit values for the 500 individual footfall waveforms generated using m₁ values equal to 1.5% and 16.0% of m_b accounted for 12% and 21% less variance, respectively, versus the original m₁=8.0% results (m₁−1.5% R²=0.83±0.16 and m₁−16.0% R²=0.74±0.21). Predictive error values for the waveforms generated using the two alternative m₁ segment values were approximately twice as large as those obtained using the original 8.0% value (m₁−1.5% RMSE=0.36±0.17W_b; m₁−16.0% RMSE=0.47±0.24W_b).

As hypothesized, we found that the force contributions of two discrete body mass components do indeed suffice to predict running vertical ground reaction force–time waveforms in full. Our two-mass, two-impulse model independently predicted nearly all of the variation in measured running ground reaction force waveforms, which differed considerably in their shape, amplitude and duration. Although our prior, best-fit analysis indicated that this outcome might be theoretically possible (Clark et al., 2014), an experimental test incorporating simultaneous motion and force data had not been previously undertaken. Our direct test here indicated that regardless of whether our 42 human subjects jogged, ran at intermediate speeds or sprinted, and whether they first contacted the running surface with the fore or rear portions of their feet, there was near-complete agreement between the model-predicted and measured force–time waveforms across the 500 footfalls we analyzed (Table 1). Thus, we conclude that the vertical ground reaction forces of human runners can be broadly understood from the motion of two discrete portions of the body: (1) the contacting lower limb and (2) the remainder of the body's mass.

A primary scientific challenge here was not knowing a priori how well the waveform contributions from the 92% of the body's mass could be predicted when modeled as a single mass component. Conceivably, the summed force contributions resulting from the motion of the head, arms, trunk, upper portion of the contacting leg and full mass of the non-contacting leg might defy accurate prediction when treated as a single mass (Bobbert et al., 1991). This large, multi-jointed mass component lacks a fixed, readily measurable center because of the non-uniform motion of the different segments that comprise it. Hence, the complexity of the within-segment and total motion of our model's mass component m₂ and its corresponding force contribution would be difficult to measure and predict from positional data. We ultimately implemented an indirect approach that allowed the force contributions of mass m₂ to be quantified without position and time data. We simply subtracted impulse J₁ from the total ground reaction impulse J_T, after quantifying the latter from body mass, gravity and the contact proportion of the total step time (Eqns 1–3). The resulting fits supported the efficacy of the approach as the agreement between model-predicted and measured waveforms was consistently excellent over the later portions of the stance period where the total waveform is predominantly determined by impulse J₂ (Figs 5A,C,E and 6A,C,E, waveform trailing edges).

A potential limitation of our impulse subtraction approach was the required assumption that the net vertical displacement of the body's center of mass over the course of one or many strides is zero. However, our analysis indicated that little to no predictive error was introduced by this assumption under our level, treadmill-running test conditions. We found similar levels of predictive accuracy for trial-averaged and individual-footfall waveforms even though non-zero vertical displacements of the center of mass, when considered on a per-step basis, could have been substantially greater over the course of an individual step versus multiple consecutive steps. However, for our across-speed comparisons, it seems likely that our slightly less accurate waveform predictions for the fastest speeds, where step-to-step mechanics are generally less consistent versus medium and slower speeds (Weyand et al., 2010), are probably attributable to marginally greater violations of this assumption for the faster trials (Table 1).

An additional technical challenge involved accurate and consistent quantification of the duration of mass component m₁'s impact-period deceleration to a zero velocity (i.e. Δt₁ in Eqn 3; see also Fig. 2) for all footfalls. Following RFS impacts, the ankle marker position–time trajectories consistently provided a discrete vertical position minimum corresponding to the zero velocity needed to quantify half-impact duration Δt₁ in our model (Fig. 3). However, for some FFS footfalls, the rates of the positional change versus time during the last fraction of the deceleration period were less consistent and more prolonged than the FFS data in Fig. 3. This led us to implement an ankle marker projection technique (see  Appendix) to avoid basing impulse J₁ predictions on Δt₁ values that, in these cases, are not representative of the overall timing of the impact-deceleration event. Upon implementation across all 500 footfalls, we found that the measured ankle marker velocities at the time that the projection technique identified a zero velocity value were very near the true zero value desired (mean v_ankle=−0.05±0.04 m s⁻¹). Given the overall mean ankle marker vertical velocity, Δv₁, of −1.59±0.02 m s⁻¹ at impact across all 500 footfalls (Δv₁, Table 2), this method, on average, captured 97% of m₁ total post-impact deceleration to zero.

The ability to consistently predict rising-edge waveform peaks that occurred from 15 to 50 ms after initial impact is a noteworthy aspect of our model validation. As the timing of rising-edge peaks on individual waveforms is determined by the overlap of the two impulses in our model, successful prediction required precisely capturing the timing of both the high- (J₁) and low-frequency (J₂) components of the waveforms. The timing of impulse J₁ was determined from motion data, and was therefore fully independent of our force data, filtering and processing routines. Conversely, the timing of impulse J₂ was directly dependent on our force data and processing routines, and was therefore fully independent of our kinematic data and processing routines. Had even a minor degree of temporal inaccuracy been present in either our kinematic or kinetic data, the predictive accuracy with which the two-mass model identified rising-edge force peaks would not have been possible.

Our experimental goal was to identify the most concise mechanical explanation that running ground reaction force waveforms might have. The multi-mass models, in contrast, were formulated to evaluate the potential influence of numerous factors on the waveform rising edge. Many of the features incorporated into the multi-mass models provide reasonable theoretical representations of the numerous, potentially influential musculoskeletal complexities present (Liu and Nigg, 2000; Ly et al., 2010; Nigg and Liu, 1999; Nikooyan and Zadpoor, 2011; Zadpoor and Nikooyan, 2010). These include mass components that vary in stiffness, that are both rigid and wobbling in nature, and that are connected with both serial and parallel elements (Fig. 1A). While the quantitative descriptions derived for these features undoubtedly include uncertainty, their basic influence on the RFS jogging waveforms thus far analyzed are plausible and largely consistent with other experimental and modeling approaches (Gruber et al., 1998; LaFortune et al., 1996; Pain and Challis, 2001; Shorten and Mientjes, 2011). Accordingly, our demonstration that a substantially more concise model can account for running ground reaction force waveforms in full should not be regarded as incompatible with the more theoretical results the multi-mass models have provided. Indeed, the most reasonable conclusion from the different approaches is that the collective influence of the many biological complexities incorporated into the multi-mass models is, in sum, accurately described by the concise impulse–momentum relationships the two-mass model provides.

In contrast to the prevailing view that the impact forces resulting from RFSs and FFSs involve different mass quantities (Lieberman et al., 2010; Nigg, 2010), our results indicate that mass quantities and force–motion relationships do not differ across strike types. This becomes fully apparent when J₁ impact impulses are quantified using measured, foot-strike-specific deceleration durations in conjunction with the invariant, anatomically based mass quantities in our model. With both factors in place, we were able to accurately predict the waveforms in full for both foot-strike types at slow, intermediate and fast running speeds alike.

As is evident from both our impulse J₁ illustrations (Figs 5 and 6) and the Δv₁ values in Table 2, rising-edge force peaks that are generally visible for RFS but not FFS waveforms can be fully attributed to the different m₁ deceleration durations we measured. Longer FFS deceleration durations reduce and delay the J₁ peak force values; they also result in greater waveform force contributions from impulse J₂ at the time of the J₁ force peak. In combination, these timing-dependent factors typically do not allow the J₁ impulse peak to introduce a localized force peak along the rising edge of the full waveform (Fig. 6B and Fig. 8B,E). We found this to be the case even though the total ground reaction impulses resulting from FFSs in our data set were just as large as those resulting from RFSs (note: J₁ impulses are ∝Δv₁ in Table 2). One noteworthy exception to these general foot-strike-specific waveform patterns has recently been documented for specialized sprinters (Clark and Weyand, 2014). These fore-foot striking athletes have waveforms that are characterized by prominent rising-edge force peaks that result from brief, large J₁ impulse peaks. The high pre-impact limb velocities (Δv₁) and brief impact deceleration periods (Δt₁) of specialized sprinters introduce conspicuous rising-edge force peaks, particularly at faster running speeds (Fig. 6E,F and Fig. 8H).

Our waveform predictions using alternative mass fractions for segments m₁ and m₂ also support the validity of describing FFS and RFS waveforms with the same fractional body mass quantities and force–motion relationships. As hypothesized for these alternative mass distributions, substantially reducing the lower limb mass value of m₁ from the anatomical fractional value of 0.08 (Plagenhoef et al., 1983) to the much lower value of 0.015 resulted in predicted rising-edge force values that fell consistently below measured values (Figs 7A,D,G and 8A,D,G). Conversely, substantially increasing the m₁ fractional mass value to 0.16 resulted in predicted rising-edge force values that fell consistently above the measured values (Figs 7C,F,I and 8C,F,I). The magnitudes of the respective predictive errors and trends observed across speed were similar for the two foot-strike types using the aforementioned alternative m₁ fractional mass values. In both cases, the increases in the magnitude of the J₁ impulse across speed introduced greater predictive errors at higher versus lower speeds. More globally, the contrast between the systematic predictive errors present in both sets of alternative-mass generated waveforms versus the near-full agreement achieved with the original values provides strong support for the validity of the anatomical values originally assigned to mass components m₁ and m₂ in our two-mass model.

The accuracy and conciseness of our mechanical explanation for the variable ground reaction force–time patterns of dozens of human runners performing across their full range of speeds raises a noteworthy possibility. Our two-mass, two-impulse approach may offer a mechanical explanation that generalizes to the ground reaction force–time patterns of other running species. However, for non-human runners, the lesser distal-limb mass segments typically present (Hildebrand, 1960; Rubenson et al., 2011) could alter both the form and applicability of the two-mass approach. Minimally, species-specific mass distributions would need to be incorporated to generate waveforms from the basic stride measures included in the model. In this regard, comparative anatomical variation constitutes both a challenge and an experimental opportunity to evaluate the basic tenets of the model. More broadly, the testable hypotheses the model offers should be tractable using a variety of direct and indirect approaches.

Finally, our empirical validation of a concise model that can fully predict running vertical ground reaction forces offers basic insights with immediate potential for application. In contrast to simplified, single-mass models (Blickhan, 1989; Blum et al., 2009; Farley and Gonzalez, 1996; McMahon and Cheng, 1990; Silder et al., 2015), which are incapable of capturing the high-frequency, impact-phase components of the waveform (Bullimore and Burn, 2007; Clark and Weyand, 2014; Shorten and Mientjes, 2011), and multi-mass models that do not account for the whole waveform and are too theoretical and complex for practical application, the two-mass model requires only body mass and very limited motion data in order to predict the waveform in full. These attributes enable indirect assessment of impact forces, limb loading rates, and other informative, force-based outcomes using video or other inexpensive motion-sensing technologies. Potential model applications include: informing the design of running robots, exoskeletons and prostheses, the customization of running shoes and orthotics from individual gait mechanics, the development of wearable ground force sensors, and the improvement of evidenced-based feedback for gait analyses, intervention and modification.

where m₁ is 8.0% of the body's mass, Δt₁ is the time interval between touchdown and the vertical velocity of m₁ slowing to zero, Δv₁ is the change in vertical velocity of m₁ during Δt₁, and F_1,avg is the average force during the total time interval (2Δt₁) of impulse J₁.

where A is the peak amplitude, B is the center time of the peak and C is the half-width time interval. Because of the simple properties of this function, peak amplitude A=2F_avg, and the area under the curve is J=AC.

where A₁=2F_1,avg using F_1,avg in Eqn A8, and B₁ and C₁ equal the time Δt₁ after touchdown for the vertical velocity of m₁ to reach zero.

where A₂=2F_2,avg using F_2,avg in Eqn A9, and B₂ and C₂ equal 0.5t_c for a symmetrical waveform.

Impulse J₁ results from the acceleration of the lower limb during the limb impact interval and is quantified by Eqn A8. Δv₁ and Δt₁ are determined using ankle marker kinematics to represent the motion of the lower limb mass m₁. The lower limb attains a relatively constant positive velocity after the impact interval (Fig. 2C), resulting in a near-zero acceleration of m₁ (Fig. 2D) and negligible force F₁=m₁g (Fig. 2E) as described in Eqn A6. This force is less than 55 N for a subject with body mass m_b=70 kg and lower limb mass m₁=5.6 kg (m₁=0.08×70 kg). Thus, J₁ can be modeled by a finite impulse during the impact interval.

Kinematic data for the ankle marker position during the impact interval for representative RFS and FFS subjects appear in Fig. 3. Δt₁ is the time interval between touchdown and the vertical velocity of m₁ slowing to zero. For some FFS runners at slower speeds, minor fluctuations in the ankle marker velocity–time profile near the end of the m₁ impact interval can create variability in the Δt₁ measurements as a result of slow m₁ impact velocities and skin marker motion artifact during the impact interval (Bobbert et al., 1991). Accordingly, Δt₁ was quantified using a technique that represented the functional end of the m₁ impact time interval. An ankle marker velocity of −0.25 m s⁻¹ was used as a threshold point, and the previous 10 ms of data were utilized for a linear projection of the ankle marker velocity to zero. For consistency in analysis, this method was applied to all ankle marker kinematics data, regardless of the speed or foot-strike mechanics of the runner (see Results).

where A₂ is the peak amplitude, B₂ is the center time of the peak, C_2L is the leading half-width time interval, and C_2T is the trailing half-width time interval. The location of peak B₂ was set at 0.47t_c as per the center of mass asymmetry originally reported by Cavagna et al. (1977) (see their table 4). With the symmetry control, C_2L=B₂ and C_2T=t_c−B₂.

This equation is specifically evaluated for each width parameter at time t where the force F₂ is equal to the threshold. A force threshold of 40 N was used for the AMTI high-speed force treadmill system. The leading width parameter C₂=C_2L is determined at the first channel (t=1 ms) and the trailing width parameter C₂=C_2T is determined at the last channel (t=t_c). Leading and trailing width parameters C_2L and C_2T are calculated using the same peak location B₂ and peak amplitude A₂.

As a result of this recalculation, the impulse is approximately equal to the original impulse and the values at the first and last channel are approximately equal to the threshold value.

The authors thank Andrew Udofa for his help with data collection and the artwork in Fig. 1B, and the participants for their time and effort.