The force response of muscles to activation and length perturbations depends on length history

Muscles produce variable forces to adjust to unexpected perturbations during locomotion on uneven terrain (Daley and Biewener, 2011). Compared with steady-state isometric contractions (Baratta et al., 1998; Kim et al., 2020; Vredenbregt and Rau, 1973) or steady-state locomotion on uniform terrain (Eng et al., 2019; Roberts et al., 1997), our ability to predict in vivo muscle force under dynamic conditions remains relatively poor (Dick et al., 2017; Lee et al., 2013). As most daily activities are performed under unsteady or perturbed conditions, understanding how muscles produce force under dynamic conditions is essential to understand the control of movement.

In the traditional view of ‘muscle as motor’ (Nishikawa, 2020; Seth et al., 2011, 2018), activation dynamics, typically derived from activation and deactivation kinetics during isometric twitch (Thelen, 2003), represents the contractile element force over time, and this force is scaled based on isometric force–length and isotonic force–velocity relationships. The phase of activation determines the onset and duration of muscle force, and the magnitude of activation is proportional to the force amplitude (Lloyd and Besier, 2003; Zajac, 1989). That is, force is mainly controlled by activation (Brenner and Eisenberg, 1986; Davies et al., 1967). However, Hill-type models based on the ‘muscle as motor’ viewpoint have failed to accurately predict muscle force in dynamic conditions (Dick et al., 2017; Lee et al., 2013). Some limitations of Hill-type models include the fact that dynamic in vivo conditions are neither isometric nor isotonic, and the models fail to include history-dependent effects (McGowan et al., 2013).

Activation dynamics also depend on muscle length history (Sandercock and Heckman, 1997; Shue and Crago, 1998). Muscles require more time to develop force when they shorten during activation and less time when they lengthen. Accordingly, Roberts and Gabaldón, (2008) noted that the relationship between muscle activation and force is not constant. The electromechanical delay varies depending on the task, and the relationship between mean electromyogram (EMG) amplitude and force amplitude depends on the phase of the gait cycle. The force response to EMG amplitude in the swing phase is lower than that in the stance phase. The ‘muscle as motor’ paradigm clearly applies to near-isometric and isotonic conditions, but its usefulness for predicting muscle force during natural movements is less clear because natural movements are characterized by time-varying applied loads that produce length perturbations in active muscles (Daley and Biewener, 2011).

An alternative to the view of ‘muscles as motors’ is the view that muscles are tunable materials (Higueras-Ruiz et al., 2021; Nishikawa and Huck, 2021). In this view, activation modulates muscle viscoelastic properties, and force develops in response to deformation by applied loads. Thus, muscle force is determined not only by activation but also by deformations in length. This view allows muscle force to be decoupled from activation, given the same length and velocity.

Several recent studies suggest that muscle force may be decoupled from activation via effects of length changes (Daley, 2018; Daley and Biewener, 2011; Daley et al., 2009; Lappin et al., 2006; Libby et al., 2020; Robertson and Sawicki, 2015; Sponberg et al., 2011b). In cockroaches, the mechanical response of a leg muscle differed depending on the motor task, such as standing or running, despite similar activation (Sponberg et al., 2011a). In addition, the shape of muscle work-loops is sensitive to small changes in velocity (Robertson and Sawicki, 2015). Furthermore, it is well known that isometric force following active lengthening is higher (Abbott and Aubert, 1952; Edman et al., 1982; Herzog et al., 2016; Nishikawa, 2016) and isometric force following shortening is lower (Abbott and Aubert, 1952; De Ruiter et al., 1998; Herzog and Leonard, 1997; Lee and Herzog, 2003; Tahir et al., 2020) than the purely isometric force at the same final length and activation. These observations demonstrate that muscle force in dynamic contractions is not well explained by activation alone.

Previous studies of in vivo muscle force–length dynamics in guinea fowl during running over potholes and obstacles have shown that length perturbations are associated with changes in muscle force (Daley and Biewener, 2011; Daley et al., 2009). Length perturbations are defined as abrupt changes in muscle length and velocity produced by applied external forces, which lead to history-dependent effects including force enhancement and depression (Libby et al., 2020). For example, the lateral gastrocnemius muscle undergoes length perturbations during foot contact following active shortening. Uneven terrain, such as obstacles and potholes, alters the phase of length perturbations compared with level running, which results in variable timing and magnitude of peak force (Fig. 1). When stepping up onto an obstacle, length perturbations occurred earlier and peak muscle force increased (Fig. 1, line ‘a’), whereas stepping down from obstacle to the treadmill delayed the phase of the length perturbation and resulted in lower peak force. When stepping down, average EMG intensity was similar to the ‘no obstacle’ condition, although EMG intensity increased on average during stepping up (Daley and Biewener, 2011). These results suggest a decoupling between muscle force and activation under dynamic conditions, and that the phase of length perturbations appears to play an important role in determining muscle force.

The muscle force produced in response to length perturbations can be increased or decreased by length history, given constant activation (Libby et al., 2020). Libby et al. (2020) demonstrated that length history influences preload, or the force at the onset of a length perturbation. Length history was varied by controlling the frequency of sinusoidal length changes. The phase determined whether length history was composed of shortening or lengthening and the frequency determined the velocity of the length change. When the preload was smaller, the force response to a length perturbation was smaller (Libby et al., 2020). In the context of the ‘muscle as tunable material’ paradigm, the length history-dependent preload is equivalent to muscle stiffness. These results suggest that the relative contribution of length perturbations to muscle force is modulated by length history.

The contribution of activation to muscle force may also be affected by length history. As noted previously, the relationship between force and activation depends on the task (Roberts and Gabaldón, 2008). If tasks differ with respect to muscle length history, then the response of muscle force to activation could be influenced by length history. Yet, how length history affects the activation–force relationship has rarely been studied.

The ‘muscle as tunable material’ paradigm predicts that the force response to deformation will be related to muscle stiffness. Because muscle force is directly proportional to stiffness (Kirsch et al., 1994) and shortening decreases muscle force as velocity increases (Julian and Sollins, 1975), we varied shortening velocity to modulate muscle stiffness. We used two protocols to test the hypotheses that length history affects the force response to changes in activation amplitude (protocol 1) or length perturbations (protocol 2) via changes in muscle stiffness. In both protocols, we designed experiments in which muscles first develop force isometrically and then are shortened at different velocities to modulate muscle stiffness.

In protocol 1, we used time-varying activation during active shortening to quantify how differences in muscle stiffness affect the relationship between activation amplitude and muscle force. In protocol 2, length perturbations with varying amplitudes were applied to muscles during active shortening to investigate the effect of muscle stiffness on the force response to length perturbations. These protocols were designed to provide insight into the inter-relationships among muscle force, length perturbations and activation.

Mice of the strain B6C3Fe a/a-Ttn^mdm/J were obtained from the Jackson Laboratory (Bar Harbor, ME, USA) and a breeding colony was established in the animal care facility at Northern Arizona University (NAU). Soleus muscles from homozygous wild-type mice (n=6 female, n=6 male; mean±s.d. mass 34.3±6.7 g; age 128– 199 days) were used in this study. The institutional Animal Care and Use Committee at NAU approved the experimental protocol and use of animals.

Mice were euthanized by isoflurane overdose (1–2% isoflurane mixed with oxygen) confirmed by cervical dislocation. Soleus muscles were isolated from the hindlimb and immersed in a mammalian Krebs–Ringer bath (in mmol l⁻¹: 137 NaCl, 5 KCl, 1 MgSO₄, 24 NaHCO₃, 1 NaH₂PO₄, 11 dextrose and 2 CaCl₂, pH 7.2 at room temperature) with 95% O₂ and 5% CO₂. Soleus muscles were attached to a dual servomotor force lever (Aurora Scientific, Inc., Series 300B, Aurora, ON, Canada) at one end and an inflexible hook at the other end. The servomotor controlled muscle length and measured position and force.

All soleus muscles were stimulated using an electrical field generated between two platinum electrodes connected to a Grass S88 stimulator. In order to determine optimal muscle length (L₀), muscles were fully activated using square-wave pulses at 60 V and 75 Hz, and length was adjusted until maximum isometric force (P₀) was established. At the end of the experiments, P₀ was measured to ensure that it was maintained throughout the experiments. Muscles used in the analysis experienced <15% drop in force during an experiment.

Two protocols were developed to evaluate the effects of active shortening on muscle stiffness, as well as the effects of muscle stiffness on the response to activation (protocol 1, n=7) or length perturbations (protocol 2, n=5). Protocol 1 measured how muscles respond to time-varying activation depending on muscle stiffness. Muscles were stimulated isometrically at submaximal voltage (mean±s.d. 3.87±0.88 V) and frequency (5 Hz) until they reached a mean (±s.d.) steady-state force of 7.0±3.0% P₀ in ∼0.7 s. This submaximal stimulation was chosen to represent a low-effort contraction in vivo. Once isometric force reached steady state, active shortening was performed at two velocities (1.0 L₀ s⁻¹, 0.5 L₀ s⁻¹) to cover the maximal power range of soleus muscle (Askew and Marsh, 1997). At the midpoint of active shortening, muscle length was held constant at L₀ (Fig. 2).

To evaluate how muscle force depends on the amplitude of activation, the amplitude of activation during active shortening was varied using an A-M Systems Model 4100 stimulator (A-M Systems, Sequim, WA, USA). Activation started increasing at the onset of active shortening. The amplitude of activation increased linearly for 100 ms, and then decreased linearly for 100 ms (Fig. 2). Peak activation occurred at the midpoint of active shortening at optimal muscle length (L₀). Six different activation amplitudes were used: 100%, 110%, 120%, 130%, 140% and 150% of the amplitude for the first submaximal isometric contraction (see above). This range of activation produced forces ranging from ∼3% to 30% P₀, similar to low-effort in vivo movements.

In protocol 2, the response of muscle force to length perturbations was measured for different values of muscle stiffness, established by varying the shortening velocity. As in protocol 1, isometric contraction and active shortening at 1.0 and 0.5 L₀ s⁻¹ were performed (Fig. 3). Length perturbations were performed at the midpoint of active shortening where muscle length was L₀ in all trials. The shape of the length perturbation was varied based on a half sinusoidal curve at a frequency of 2.5 Hz. For the trials with no active shortening, the amplitudes of the sinusoidal curves were 0.2%, 0.4%, 0.6%, 0.8%, 1% and 3% of L₀ (Fig. 3B), and for the trials with active shortening, the amplitudes were 1%, 3%, 5%, 7%, 9% and 11% of L₀ (Fig. 3A). The amplitude of length perturbations was chosen to produce approximately the same force for all perturbations, as well as to fall within the operating range of most muscles (±20% L₀; Burkholder and Lieber, 2001). Activation was constant at ∼35% P₀ (Fig. 3A,B).

To measure instantaneous muscle stiffness in protocols 1 and 2, quick length transients were performed at L₀ during active shortening (Fig. 4). In a preliminary study, we found that the smallest transient that produced a reproduceable change in force was 0.9% L₀. The speed of the transient (3 L₀ s⁻¹) was based on a previous study (Pinnell et al., 2019). Instantaneous muscle stiffness was obtained from the linear portion of the force–time trajectory during a quick transient. Muscle stiffness was calculated by the force difference between the linear portion divided by the amplitude of the transient in the linear portion.

Data were sampled at 4000 Hz and collected using a DAQ box (National Instruments, Austin, TX, USA). Muscle stiffness, peak muscle force and peak muscle length during active shortening following submaximal isometric contraction were normalized to values measured at P₀. The relative values were used for statistical analysis. The sample sizes were sufficient to yield power estimates of 0.8 given the observed variance and effect sizes. Differences were considered significant if P<0.05 and Tukey's honestly significant difference (HSD) test was used to evaluate post hoc differences among means. All data are reported as means±s.d.

In protocol 1, two-way ANOVA was used to test whether active shortening and activation amplitude induce changes in muscle stiffness. To evaluate the relationship between variable muscle stiffness and force resulting from active shortening at different activation amplitudes, the correlation coefficient between muscle stiffness and force was calculated using standard linear regression. In order remove the effect of shortening velocity, muscle force with time-varying activation was subtracted from muscle force with constant activation (Fig. 5A). The peak force after subtraction (Fig. 5B) was used to analyze the effect of activation on muscle force during active shortening. Analysis of covariance (ANCOVA) was performed to identify changes in mean peak force, with shortening velocity as the covariable and amplitude of activation as the main effect. In protocol 2, one-way ANOVA was performed to investigate the effect of shortening velocity on muscle stiffness. ANCOVA was used to test the effects of length perturbations and shortening velocity on muscle force, with shortening velocity as the covariable and length perturbations as the main effect. The data are available from the Dryad digital repository (https://doi.org/10.5061/dryad.5hqbzkh88).

In protocol 1, relative muscle stiffness increased significantly with increasing activation amplitude (ANOVA, F=13.28, P<0.001) and decreasing shortening velocity (F=102.57, P<0.001) as expected. The interaction between activation amplitude and shortening velocity on muscle stiffness was not significant (F=1.18, P=0.32). Relative muscle stiffness at the highest activation amplitude (150% of activation amplitude during the first submaximal isometric contraction) was significantly higher than that at 130% (HSD, P<0.05), 120% (HSD, P<0.01) and 110% (HDS, P<0.001) (Fig. 6A). Stiffness at 140% and 130% was significantly higher than that at 110% (HSD, P<0.001) (Fig. 6A). All shortening group means were significantly different from each other (HSD, P<0.001) (Fig. 6B). Linear regression analysis showed a significant correlation between muscle stiffness and force (r²=0.82, P<0.001; Fig. 7). The slope of the relationship between maximum force and activation amplitude differed significantly among shortening velocities (ANCOVA, F=16.96, P<0.001). The slope for isometric contractions was greater than the slope at 0.5 L₀ s⁻¹ (HSD, P<0.001) and 1.0 L₀ s⁻¹ (HSD, P<0.001) (Fig. 8A). The intercept did not differ significantly among shortening velocities.

In protocol 2, ANOVA demonstrated a significant effect of shortening velocity on muscle stiffness (F=199.17, P<0.001). All shortening velocity group means were significantly different from each other (HSD, P<0.001; Fig. 6C). As in protocol 1, muscle stiffness decreased significantly with increasing shortening velocity. The intercept (ANCOVA, F=33.29, P<0.001) and slope (ANCOVA, F=4.13, P<0.05) of the relationship between length perturbation amplitude and maximum force differed significantly among shortening velocities. For the isometric condition, the maximum force increased faster with length perturbation amplitude than for 0.5 L₀ s⁻¹ (HSD, P<0.05) and 1.0 L₀ s⁻¹ (HSD, P<0.05) (Fig. 8B). Additionally, the isometric condition had a greater y-intercept than for 0.5 L₀ s⁻¹ (HSD, P<0.05) and 1.0 L₀ s⁻¹ (HSD, P<0.001) (Fig. 8B). There was no significant difference in slope or y-intercept between 0.5 L₀ s⁻¹ and 1.0 L₀ s⁻¹.

Under isometric conditions, activation is primarily responsible for the production of muscle force (Inman et al., 1952; Seyfert and Künkel, 1974; Smith and Stokes, 1993), whereas under dynamic conditions, the contribution of activation to force depends on length history (Daley and Biewener, 2011; Daley et al., 2009; Sponberg et al., 2011a,b). Additionally, muscles respond differently to length perturbations depending on length history (Libby et al., 2020). In the force–velocity relationship, as active shortening velocity increases, force decreases. In the present study, we demonstrated that active shortening reduces muscle stiffness, and further showed that the force response of muscles to activation and length perturbations depends on muscle stiffness.

In order to evaluate whether the contribution of activation to muscle force depends on muscle stiffness, we varied muscle stiffness by shortening muscles actively at different velocities following submaximal isometric contraction. Under isometric conditions, activation was strongly correlated with muscle force and stiffness. However, the contribution of activation to muscle force was less clear under dynamic conditions. If variation in activation could explain variable forces during locomotion, then muscle stiffness should be modulated by varied activation history, rather than length history. However, during locomotion, the activation pattern is periodical with relatively little variation. Previous studies showed that the variation in activation does not explain the variation in force (Daley and Biewener, 2011; Daley et al., 2009). Therefore, to understand how muscle produces variable forces against unexpected perturbations during locomotion in this study, muscle stiffness was modulated by varying length history rather than through variation in activation history.

As expected, protocol 1 showed that muscle stiffness decreased significantly with increasing shortening velocity (Fig. 6B) and decreasing activation amplitude (Fig. 6A). Reduced muscle stiffness was associated with a decrease in muscle force (r²=0.90, P<0.001; Fig. 7) in agreement with previous studies (Kirsch et al., 1994; Sugi and Tsuchiya, 1988). The results suggest that muscle stiffness is regulated by both activation and shortening velocity. Given different length and activation history, muscles can have different stiffness even at the same length.

Based on the relationship between activation amplitude and muscle stiffness, we evaluated how active shortening at different velocities affects the force response of muscles to activation. We observed that muscles were more sensitive to activation amplitude in isometric conditions than during active shortening (Fig. 8A). As muscle stiffness decreases as shortening velocity increases, the results support the hypothesis that the relative contribution of activation to muscle force decreases with decreasing muscle stiffness.

In protocol 2, we investigated how muscle stiffness modulates the effect of length perturbations on muscle force. Length perturbations were performed following active shortening at different velocities to modulate muscle stiffness. The effect of length perturbations on muscle force also varied with shortening velocity (Fig. 8B). Muscle force increased faster with increasing amplitude of length perturbations during isometric contractions than during active shortening (Fig. 8B). When muscles shortened actively, a greater length perturbation was needed to produce the same force than for isometric contraction. Similar to Libby et al. (2020), these results demonstrate that the effect of length perturbations on muscle force depends on muscle stiffness at the onset of the perturbation. Greater muscle stiffness produces greater peak force during length perturbations for the same length and activation.

Although it has been demonstrated that muscles can perform a variety of roles depending on the phase of activation (Ahn et al., 2006; Daley and Biewener, 2003; Dickinson et al., 2000; Roberts et al., 1997), how muscle force adjusts to unsteady conditions that pertain to natural movements remains less well understood. Our findings in this study provide insight into muscle force during unsteady movement by investigating how the response to length perturbations depends on stiffness and activation.

Recent in vivo data from the lateral gastrocnemius (LG) of guinea fowl running on a treadmill with obstacle perturbations (Daley and Biewener, 2011) offer a unique opportunity to evaluate how muscle force is controlled during unsteady movement. Daley and Biewener (2011) measured in vivo LG force, strain and activation during treadmill running over obstacles using tendon buckles, sonomicrometry and EMG. Peak LG force during stance was varied to negotiate obstacles. When stepping up onto obstacles, the average peak LG force was greater, but when stepping down, the force was less than during level locomotion. Interestingly, the peak force was more highly correlated with LG length at foot contact (r²=0.63; Fig. 9A) than with peak EMG intensity (r²=0.34; Fig. 9B).

The effect of earlier foot contact for stepping up is to reduce the LG shortening velocity earlier than in level steps. When stepping up, foot contact occurs earlier (Fig. 1, line ‘a’) than in level steps (Fig. 1, line ‘b’). The earlier reduction in shortening velocity allows for a larger preload at the onset of length perturbations in level steps, which results in greater peak force. These guinea fowl LG data are in accordance with Libby et al.’s (2020) study which also showed that work during length perturbations depends on preload (Libby et al., 2020). The results are also consistent with the results of the present study, which demonstrates that the force response to length perturbation depends on length history-dependent muscle stiffness (Fig. 8B) because muscle stiffness is highly correlated with preload, i.e. muscle force at the onset of the perturbation (Libby et al., 2020). In addition, the low r² value between peak LG muscle force and peak EMG intensity (Fig. 9B) is also consistent with our results demonstrating that activation has a relatively weak effect on muscle force during active shortening (Fig. 8A, 1.0 L₀ s⁻¹). Taken together, our results and the in vivo guinea fowl data (Daley and Biewener, 2011) suggest that peak force during active shortening is determined largely by the length history-dependent force response to length perturbations, rather than activation amplitude.

In the traditional view, muscles function as a motor (Nishikawa, 2020). Force is mainly produced by activation, but can be modified by length and velocity (Hof and Van den Berg, 1981; Olney and Winter, 1985; Seth et al., 2011; van Ruijven and Weijs, 1990). For example, ‘Hill-type’ muscle models predict muscle force on the basis of activation dynamics, isometric force–length and isotonic force–velocity relationships (Dick et al., 2017; Wakeling et al., 2021). However, Hill models fail to predict the history dependence of muscle force (Perreault et al., 2003; Sandercock and Heckman, 1997) because they must produce the same force given the same length and activation (Libby et al., 2020). It is difficult to reconcile these models with our findings that muscles produce different forces depending on their length history and stiffness given the same length and activation.

In contrast to the view of ‘muscle as motor’, the view that muscle is a tunable material emphasizes the interaction between activation and length changes in producing muscle force. The view that muscle is a tunable material predicts the dependence of muscle force on length history. Libby et al. (2020) demonstrated that the force response of muscle to length perturbations is history dependent, suggesting that muscles retain a viscoelastic memory of their length history. They developed a simple viscoelastic model that used a spring, damper and equilibrium length to predict the muscle force response to perturbations. When the equilibrium length of the spring was proportional to the preload, the model predicted the muscle work during length perturbations with high accuracy (r²>0.99).

We developed a three-parameter viscoelastic model with constant stiffness and a two-parameter viscoelastic model with history-dependent stiffness equal to measured muscle stiffness. Both models included a linear spring and a linear damper in parallel (Fig. 10A). The three-parameter model included a linear spring constant k, a linear damping coefficient c, and equilibrium position x*. The two-parameter model included a linear damping coefficient c and equilibrium position x* but the linear spring constant was replaced with the measured history-dependent muscle stiffness. The parameters were optimized for the best r² value for length perturbations using data from all subjects in protocol 2, and then the model was used to predict each subject's data using the optimized parameter set. Work during length perturbations (Fig. 10B,C) was predicted accurately by both the model with constant spring stiffness (mean r²=0.91±0.08) and the model with history-dependent stiffness (mean r²=0.92±0.08) as in the previous study (Libby et al., 2020). However, the model that included history-dependent stiffness predicted peak force with higher accuracy (mean r²=0.77±0.05) than the model with constant stiffness (mean r²=0.46±0.25) (Fig. 10D,E). The simulation results support the hypothesis that muscle force under dynamic conditions can be accurately predicted as the force response of a history-dependent viscoelastic material with tunable stiffness to length perturbations.

Most muscle models have not considered muscle as a ‘tunable material’. Hill models predict muscle force well under near-isometric conditions. However, they would perform poorly at predicting the force response to perturbations in both protocols in our study, largely because the activation dynamics in Hill models lack the length and length history dependence of activation dynamics (Sandercock and Heckman, 1997; Shue and Crago, 1998). In principle, a titin-inspired model with a pulley representing titin–actin interactions (Dutta et al., 2018) could predict length-dependent activation dynamics because the pulley can rotate and translate simultaneously (Nishikawa and Huck, 2021).

In addition to tunable stiffness, tunable damping also likely contributes to the muscle force response to length perturbations and activation, and hence to in vivo locomotion. In the present study, muscle stiffness was emphasized and measured using quick transients. The method reduces or eliminates viscous effects. Therefore, no conclusions about damping can be drawn from this study. Future studies should explore the effect of tunable damping on muscle force response to length perturbations and activation.

Intrinsic muscle properties contribute importantly to control of movement, along with central nervous system feedforward and feedback mechanisms (Nishikawa et al., 2013). Reflexes act over longer time scales by modulating muscle activation in response to muscle stretch or unloading, although intrinsic muscle properties act on shorter time scales and appear to be important for maintaining stability in the face of perturbations (Nichols and Houk, 1976). The observation of decoupling between muscle activation and force during in vivo locomotion suggests the existence of sensitive periods for control, likely when conditions are closer to isometric. However, our findings do not directly address the roles of feedforward and feedback control in locomotion. Future ex vivo work loop experiments could explore this question by systematically varying stimulation at different phases of a strain cycle.

We thank Philip Service and the anonymous reviewers for comments that improved the manuscript.