Vol. 85, Issue 6, 2176-2189, December 1998
Two-step, predictive, isometric force model tested on data
from human and rat muscles
Jun
Ding1,
Stuart A.
Binder-Macleod1,2, and
Anthony S.
Wexler1,3
1 Interdisciplinary Graduate
Program in Biomechanics and Movement
Science, 2 Department of
Physical Therapy, and 3 Department
of Mechanical Engineering, University of Delaware, Newark, Delaware
19716
 |
ABSTRACT |
Functional
electrical stimulation can assist paralyzed individuals to perform
functional movements, but muscle fatigue is a major limitation to its
practical use. An accurate and predictive mathematical model can
facilitate the design of stimulation patterns that optimize aspects of
the force transient while minimizing fatigue. Solution nonuniqueness, a
major shortcoming in previous work, was overcome with a simpler model.
The model was tested on data collected during isometric contractions of
rat gastrocnemius muscles and human quadriceps femoris muscles under
various physiological conditions. For each condition tested, parameter
values were identified using the force response to one or two
stimulation trains. The parameterized model was then used to predict
forces in response to other stimulation patterns. The predicted forces
closely matched the measured forces. The model was not sensitive to
initial parameter estimates, demonstrating solution uniqueness. By
predicting the force that develops in response to an arbitrary pattern
of stimulation, we envision the present model helping identify optimal
stimulation patterns for activation of skeletal muscle during
functional electrical stimulation.
functional electrical stimulation; mathematical model; catchlike
property; fatigue; muscle length
| |
INTRODUCTION |
FUNCTIONAL ELECTRICAL stimulation (FES) can assist
paralyzed individuals to perform functional movement (18, 20). Muscle fatigue, however, is a major limitation to its practical implementation (18, 20). It has long been known that the stimulation frequency affects
the forces produced and the rate of muscle fatigue (8). More recently,
however, it has been shown that, in addition to the mean frequency, the
pattern of the pulses can markedly affect force production and fatigue
(4, 5). Studies showed that the activation pattern that produces the
greatest force varies from person to person (19) and, even for the same
individual, changes for different physiological conditions [e.g.,
muscle length (21)]. Therefore, it is important and difficult to
identify activation patterns that maximize force and minimize fatigue
for individual patients.
The development of a mathematical model that predicts the force
responses of muscles to electrical stimulation would enable clinicians
to identify the optimal activation patterns for individual subjects.
Several mathematical models of muscle force have been developed,
including pure mathematical models (9, 12), Huxley-type models (17, 26,
27, 32), and Hill-type models (11, 13, 16, 22, 31). Different aspects
of the muscle have been modeled. Our long-term goal is to develop a
model that can predict forces of muscle during FES. A reliable
isometric force model is a crucial step toward this goal. The only
model able to predict isometric force response to a wide range of
stimulation patterns incorporates the Hill paradigm (31). Hill-type
models are phenomenological, originating from experimental observations
of the muscle's mechanical behavior, and describe the whole muscle as
a contractile element in series with a spring, together with various
modifications such as dampers placed in parallel or series with other components.
In our previous work we decomposed the contractile response of skeletal
muscle into three steps: 1)
Ca2+ release and uptake by the
sarcoplasmic reticulum, 2)
Ca2+ reaction with troponin, and
3) force mechanics modeled by a
linear spring, damper, and motor in series (31). Steps
1 and 2 make the model
less phenomenological than Hill's original model (16). This three-step
model consisted of three ordinary differential equations and nine
parameters that were identified by fitting the model to forces measured
in response to a single constant-frequency train (CFT) stimulation.
Then, using the parameter values identified from this fit, we tested
the fidelity of the model by comparing its predictions with responses
measured from other stimulation patterns. For rat gastrocnemius and
soleus muscles, the model accurately fit and predicted the force
responses during stimulation with a wide range of stimulation
frequencies and patterns. The primary shortcoming of the model was its
sensitivity to different initial parameter estimates. Numerous fitting
trials were necessary to identify the initial parameter values, leading
to identified values that reliably predicted force responses. In
addition, subsequent analysis of these previously published results
showed that the model was only able to predict the stimulation pattern
that produced the greatest force-time integral for 9 of the 15 gastrocnemius muscles studied. All the muscles for which the model was
successful had twitch contraction times <25 ms. The purpose of this
study, therefore, was to develop a model with greater reliability and range of applicability.
The simpler model presented here consists of only two ordinary
differential equations and four parameters. It has been tested on
isometric contractions from rat and human muscles in response to a wide
range of brief, physiologically relevant stimulation patterns. Forces
used to test the model were collected from fresh and fatigued muscles
and at long (near-optimal) and short (suboptimal) muscle length.
| |
FORMULATION OF THE MODEL |
In the new muscle model the contractile response is decomposed into two
distinct physiological steps: cross-bridge activation and force generation.
Cross-bridge activation.
Pilot studies with our previous model (31) showed that the ability of
the model to predict forces is relatively insensitive to the specific
shape and magnitude of the Ca2+
and Ca2+-troponin complex
transient. This implies that steps 1 and 2 in our previous model could be
combined into one step. Therefore, we use one unitless factor,
CN, the dynamics
of which are governed by a time constant,
c, to qualitatively describe
the rate-limiting step before the myofilaments mechanically slide across each other and generate force (1; see
APPENDIX for details). The
differential equation describing this dynamic is
|
(1)
|
and
the chosen analytic solution is
|
(2)
|
where
the transient behavior of
CN characterizes
the rate-limiting step in force generation,
t (ms) is the time since the beginning
of the stimulation train,
ti (ms) is the
time of the ith stimulation pulse
since the beginning of the stimulation train, n is the number of stimulation pulses
before time t in the train (its
maximum value is 6 in this study), and
c (ms) is the time constant
controlling the transient shape of
CN.
Force generation.
Similar to our previous model (31), the force generation is modeled
mechanically by a linear spring, damper, and motor in series, which
represent the elastic, viscous, and contractile elements of the muscle,
respectively (Fig. 1). The damper
represents the viscous resistance of the contractile and connective
tissue (22). The force (F) exerted by the damper is given by
|
(3)
|
where
b is the damping coefficient,
x is the length of the spring, and
V is the contractile velocity of the
motor.
The motor represents the contractile element, or the sliding of actin
and myosin filaments, and is driven by strongly bound cross bridges
(1). Because there is a sigmoidal force-pCa relationship (2), we tested
for a similar relationship between
CN and the strongly bound cross bridges and found that a simple Michaelis-Menten function, CN/(1 + CN),
sufficed. The contractile capability of the muscle,
V, therefore, is
|
(4)
|
in
which B is a proportionality constant.
A linear spring represents the series elastic components of the muscle
and tendon. The force exerted by the spring is given by
|
(5)
|
where
K is the spring constant.
Differentiating Eq. 5 with respect to
time and using Eq. 3 to eliminate
dx/dt
and Eq. 4 to eliminate
V gives
|
(6)
|
The term b/K represents the time
constant over which the force decays. The friction between the actin
and myosin filaments increases with an increase in bound cross bridges,
so we set b/K = 1 + 2CN/(1 + CN), where
1 is the value of the time
constant in the absence of
CN and
2 is the additional frictional
component due to the cross-bridge chemical bonds. Substituting for
b/K and replacing
KB with a new constant,
A, in Eq. 6 gives
|
(7)
|
The
isometric force dynamics are governed by Eqs.
1 and 7, which
describe the transient behavior of the two state variables, CN and F,
subject to the four parameters
c,
A,
1, and 2 (Table
1).
| |
METHODS |
Experimental procedure.
Rat and human data were taken from previously published or unpublished
studies in our laboratory, which were approved by the University of
Delaware's Animal and Human Subjects Review Board. Data were collected
from 16 nonfatigued and 8 fatigued rat gastrocnemius muscles. The
experimental methods are similar to those previously described (31).
Briefly, the rats were deeply anesthetized, and a hook stimulating
electrode was positioned around the nerve leading to the gastrocnemius
muscle. Force was recorded directly from the Achilles' tendon. The
gastrocnemius muscle was stimulated with a series of nine six-pulse
trains. During traditional FES, skeletal muscle is activated with short
CFTs, where the pulses within each train are separated by regular
intervals. Five of the trains in our study, therefore, were CFTs with
interpulse intervals (IPIs) of 10, 20, 30, 40, or 50 ms (Fig.
2A):
CFT10, CFT20, CFT30, CFT40, and CFT50, respectively. The other four
trains were variable-frequency trains (VFTs), with an initial IPI of 10 ms, and the remaining four IPIs within the train were 20, 30, 40, or 50 ms (Fig. 2B): VFT20, VFT30, VFT40,
and VFT50, respectively. The VFTs were chosen, because recently it has
been shown that VFTs that take advantage of the catchlike property of
skeletal muscle can markedly increase the force output from fatigued
muscle compared with CFT stimulation (4, 5, 8). The catchlike property
is the tension enhancement seen when an initial brief high-frequency
burst of pulses (2-4 pulses) at the onset of a subtetanic CFT is
used to activate the muscle (6, 10) and is a fundamental muscle
property that is not attributable to the motor axon or neuromuscular
junction (4, 10). To collect the responses to nonfatiguing stimulation,
~10 s separated each train. To fatigue the muscle, 3 random sequences
of these 9 trains were used to form a cycle of 27 trains, and 8 such
cycles (216 trains) were used. Each train was separated by ~1 s. The
force responses to the last cycle of trains were used to represent the force responses from the fatigued muscles.
View larger version (12K):
[in this window]
[in a new window]
|
Fig. 2.
Schematic representation of constant-frequency trains (CFTs) and
variable-frequency trains (VFTs) used to activate rat gastrocnemius
muscles. Each vertical line represents a stimulation pulse.
A: CFTs with interpulse interval (IPI)
of 10-50 ms in 10-ms increments
(top to
bottom).
B: comparable VFTs with an initial IPI
of 10 ms.
|
|
The experimental setup for collection of human data was similar to that
previously described (8). Briefly, subjects were seated on a force
dynamometer with their hips flexed to ~75° and their knees flexed
to 90° for the long-muscle-length study or 15° for the
short-muscle-length study. Two stimulating electrodes were placed over
the quadriceps femoris muscle. The force transducer was positioned
against the anterior aspect of the leg, proximal to the lateral
malleolus. Nonfatiguing and fatiguing protocols, similar to those used
for rat muscles, were used for human subjects, except the initial IPI
of the VFTs was 5 ms, and the IPIs of the CFTs were 10-120 ms in
steps of 10 ms, for a total of 24 stimulation patterns. Nonfatiguing
force responses were collected for long-muscle-length (n = 12) and short-muscle-length
(n = 11) studies. Fatigued force data
were collected only at long muscle length
(n = 11).
Simulations.
Reliability of the model was first investigated for rat and human
muscles by using the model to fit the force data for each animal or
subject by incorporation of a wide range of initial parameter values.
Convergence to the same parameter values for each rat or human subject
under each condition would suggest a unique solution of the model
within the range of initial parameter values tested.
Parameter identification.
SAAM II, a modeling system developed by the University of Washington
and the National Institutes of Health, was used to identify the optimum
parameter values and predict forces for other stimulation trains (24).
For each muscle, 1 was first
calculated using the numerical module of SAAM II. The procedure to
calculate 1 was described in
detail previously (31). Briefly, it was found by setting
CN = 0, so that
Eq. 7 became
dF/dt = 
F/1. By taking the force
decay at the end of the tetany and performing a linear regression of
lnF vs. t, we obtained the value of
1 from the slope of the fit.
The other three parameter values were then identified using the
compartment module of SAAM II, which accepted as input 1) Eqs.
1 and 7 expressed
compartmentally, 2) force measured
for a contraction, 3) the
stimulation pattern, and 4) initial
parameter estimates. SAAM II performs a multidimensional fit between
the measured and predicted force data to identify the best parameter values.
Preliminary studies were conducted to find the best train(s) for
parameter identification. For rat and human muscles, trains that
produce nonfused tetanic force responses were tested first, because
they exhibited a clear rise and fall of force for each pulse and a
buildup of force over multiple pulses. CFTs with IPIs of 30 ms (CFT30,
Fig. 3C)
and 70 ms (CFT70, see Fig. 6C) were initially used to fit the model for rat gastrocnemius and human quadriceps femoris muscles, respectively.
View larger version (23K):
[in this window]
[in a new window]
|
Fig. 3.
Measured (exp) and predicted (mod) forces from a typical nonfatigued
rat gastrocnemius muscle. A-E:
6-pulse CFTs of 10, 20, 30, 40, and 50 ms.
F-I: comparable VFTs. Model
parameter values were obtained by fitting model to responses to CFT20
(B) and CFT50
(E).
|
|
For rat nonfatigued or fatigued muscles, the 30-ms CFT did not always
lead to a set of parameter values that successfully predicted the
forces from other stimulation patterns. In fact, no other train was
found that could satisfactorily predict forces from all other
stimulation patterns. Further testing showed that only a combination of
two stimulation trains, one that produced a nearly fused tetanic
response (i.e., CFT20) and another that produced an unfused tetanic
response (i.e., CFT50), identified a reliable set of parameter values
(see Fig. 3 for examples of force responses to CFT20 and CFT50). The
CFT20 and CFT50, therefore, were used to parameterize the model for the
rat muscles.
Unlike the rat muscles, parameter values for human quadriceps femoris
muscles obtained by fitting the combination of a high- and a
low-frequency CFT (e.g., CFT30 and CFT100) did not improve the accuracy
of the predicted responses. When muscles were at a long length and not
fatigued, the CFT with IPI of 70 ms (see Fig.
6C) yielded parameter values that
best predicted the forces. However, when muscles were at long length
and fatigued or at short length and not fatigued, we found that the
CFT70 could only be used to identify parameter values to predict forces
to other CFTs and that the VFT70 (see Fig.
7H) best identified parameter values for the VFTs.
Force predictions.
For each muscle, with use of the parameterized model, the state
variable CN and
the force responses to other stimulation patterns were predicted by
simply changing the simulation pattern input to the parameterized
compartment module. We then tested the fidelity of the model by
comparing the predicted forces with the measured forces.
Data analysis.
To validate the model, several comparisons were made between the
predicted and measured data. To compare the shape of the forces, a
correlation coefficient was calculated by comparing the measured and
predicted forces at 2-ms-interval time steps. Any phase shift or
disagreement in shape between the predicted and measured data lowers
the correlation coefficient. Next, for each animal or subject, the
force-time integral was compared for the measured and the predicted
responses to each stimulation pattern. For each condition, the measured
force-time integrals were plotted against the predicted ones. Linear
regression trend lines, with intercepts of zero, were used to compare
force response across all the stimulation patterns and determine how
well the model predicts the force-time integrals. A perfectly accurate
model would have a slope and correlation coefficient of 1. Also within each condition, percent differences were calculated between the average
measured and predicted force-time integral for each stimulation pattern. If the percent differences were greater than a chosen criterion, 5%, a paired t-test was
used to check whether the difference was significant across the
subjects. We chose 5%, because we believed that this was the minimal
difference that would be clinically significant.
The best-fit parameter values were pooled for all subjects within each
condition, and means ± SD were calculated. For each parameter, a
paired t-test was used to compare the
parameter values for the same animal or subject across conditions, with
nonfatigued and long-muscle-length-nonfatigued data as the control
condition for rat and human studies, respectively.
| |
RESULTS |
Pilot studies on reliability of the model showed convergence of the
parameter values within the range of values tested for each of the four
rat and six human muscles tested, suggesting that the model was able to
fit the data with a unique solution. Subjective evaluation of the force
responses showed that the model accurately fit and predicted the
responses to a wide range of stimulation patterns for the nonfatigued
and fatigued rat muscles (Figs. 3 and 4,
respectively). The median correlation coefficients were 0.98 for the
nonfatigued and fatigued muscles (Table 2), which is consistent with the subjective impression that the predicted response accurately matched the shape of the measured data. Typical trend lines for rat gastrocnemius force-time integrals had slope and
correlation coefficient close to 1 (Fig. 5,
A and
D), which also suggested an accurate
fit and prediction with the model. For nonfatigued rat muscles, the
model significantly overestimated the CFT10 force-time integral by 16%
(Fig. 5B). For fatigued rat muscles,
although the model significantly overestimated the VFT30, VFT40, and
VFT50, the differences between the measured and predicted forces were
small (8.1, 6.7, and 8.5%, respectively; Fig.
5F). The model successfully
predicted the activation patterns that produced the greatest force-time
integrals for seven of the nine rat muscles with twitch contraction
times <25 ms and for seven of the eight muscles with twitch
contraction times >25 ms.
View larger version (23K):
[in this window]
[in a new window]
|
Fig. 4.
Measured and predicted forces from a typical fatigued rat gastrocnemius
muscle (same muscle presented in Fig. 3).
A-E: 6-pulse CFTs of 10, 20, 30, 40, and 50 ms. F-I: comparable
VFTs. Model parameter values were obtained by fitting model to
responses to CFT20 (B) and CFT50
(E).
|
|
View larger version (25K):
[in this window]
[in a new window]
|
Fig. 5.
Comparison of force-time integrals (N · ms) between
measured and predicted forces of 16 nonfatigued
(A-C) and 8 fatigued
(D-F) rat gastrocnemius
muscles. A and
D: plots of predicted vs. measured
force-time integrals for 9 CFTs and VFTs.
B and
E: force-time integrals (means ± SE) vs. stimulation type for 5 CFTs. C
and F: force-time integrals (means ± SE) vs. stimulation type for 4 VFTs. Significant difference
between predicted and measured time-force integrals for each
stimulation pattern: * P < 0.05, ** P < 0.01.
|
|
For human studies, subjective evaluation of the force responses also
showed that, in general, the model accurately fit and predicted the responses to a wide range of stimulation patterns and
conditions (Figs.
6-8).
The correlation coefficients between the measured and predicted forces,
averaged across subjects, had a median of 0.96, 0.96, and 0.97 (Table
2) for nonfatigued long-length, fatigued long-length, and nonfatigued
short-length conditions, respectively, suggesting good agreement
between the shape of the measured and predicted forces.
View larger version (23K):
[in this window]
[in a new window]
|
Fig. 6.
Measured and predicted forces from a typical subject's nonfatigued
quadriceps femoris muscle held at long length (90° knee flexion).
Model parameter values were obtained by fitting model to CFT70
(C).
|
|
View larger version (23K):
[in this window]
[in a new window]
|
Fig. 7.
Measured and predicted forces from a fatigued human quadriceps femoris
muscle held at long length (same subject presented in Fig. 6). Model
parameter values obtained by fitting model to CFT70
(C) were applied to other CFT
stimulation patterns; parameter values obtained by fitting VFT70
(H) were applied to other VFT
stimulation patterns.
|
|
View larger version (24K):
[in this window]
[in a new window]
|
Fig. 8.
Measured and predicted forces from a human quadriceps femoris muscle
held at short length (same subject presented in Fig. 6). Model
parameter values obtained by fitting model to CFT70
(C) were applied to other CFT
stimulation patterns; parameter values obtained by fitting model to
VFT70 (H) were applied to other VFT
stimulation patterns.
|
|
When muscles were not fatigued and were held at long length, the trend
lines for force-time integral had slopes close to 1, and the model
accounted for 65-68% of the variability in the measured data
(Fig. 9). The model significantly
overestimated the force produced in response to the CFTs and VFTs with
short IPIs (10-60 ms) and underestimated those with long IPIs
(80-110 ms).
View larger version (39K):
[in this window]
[in a new window]
|
Fig. 9.
Comparison of force-time integrals between measured and predicted
forces of 12 nonfatigued human quadriceps femoris muscles held at long
length. A and
C: plots of predicted vs. measured
force-time integral for 12 CFTs and 12 VFTs, respectively, analyzed for
each of 12 muscles. B: force-time
integrals (means ± SE) vs. stimulation type for CFTs
(CFT10-CFT120). D: force-time
integrals (means ± SE) vs. stimulation type for VFTs
(VFT10-VFT120). Significant difference between measured and
predicted force-time integrals for each stimulation pattern:
* P < 0.05, ** P < 0.01.
|
|
Surprisingly, when muscles were fatigued the model did a better job of
predicting the force output than when muscles were not fatigued. The
slopes of force-time integral trend lines were only slightly <1, and
the model accounted for 90-96% of the variability in the measured
data (Fig. 10). The model tended to
underestimate the force-time integrals of CFTs with short IPIs
(20-40 ms), the percent error being worst (25%) for CFT10. The
model slightly overestimated the force-time integrals of the VFTs with
short IPIs (20-40 ms) and underestimated the VFTs with long IPIs
(90-110 ms).
View larger version (38K):
[in this window]
[in a new window]
|
Fig. 10.
Comparison of force-time integrals between measured and predicted
forces of 11 fatigued human quadriceps femoris muscles held at long
length. A and
C: plots of predicted vs. measured
force-time integral for 12 CFTs and 12 VFTs, respectively, analyzed for
each of 11 muscles. B: force-time
integrals (means ± SE) vs. stimulation type for CFTs
(CFT10-CFT120). D: force-time
integrals (means ± SE) vs. stimulation type for VFTs
(VFT10-VFT120). Significant difference between measured and
predicted force-time integrals for each stimulation pattern:
* P < 0.05, ** P < 0.01.
|
|
When muscles were not fatigued and were at short length, trend lines of
the force-time integral had slope values close to 1 and correlation
coefficients of 0.89-0.94 (Fig. 11).
The model significantly overestimated the force-time integral of VFTs
with short IPIs (20-40 ms; Fig. 11). However, these differences
were <8%.
View larger version (37K):
[in this window]
[in a new window]
|
Fig. 11.
Comparison of force-time integrals between measured and predicted
forces of 11 nonfatigued human quadriceps femoris muscles held at short
length. A and
C: plots of predicted vs. measured
force-time integral for 12 CFTs and 12 VFTs, respectively, analyzed for
each of 11 muscles. B: force-time
integrals (means ± SE) vs. stimulation type for CFTs
(CFT10-CFT120). D: force-time
integrals (means ± SE) vs. stimulation type for VFTs
(VFT10-VFT120). Significant difference between measured and
predicted force-time integrals for each stimulation pattern:
* P < 0.05.
|
|
Overall, the model did an excellent job of predicting the pattern of
activation that produced the maximum force-time integral. The only
exceptions were the human nonfatigued data collected at long muscle
length (Fig. 9). In this condition the error between the maximum
predicted and measured force-time integrals was, however, only ~5%.
The model predicted that the maximum force-time integral is obtained
with IPIs that are short compared with the measurements, leading to an
actual force-time integral that is 10% lower than the measured maximum.
Fatigue produced changes in all four parameter values that were
consistent for rat and human muscles (Table 3). Paired
t-test showed a significant increase
in 1 and a significant decrease in A when muscles were fatigued or
held at short length for rats and humans. Only for humans, however,
were significant differences seen for
2 and
c across conditions.
| |
DISCUSSION |
The current model predicted isometric forces for rat gastrocnemius and
human quadriceps femoris muscles with reasonable accuracy during brief
subtetanic and tetanic isomeric contractions by using CFTs and VFTs. It
successfully predicted nonfatigued force responses when muscles were
held at long or short length and fatigued force responses at long
length. Subjective and objective evaluations of the force response
showed that the shape of the predicted force responses closely matched
the responses that were measured. Comparisons between the predicted and
measured force-time integrals also suggested excellent agreement
between the predicted and measured forces. Although significant
differences between the predicted and measured force-time integrals
were observed with some stimulation patterns, most differences were small.
Comparison with previous work.
The present model differs from most previously developed models, in
that it is simpler in structure and has fewer parameters (11, 17, 26,
31). With only four free parameters, the model allows precise and
unique parameter estimation from the force data. This conclusion was
based on the convergence to the same parameter values starting from
different initial estimates. In addition, this model is the first to be
tested on muscles under different physiological conditions.
The process of muscular contraction at the molecular level is very
complicated. Many steps, however, may not contribute significantly to
force generation and would not be detectable from analysis of force
data alone (27). Previous investigators have tried other approaches to
a simpler model. Schultz and colleagues (25) developed a simple
Hill-type model, also with four independent parameters. Their model was
designed to predict the time course of the muscle force during
isometric tetanic and isovelocity shortening contractions. Unlike the
model presented in this study, however, the model developed by Schultz
and colleagues cannot predict the force-time response to unfused
tetanic, subtetanic, or variable-frequency stimulations. Recently,
Shames and colleagues (26) developed a simple version of a Huxley-type
model, also with four parameters, to explore the
Ca2+-force relationship during
isometric tetanic contractions. With use of the
Ca2+ transient as input, their
model successfully fit the forces produced during isometric tetanic
contractions of various durations and isometric twitch contractions.
However, this model has limited predictability, because the
Ca2+ transient is not typically available.
An advantage of the present model over our previously reported model
(31) is that it can predict force responses from fatigued muscles. The
previous model failed to predict the forces from fatigued muscles by
fitting the force in response to a single CFT. Subsequent studies with
our previous model also showed that using the responses to more than
one frequency of stimulation (e.g., the force produced in response to a
high- and a low-frequency train, see
METHODS,
Stimulations) to identify the
parameter values of the model did not improve the accuracy of the model
to predict forces. The present model was also successful with human
muscles at long and short lengths, although the short-muscle-length
condition was not tested on the previous model. The quality of the
present model's predictions was independent of the twitch contraction times of the muscles. In contrast, the previous model failed to predict
forces of all six nonfatigued rat gastrocnemius muscles with
twitch-contraction times >25 ms.
Parameter interpretation.
The advantage of Hill- and Huxley-type models over pure mathematical
models is that the parameter values may be physiologically meaningful
and the interpretation of the changes of those values with muscle
condition can provide useful insight into the mechanisms of muscle
contraction and fatigue (12). The parameters of this model, derived
from Hill's original model and based on the mechanisms of muscle
contraction, may also be helpful in elucidating characteristics of
muscle physiology.
Of the four model parameters, A is
directly responsible for the magnitude of the force. When the muscle is
fatigued or at short length, the force declines. Thus a decrease in the
values of A for fatigued muscle or
muscle placed at short length is expected. Our results are consistent
with this expectation (Table 3).
The relaxation of the force is primarily determined by
1, because the contribution of
2 diminishes rapidly with the
decline of CN,
which always precedes the force relaxation (23). Fitch and McComas (14)
observed a slower muscle relaxation in fatigued muscles, whereas
Gandevia and McKenzie (15) observed faster relaxation of nonfatigued
shortened human muscle. The increases and decreases in the fitted
values of 1 (Table 3) for
fatigued and short-length muscles, respectively, are consistent with
these previous findings.
The parameter c is the
characteristic time for
CN to reach its
peak and for it to decay to zero. The interpretation of the changes in
c across conditions is
difficult. Many underlying physiological steps, such as
Ca2+ release from and uptake by
the sarcoplasmic reticulum, Ca2+
binding to and dissociation from troponin, and cycling between the
strongly bound and weakly bound cross bridges, may affect the time
course of CN.
Also, it is unclear which of the steps is rate limiting (27). Until we
have a more complete understanding of cross-bridge activation, we
cannot properly interpret c. However, some previous experimental findings could, in part, explain our results. The prolonged Ca2+
transient with fatigued muscles observed by Westerblad and Allen (30)
could be one of the causes for an increased
c for fatigued muscles (Table
3). Reduced Ca2+ sensitivity to
troponin (3, 28) when muscles are at short length could cause less
troponin-Ca2+ complex to be formed
and, consequently, lead to activation of fewer cross bridges. Thus the
decrease in c with muscles at
short length could also be explained.
As the characteristic time in the force equation (Eq. 7), the sum of
1 and some fraction of
2 determines the overall rate of rise of force. When muscles were activated with CFTs, Binder-Macleod et al. (8) observed a slowing in the rate of rise of force as the
muscles fatigued. In contrast, there was no difference in the rise time
between the fatigued and the nonfatigued muscles when they were
stimulated with VFTs (8). The increases of
1 and
2 for CFTs and increase of
1 and decrease of
2 for VFTs as the muscles
become fatigued (Table 3) are, therefore,
consistent with previous findings.
Model limitation.
The major limitation of this model is that, for human muscles at short
length and fatigued conditions, separate sets of parameter values were
needed for CFTs and VFTs. In contrast, when human muscles were held at
long lengths and not fatigued, only one set of parameter values was
needed. However, the model did not always identify the stimulation
pattern that produced the greatest force-time integral for this
condition (Fig. 9). For this condition, CFT60 and VFT60 were identified
as producing the greatest force-time integral. In contrast, the
measured data showed that CFT90 and VFT100 produced the greatest
force-time integrals among CFTs and VFTs, respectively. Similarly, the
model consistently underestimated the forces at low frequencies (long
IPIs) and overestimated the forces at high frequencies (short IPIs;
Figs. 6 and 9). All these shortcomings could be due to some
physiological steps not being taken into consideration, such as unequal
amounts of Ca2+ released in
response to each pulse of the train.
Conclusion.
The mathematical model presented here is unique in its ability to
predict isometric forces from rat and human skeletal muscles in
nonfatigued and fatigued conditions at various muscle lengths. Modification of the present model are needed to allow a single set of
parameter values to predict forces produced by CFTs and VFTs and to
improve the model's accuracy under various physiological conditions
such as fatigued or nonfatigued muscle and short or long muscle length.
In addition, the model will need to be modified to predict forces
during repetitive stimulation of muscles that are allowed to shorten or
lengthen if the model is to be used to predict forces during FES. By
predicting the force that develops in response to an arbitrary pattern
of stimulation, we envision mathematical models similar to the present
one helping identify the optimal stimulation patterns for activation of
skeletal muscle during FES.
| |
APPENDIX |
Muscle force has been assumed to be proportional to the amount of
Ca2+-troponin complex in some
models (11, 27, 31) or the amount of strongly bound cross bridges
in others (26, 27, 29), suggesting that these factors play a similar
role in muscle activation. Therefore, we used one unitless factor,
C, to represent the single rate-limiting step before the myofilaments mechanically slide across
each other. No experimental data have elucidated the dynamics of the
Ca2+-troponin complex or the
strongly bound cross-bridge transients during muscle contraction.
However, studies of Chou and Hannaford (11) and Wexler et al. (31)
suggest a Ca2+-transientlike curve
for the dynamics of the
Ca2+-troponin complex. Therefore,
we assumed a similar transient shape for
C, which can be expressed as
|
(A1)
|
where
n is the number of pulses before
time t in the train and has a
maximum value of 6 in this study,
ti represents
the times that the electrical stimulation was applied, and
c1 and
c2 are the time
constants controlling the shape and height of
C.
Pilot studies with this model showed that the force production was not
sensitive to the magnitude of C,
because the effect of increased height of
C could be compensated for by
adjustment of parameter values in the force equation. Therefore, we
normalized Eq. A1 by dividing the
right and left side by
c2/(c1 + c2), which gives
|
(A2)
|
where
CN represents
the normalized C and is equal to
C(c1 + c2)/c2.
Pilot studies showed a relationship between
c1 and
c2 and that the
model worked best when
c1 was much larger
than c2. Therefore,
Eq. A2 can be further simplified to
|
(A3)
|
There
are a number of linearly dependent analytic solutions to
Eq. A3. The one chosen to represent
the transient behavior of
CN is
|
(A4)
|
because,
among the choices tested, this solution provided the best fits to the
force data. Examples of the changes in
CN are shown in
Fig. 12.
View larger version (16K):
[in this window]
[in a new window]
|
Fig. 12.
Typical predicted responses of
CN to 3 types of
stimulation patterns. Thin solid line, CFT100 (CFT with 100-ms IPIs);
dashed line, VFT100 (VFT with initial IPI of 5 ms and remaining pulses
equally spaced by 100 ms); thick solid line, CFT20 (CFT with 20-ms
IPIs). c = 20 for all 3 responses.
|
|
| |
ACKNOWLEDGEMENTS |
The authors thank Samuel C. K. Lee, Lisa Landis, April Fritz,
Lorin Kucharki, and Michelle Gerdom for assistance in data collection.
| |
FOOTNOTES |
This study was supported by National Institutes of Health Grant
AR/HD-41254.
The costs of publication of this
article were defrayed in part by the
payment of page charges. The article
must therefore be hereby marked
"advertisement"
in accordance with 18 U.S.C. §1734 solely to indicate this fact.
Address for reprint requests: S. A. Binder-Macleod, 315 McKinly
Laboratory, University of Delaware, Newark, DE 19716.
Received 5 May 1998; accepted in final form 31 July 1998.
| |
REFERENCES |
1.
Aidley, D. J.
The Physiology of Excitable Cells. New York: Cambridge University Press, 1989.
2.
Allen, D. G.,
and
J. R. Blinks.
The interpretation of light signals from aequorin-injected skeletal and cardiac muscle cells: a new method of calibration.
In: Detection and Measurement of Free Ca2+ in Cells, edited by C. C. Ashley,
and A. K. Campbell. Amsterdam: Elsevier/North-Holland, 1979, p. 159-164.
3.
Balnave, C. D.,
and
D. G. Allen.
The effect of muscle length on intracellular calcium and force in single fibers from mouse skeletal muscle.
J. Physiol. (Lond.)
492:
705-713,
1996[Medline].
4.
Bevan, L.,
Y. Laouris,
R. M. Reinking,
and
D. G. Stuart.
The effect of the stimulation pattern on the fatigue of single motor units in adult cats.
J. Physiol. (Lond.)
449:
85-108,
1992[Abstract/Free Full Text].
5.
Binder-Macleod, S. A.,
and
C. B. Barker.
Use of a catchlike property of human skeletal muscle to reduce fatigue.
Muscle Nerve
13:
850-857,
1991.
6.
Binder-Macleod, S. A.,
and
H. P. Clamann.
Force output of cat motor units stimulated with trains of linearly varying frequency.
J. Neurophysiol.
61:
208-216,
1989[Abstract/Free Full Text].
7.
Binder-Macleod, S. A.,
E. E. Halden,
and
K. A. Jungles.
Effects of stimulation intensity on the physiological responses of human motor units.
Med. Sci. Sports Exerc.
27:
556-565,
1995[Medline].
8.
Binder-Macleod, S. A.,
S. C. K. Lee,
and
S. A. Baadte.
Reduction of the fatigue-induced force decline in human skeletal muscle by optimized stimulation trains.
Arch. Phys. Med. Rehabil.
78:
1129-1137,
1997[Medline].
9.
Bobet, J.,
R. B. Stein,
and
M. N. Oguztöreli.
A linear time-varying model of force generation in skeletal muscle.
IEEE Trans. Biomed. Eng.
40:
1000-1006,
1993[Medline].
10.
Burke, R. E.,
P. Rudomin,
and
F. E. Zajac.
Catch property in single mammalian motor units.
Science
168:
122-124,
1970[Abstract/Free Full Text].
11.
Chou, C. P.,
and
B. Hannaford.
Dual stable point model of muscle activation and deactivation.
Biol. Cybern.
66:
511-523,
1992[Medline].
12.
Donaldson, N. N.,
H. Gollee,
K. J. Hunt,
J. C. Jarvis,
and
M. K. N. Kwende.
A radial basis function model of muscle stimulated with irregular inter-pulse intervals.
Med. Eng. Phys.
16:
431-441,
1995.
13.
Durfee, W. K.,
and
K. I. Palmer.
Estimation of force-activation, force-length, and force-velocity properties in isolated, electrically stimulated muscle.
IEEE Trans. Biomed. Eng.
41:
205-216,
1994[Medline].
14.
Fitch, S.,
and
A. McComas.
Influence of human muscle length of fatigue.
J. Physiol. (Lond.)
362:
205-213,
1985[Abstract/Free Full Text].
15.
Gandevia, S. C.,
and
D. K. McKenzie.
Activation of human muscles at short muscle lengths during maximal static efforts.
J. Physiol. (Lond.)
407:
599-613,
1988[Abstract/Free Full Text].
16.
Hill, A. V.
Heat and shortening and the dynamic constants of muscle.
Proc. R. Soc. Lond. B Biol. Sci.
125:
136-195,
1938.
17.
Huxley, A. F.
Muscle structure and theories of contraction.
Prog. Biophys. Mol. Biol.
7:
257-318,
1957.
18.
Isakov, E.,
J. Mizrahi,
and
J. Najenson.
Biomechanical and physiological evaluation of FES-activated paraplegic patients.
J. Rehabil. Res. Dev.
23:
9-19,
1986[Medline].
19.
Karu, Z. Z.,
W. K. Durfee,
and
A. M. Barzilai.
Reducing muscle fatigue in FES applications by stimulation with N-let pulse trains.
IEEE Trans. Biomed. Eng.
42:
809-816,
1995[Medline].
20.
Kralj, A.,
T. Bajd,
and
R. Turk.
Enhancement of gait restoration in spinal injured patients by functional electrical stimulation.
Clin. Orthop.
233:
34-43,
1988.
21.
Lee, S. C. K., M. L. Cullen, and S. A. Binder-Macleod. Effects of muscle length on the catchlike
property in fresh and fatigued human muscle (Abstract).
Med. Sci. Sports Exerc. 28, Suppl. 5: S130, 1996.
22.
Martin, A.,
L. Martin,
and
B. Morton.
Theoretical and experimental behavior of the muscle viscosity coefficient during maximal concentric actions.
J. Appl. Physiol.
69:
539-544,
1994.
23.
Rüegg, J. C.
Calcium in Muscle Activation: A Comparative Approach. New York: Springer-Verlag, 1986.
24.
SAAM II Institute.,
University of Washington,
and
National Institutes of Health.
SAAM II Seattle: University of Washington, 1994.
25.
Schultz, A. B.,
J. A. Faulkner,
and
V. A. Kadhiresan.
A simple Hill element-nonlinear spring model of muscle contraction biomechanics.
J. Appl. Physiol.
70:
803-812,
1991[Abstract/Free Full Text].
26.
Shames, D. M.,
A. J. Baker,
M. W. Weiner,
and
S. A. Camacho.
Ca2+-force relationship of frog skeletal muscle: a dynamic model for parameter estimation.
Am. J. Physiol.
271 (Cell Physiol. 40):
C2061-C2071,
1996.
27.
Stein, R. B.,
J. Bobet,
M. N. Oguztöreli,
and
M. Fryer.
The kinetics relating calcium and force in skeletal muscle.
Biophys. J.
54:
705-716,
1988[Abstract/Free Full Text].
28.
Stephenson, D. G.,
and
D. A. Williams.
Effects of sarcomere length on the force-pCa relation in fast- and slow-twitch skinned muscle fibers from the rat.
J. Physiol. (Lond.)
333:
637-653,
1982[Abstract/Free Full Text].
29.
Taylor, C. P. S.
Isometric muscle contraction and the active state: an analog computer study.
Biophys. J.
9:
759-780,
1969.
30.
Westerblad, H.,
and
D. G. Allen.
Changes in myoplasmic calcium concentration during fatigue in single mouse muscle fibers.
J. Gen. Physiol.
98:
615-635,
1991[Abstract/Free Full Text].
31.
Wexler, A. S.,
J. Ding,
and
S. A. Binder-Macleod.
A mathematical model that predicts skeletal muscle force.
IEEE Trans. Biomed. Eng.
44:
337-348,
1997[Medline].
32.
Zahalak, G. I.,
and
S. P. Ma.
Muscle activation and contraction: constitutive relations based directly on cross-bridge kinetics.
J. Biomech. Eng.
112:
52-62,
1990[Medline].
J APPL PHYSIOL 85(6):2176-2189
8570-7587/98 $5.00
Copyright © 1998 the American Physiological Society
Top
Abstract
Introduction
