A simple Hill element-nonlinear spring model of muscle contraction biomechanics ALBERT B. SCHULTZ, JOHN A. FAULKNER, AND V. A. KADHIRESAN Departments of Mechanical Engineering and Applied Mechanics and Physiology Program, University of Michigan, Ann Arbor, Michigan 48109-2125 SCHULTZ,ALBERTB., JOHN A. FAULKNER,AND V.A. KADHIRESAN. A simple Hill element-nonlinear spring model of muscle contraction biomechanics. J. Appl. Physiol. 70(Z): 803-812, 1991.-The purpose of this study was to develop a model to predict the mechanical responseof musclesduring isometric tetanic, afterloaded isotonic and isovelocity shortening contractions. Two versionsof the modelwere developed.Both incorporated a contractile element that obeyed a Hill force-velocity relationship and a serieselastic element. In a quadratic spring version, the serieselastic element force wasrepresentedasproportional to the squareof the stretch; in a cubic spring version, it was represented as proportional to the cube of the stretch. Both versions provided closed-form equations for response predictions that involved four independent parameters. Once the four parameterswere chosen,eachof these responsescould be predicted. Model validity was establishedby comparing predicted and observed responsesin slow and fast hindlimb muscles of rodents. Significant model-predicted responsesseldom differed by more than 15% from experimental values. The modelcan provide insights into how changesin individual properties affect the overall mechanical behavior of musclesin a variety of circumstancesand reduce the need for collection of experimental data. isometric contractions; isotonic contractions; isovelocity contractions; power; stiffness

MATHEMATICAL MODELS have the potential to predict the biomechanical behaviors on contracting skeletal muscles in a variety of circumstances. Models that represent biomechanical behaviors in vivo and in vitro with reasonable accuracy can provide insights into what factors determine those behaviors. When provided with experimental data from just a few tests, a model can predict what will happen in a array of other tests so that those tests need not be actually conducted. Consequently, use of a model can diminish the need for experimental measurements. Many investigators have constructed biomechanical models to describe various aspects of muscle behavior. MacMahon (14) reviews a number of these models. They have been designed primarily to describe cross-bridge mechanics (12, 17) or complex whole body segment biomechanics (l&20). Few studies have aimed explicitly at the development of models to explain biomechanical behaviors of single whole muscles during different types of contractions. The purpose of the present study was to develop and to test the validity of a mathematical model to predict biomechanical behaviors of single muscles 0161-7567/91

$1.50

and Bioengineering

studied in vitro and in situ during isometric, afterloaded isotonic, and isovelocity shortening contractions. Abilities to predict small-perturbation muscle stiffnesses and power generation capabilities were also examined. Model validity was tested by comparing model-predicted behaviors with experimentally observed behaviors of extensor digitorum longus (EDL) and soleus (Sol) muscles of mice and rats. METHODS

Construction of Biomechanical Model The model developed (Fig. 1) consisted of a contractile element (CE) is series with an elastic element (EE). For simplicity, the behavior of the CE was assumed to be that of the Hill representation of the force-velocity relationship (8), although this is a simplification of observed force-velocity relationships (3, 4). The Hill equation can be written in the form VW,

= (alF,)(l

- FIF,)I(aIF,

+ F/F,)

(0

where F is the contraction force, V is the velocity of the nonfixed endpoint of the CE, F, and Vm are the largest values of F and V measured under a given set of conditions for muscle activation, and a/F, is the measure of curvature of the force-velocity relationship. The elastic element was a nonlinear spring. Governing equations are derived here for the case of a quadratic elastic spring which obeys the force-displacement relationship F/F,

=

[W - W&l2

(2)

where X and Y are the displacement of the two ends of the spring (Fig. 1) and D, is the stretch of that spring when subject to a tension force of F,. In the APPENDIX, analogous governing equations are derived for a cubic elastic spring. We also examined use of both linear and square-root springs for the model EE but found neither of these models to predict behaviors accurately enough (We DISCUSSION).

The velocity of the nonfixed end of the EE and the contraction velocity of the CE, which are both equal to V, is the time derivative of X. The tangent stiffness of a muscle is the ratio of force change to length change when small perturbations in muscle length are made. The stiffness of the model muscle is a function of the stretch produced in the EE by the muscle contraction. The relation of the tangent stiffness

Copyright 0 1991 the American Physiological

Society

803

Downloaded from www.physiology.org/journal/jappl by ${individualUser.givenNames} ${individualUser.surname} (163.015.154.053) on October 4, 2018. Copyright © 1991 American Physiological Society. All rights reserved.

804

CONTRACTION

BIOMECHANICS

X

Y

r

Elastic Element

Contractile Element

I

FIG. 1. Schematic of model showing contractile and elastic elements in series.

to the contraction force is found to taking an X derivative of Eq. 2 to obtain dFldX

= 2(F,lD,)

[ (FlFm)o*5]

(3)

MODEL

where q2 = (F/F,). Equation 7 relates the muscle contraction force F to time t. Choosing a value for F and solving Eq. 7 for the corresponding t is more convenient than choosing t and solving for F. The rate of force development at any time in an isometric contraction can be calculated by finding the value of F corresponding to that time from Eq. 7 and then finding the corresponding dF/dt from Eq. 6. The tangent stiffness at any time in an isometric response can be calculated by finding the value of F corresponding to that time from Eq. 7 and then finding the corresponding dFldX from Eq. 3.

Isometric Contraction

Afterloaded Isotonic Contraction

Consider the response of this model in an isometric tetanic contraction. Both origin and insertion are held fixed in an isometric test, so the Y displacement of the spring is zero in this circumstance. Equation 2 then becomes

Two phases of responses are observed during an afterloaded isotonic contraction. The first phase ends at the time, defined as TB,at which the contraction force rises to the value of the afterload. Before time TBis reached the contraction is isometric. Consequently, the response during this phase is given by Eq. 7. If the afterload F is substituted into that equation, the value of TB can be calculated. During the second phase of the afterloaded isotonic contraction, beyond T,, the muscle insertion point is assumed to move at a constant velocity. The velocity is a function of the magnitude of the load and is derived from the Hill equation. Therefore the model-predicted quantity of chief interest is the start time for displacement, T*.

F/F,

Taking eliminate

= (XIDJ2

(4)

a time derivative of Eq. 4 and using Eq. 4 to from the result produces

X/D,

VW, = (T,/2F,)(dFldt)[(FlF~)‘-0~5’] (5) where a characteristic response time, T,, has been defined as equal to 0,/V,. T, is, for example, the time required for the force to rise to approximately one-third its maximum value during an isometric tetantic concentration. Replacing the left-hand side of Eq. 5 with Eq. 1 yields the differential equation governing isometric response dFldt = 2(F,IT,)(a/F~)[(FIF,)“*5](1 - FIFJ (a/F,

+ FIF,)

(6)

This expression can also be used to predict the time rate of force development dFldt and the largest value of the rate, (dF/dt), obtained under the same conditions as were present for the measurement of F, and Vm. Separating variables, integrating, and noting that F = 0 when t = 0 yield the isometric response WF,)(tlTJ = 0.5(1 + alF,)[ln TABLE

((1 + q)/(l - q)}] - q (7)

Isovelocity Shortening Contractions

During an isovelocity shortening contraction the velocity of the insertion, dY/dt, is constant and equal to the imposed velocity, Vi. If a time derivative of Eq. 2 is taken, then upon rearrangement and substitution of Eq. 1, the differential equation governing an isovelocity response is obtained (T,lF,)(dF/dt)

= 2[(FIF,)“*5][(a/Fm)(1 (a/F,

-t FIF,)

- FIFJ - VilV,]

(8)

Upon separation of variables, integration and satisfaction of the initial condition that F = 0 when t = 0, the isovelocity response is obtained

1. Morphological and physiological characteristics of six mouse and three rat muscles L

mm

Mouse muscles (25°C) EDL-1 EDL-2 EDL-3 Sol-l Sol-2 Sol-3 Rat muscles (35°C) EDL-1 EDL-2 EDL-3

5.5 6.3 5.6 7.7 8.4

8.0 12.0

12.8 11.6

Area, mm2

Mass, mg

2.16 1.81 2.12

12.7 12.1 12.6

1.15 0.69 0.67

6.1 5.7

6.85 9.43 15.91

9.4

87 128 196

Stim

Freq, HZ

FOY N/cm2

F mK

V LfZ

V mZ3

a/F,

100 100 150

22.1 26.5 22.2

450 434 456

10.0 8.6 9.6

55.3 54.0 54.0

0.332 0.351 0.362

80 120 100

21.8 28.0 28.7

245 191 173

3.6 4.3 4.9

27.4 36.0 39.5

0.188 0.197 0.144

200 120 120

25.1 25.5 21.2

1608 2197 3029

10.9 12.8 14.2

131 164 164

0.305 0.278 0.297

Lf, fiber length; area, total fiber cross-sectional area; mass, whole muscle wet mass; F,, force at 90% F, (see text); F,, maximum specific force; Vm, velocity extrapolated to 0 load for contractions at F,; a/F,, curvature of force-velocity relation.

Downloaded from www.physiology.org/journal/jappl by ${individualUser.givenNames} ${individualUser.surname} (163.015.154.053) on October 4, 2018. Copyright © 1991 American Physiological Society. All rights reserved.

CONTRACTION

TABLE

BIOMECHANICS

805

MODEL

2. Data regarding characteristic times T, and elastic element stretches D, at F, Values of

To Needed to Match Observations,

Quadr model

ms

Cubic model

IM

IT

IM

IT

Value Used*

Corresp Q,,t % Lf

EDL-1 EDL-2 EDL-3

15.7 11.8 10.9

15.9 13.0 9.5

14.2 10.5 9.6

14.3 18.0 8.1

15.0 11.2 10.3

15.0 9.6 9.9

Sol-l Sol-2 Sol-3

40.4 26.1 33.2

44.7 22.5 27.2

37.1 23.9 31.1

40.6 20.0 24.6

38.8 25.0 32.1

13.8 10.7 15.8

14.7

16.7

13.1 16.0 21.6

14.7

13.9 16.8 22.7

15.1 21.6 32.1

Mouse muscle

Rat muscles EDL-1 EDL-2 EDL-3

17.6 23.7

Quadr, quadratic; IM, value of To needed to match model and experimental time at which a force of 0.3 F, is reached in an isometric tetanus; IT, value of T,, needed to match model and experimental time at which a force of 0.9 times afterload is reached in an afterloaded isotonic test for an afterload of ~0.3 F,. * Mean of quadratic and cubic model isometric test-matched To was used for all subsequent model predictions. t Stretch of elastic element under tension F,, as percentage of fiber length.

t/T, = -q/B

+ (1/2C)(AlB

+ aIF,) x ln [(A + w(A

- WI

(9)

where q2 = (FIF,), A = (alF,)(l - Vi/V,), B = alF, + VilVm, and C2 = (AB). This equation is conveniently solved by assuming a value of contraction force F and solving Eq. 9 for the corresponding time t. Note that during an isovelocity shortening contraction the free end of the CE starts out moving at speed V,, after which its speed decreases. If enough time is allowed, the CE velocity will decay to the imposed velocity, Vi* Since these tests are terminated before excessive length changes occur, the CE velocity may or may not have decayed to Vi before test termination. Isovelocity shortening contractions are often used to measure power development capabilities of muscles. A prescribed velocity of shortening is imposed at time 0, and the test is terminated after a prescribed length change has occurred. The contraction force developed is zero at first and rises with time until the termination of the test. Common practice is to observe the force-time relationship and calculate a time-averaged force value. Power output is computed as the product of this mean force with the imposed velocity. Since Eq. 9 predicts the force-time relation, the model can be used to predict what power output can be expected when power output is calculated in this way. Model Parameters and Model Validity

The three types of contractions described above involve five parameters. Three of these, F,, Vm, and alF,, describe the force-velocity relationship of the contractile element. The other two are T,, the characteristic time for a response, and D,, which is the stretch in the spring element when F = F,. D, and T, are not independent, since To was defined so that D, = V,T,. Thus four of the five parameters are independent. Once a set of four independent parameters is assigned to the model, each of the biomechanical responses at issue can be predicted. The validity of the model can be tested by various

methods. Each of the methods involves first obtaining the four parameters from a set of experimental data under consistent conditions of muscle temperature, muscle activation, and measurement and then predicting what responses should occur in other experiments. Finally, the predicted responses are compared with the measured responses. We chose to measure F, directly during an isometric tetanic contraction and Vm and a/F, during afterloaded isotonic contractions. We then selected a value of T, based on data obtained during the isometric contraction. The values of T, needed to match quadratic and cubic model-predicted and experimental times at which a force of 0.3 F, was reached were calculated. The mean of these two isometric test-matched T, values was used for all subsequent model predictions. Subsequently, we used the model to predict the following biomechanical behaviors: 1) during isometric tetanic contractions, the time course of the contraction force increase, its maximum rate, (dFldt),, and the time at which that rate occurs, the time course of stiffness increase, and D,, the stretch of the elastic element at F,; and 2) during isovelocity shortening contractions, the time course of the contraction force for various shortening velocities and the power output for various shortening velocities and percentages of the change in total fiber length. On all but two of these response features, we collected experimental data for the comparison of predicted and measured responses. We did not gather data on small perturbation stiffnesses or values of D,. Stiffness was not measured because an apparatus to measure stiffness was not available to us. Instead, the model predictions of stiffness were compared qualitatively with data published by others on whole EDL and Sol muscles of mice and rats. Data regarding values of D, were not collected because, as will be detailed in DISCUSSION, no corresponding easily measurable estimate of D, exists. Experimental Methods

Experimental data were collected on young (4-mo-old) female CD-l mice or albino rats. Animals were anesthe-

Downloaded from www.physiology.org/journal/jappl by ${individualUser.givenNames} ${individualUser.surname} (163.015.154.053) on October 4, 2018. Copyright © 1991 American Physiological Society. All rights reserved.

806

CONTRACTION

A ul cn

125

E

Quad

BIOMECHANICS

a/Fm=0.30

Model

r

‘0 Ll sL

025

9 L

0.00 - 0

1

2

3

5

4

Normalized

B

Cubic

7

6

Time a/Fm=0.30

Model

MODEL

by an electric field generated between two large platinum electrodes situated in the bath on either side of the muscle (1). Operation for the in situ muscle preparation. The data on the EDL muscles of rats were collected in situ because of the large mass of those muscles (450 mg) and the demanding nature of the protocol (18). For in situ measurements, the EDL muscle of the rat was isolated, a suture (5-O) was tied around the distal tendon, and the distal tendon was cut. The knee and foot were immobilized on a baseplate and the distal tendon was tied by the suture to the lever arm of the servomotor (Cambridge Technology model 305). Needle electrodes were inserted on either side of the peroneal nerve. The muscle was kept moist with buffered mammalian Ringer solution at a temperature of 35*C. The temperature of the rat was maintained at 35*C by a heated platform (2). Measurement of whole muscle contractile properties. To stimulate the muscle, pulses of 0.2-ms duration at supramaximal voltage were used. The servomotor displaced the muscle and measured the displacement and force. The force response was displayed on a storage oscilloscope and sampled by microcomputer. Muscle length was adjusted to the length (I+,) at which isometric twitch and tetanic force were maximum (1). 500 -

Q) 0 b L

A

.

0.00 0

1

1

1

2

L 3

Normalized

1

I

L

I

4

5

6

7

Time

FIG. 2. Universal curves showing model-predicted responses in isometric tetanus. A: quadratic spring model; B: cubic spring model. For both A and B, force-time response is given by solid curve and stiffnesstime curve by dashed line.

1 0

tized with pentobarbital sodium (80 mg/kg) with supplemental doses given as required to maintain an adequate level of anesthesia. All experiments were carried out in accordance with the policy statement published by the American Physiological Society outlining the care and use of animals. Operation for the in vitro muscle preparation. The small mass (8-10 mg) of the EDL and Sol muscles of mice permitted the measurements of contractile properties to be made in vitro (18). For in vitro measurements, the EDL or Sol muscle of the mouse was isolated and the nerves to the muscle were dissected and cut. Ties were placed around the proximal and distal tendons of the muscle and the muscle was removed. The muscle was placed in buffered mammalian Ringer solution. After removal of muscles, animals were killed with an overdose of the anesthetic. The solution was maintained at 25*C and gassed continuously with a mixture of 95% O,-5% CO,. One tendon was tied to the lever of a servometer (Cambridge Technology model 300H) and the other tendon to a fixed point. The muscle was stimulated directly

30

60 Time

90

120

150

120

160

200

(ms)

200

f loo ii LL

0 0

40

80 Time

(ms)

3. Force-time responses of mouse EDL-1 muscle (B) in isometric tetanus. Open symbols measured values. Predicted responses are shown model by solid curves and for cubic spring model FIG.

muscle (A) and Sol-2 are experimentally for quadratic spring by dashed curves.

Downloaded from www.physiology.org/journal/jappl by ${individualUser.givenNames} ${individualUser.surname} (163.015.154.053) on October 4, 2018. Copyright © 1991 American Physiological Society. All rights reserved.

CONTRACTION

BIOMECHANICS

3. Measured and predicted values of (dFldt),

TABLE

Values of (dF/dt),,

Cubic model

13 15 17

14 18 21

6

2.5 2.8 1.6

Rat muscle EDL-1 EDL-2 EDL-3 rate

9 10

2.3 2.8 1.8

53 44 54 maximum

Times of Occurrence,

Quadr model

Sol-l Sol-2 Sol-3

WW,,

mN/ms

Measd Mouse muscles EDL-1 EDL-2 EDL-3

1.0 1.2 0.7

51 56 58

of increase

in contraction

24 25 26

force

over

4. Predicted stiffnesses at F, Stiffness,

Mouse muscles EDL-1 EDL-2 EDL-3 Sol-l Sol-2 Sol-3

largest

value

Quadr Model

Cubic Model

1,079 1,451 1,649

1,628 2,177 2,474

461 422 273

696 637

Measd

Quadr model

Cubic model

12.5 7.5 7.5

7.1 5.5 5.0

19.4 14.5 13.2

27.5 17.5 25.0

16.1 10.3 13.0

36.4 51.0

8.0 5.0 20.0

6.6 7.8 10.9

18.4 23.0 30.4

1,765 1,589 1,628 for contraction

force.

59.0

the force-velocity curve data and extrapolating to the velocity at zero load. Zsovelocity shortening contractions. The development of force and power were measured during single isovelocity shortening contractions. The L, was measured and fiber length (I,,) was determined by multiplying the L, by a previously determined L,IL, of 0.45 for mice and 0.40 for rats (2). The shortening, through 5,10, and 20% of L, occurred from 2.5,5, or 10% above L, to the same percentage below L,. Velocity of shortening was controlled by the servomotor and held constant throughout the contraction. The integrated area under the force50

40

75 g c1 5 3;

30

.-6 ij t 5

20

?i i= 10

0 0.00

I

I

0.20

1

0.40

0.60

409 Normalized

Rat muscle EDL-1 EDL-2 EDL-3 F,,

mN /mm

ms

time.

Isometric contractions. Five minutes after muscle length was set at L,, maximum twitch force was recorded. The frequency-force curve was determined by stimulating each muscle for 200 ms at a frequency of stimulation increasing in increments of 50 Hz until the force development plateaued at a maximum value (F,). One minute was taken between each stimulation train. The stimulation frequency that produces 90% of F, has been identified as the optimum frequency of stimulation for sustained power (2). The maximum force measured at that frequency of stimulation was designated F,, and the maximum rate for that development was designated (dFl dt),. All subsequent contractions were stimulated at this frequency. Afterloaded isotonic contractions. The velocity of shortening at zero load (V,) and the curvature of the Hill hyperbola (a/F,) were determined as described by Brooks and Faulkner (1). When recording velocities of shortening, the servomotor was loaded with forces equivalent to a range of loads between 0.05 and 0.50 of F,. Shortening velocities were measured during afterloaded isotonic contractions at 13 different loads. All velocities were calculated by microcomputer from the sampled displacement signals. For each preparation the values of Vm and alF, were estimated by fitting the Hill equation to TABLE

807

MODEL

2,648 2,393 2,441

Force,

F/Fm

4. Times at which mouse EDL-1 muscle (bottom) and Sol-2 muscle (top) begin shortening in afterloaded isotonic tests. Open symbols are experimentally measured values. Predicted responses are shown for quadratic spring model by solid curves and for cubic spring model by dashed curves. FIG.

Downloaded from www.physiology.org/journal/jappl by ${individualUser.givenNames} ${individualUser.surname} (163.015.154.053) on October 4, 2018. Copyright © 1991 American Physiological Society. All rights reserved.

CONTRACTION

BIOMECHANICS

MODEL

by division of the absolute value for each parameter the L, cross-sectional area, and mass, respectively.

by

RESULTS

0

10

20

30

40

Time

50

60

70

(md

80

60

8 b IL

40

20

The values of F,, VI, and alF, obtained from the measurements for each of the nine muscles are presented in Table 1. The values of T, and II, that matched the model-predicted and experimentally measured behaviors for both the quadratic and cubic spring models in both isometric and afterloaded isotonic tests showed some variation (Table 2). Since the measured force values usually lay between the values predicted by the quadratic and cubic spring models, the value of T, selected for subsequent calculations was the mean of the quadratic and cubic spring model values needed to match observed isometric responses at a force level of 0.3 F,. Values of D, corresponding to the selected values of T, ranged from 9.6 to 21.6% of L, (Table 2). When plotted in terms of the normalized ratios FIF, against t/T,, the isometric response of any muscle modeled will depend only on the parameter alF,. The predicted isometric response in this universal form for one value of a/Fm is shown for both the quadratic (Fig. 2A) and cubic (Fig. 2B) spring models. For this value of a/F,, normalized stiffness (the ratio of stiffness to the stiffness at F = F,) is also plotted against normalized time for each of the two models (Fig. 2A and B).

0 0

40

20 Time

60

80

(ms)

5. Force-time responses of mouse EDL-1 muscle (A) and Sol-2 muscle (B) in isovelocity shortening tests. Open symbols are experimentally measured values. Predicted responses are shown for quadratic spring model by solid curves and for cubic spring model by dashed curves. Responses for velocities of 5, 12, and 20 mm/s are shown for EDL-1 and for velocities of 6,9, and 12 mm/s for Sol-2. FIG.

time curve was used to determine the average force during the contraction. Power was calculated as the product of the velocity and the average force. Additional details. The total experimental protocol necessary to provide the input data for the mathematical model and to evaluate its subsequent predictions was unusually demanding for in vitro muscles, particularly for the large fast-fibered EDL muscles of rats (18). Muscles that did not maintain 85% of F, throughout the protocol were excluded from the data analysis. A total of one EDL muscle from a mouse and five EDL muscles of rats failed to meet this criterion. Upon completion of the experiments, the muscles were isolated, removed from the animal, and pinned at I,, in mammalian Ringer solution. The animal was then killed with an overdose of the anesthetic. Muscle length and I,, were measured under a dissecting microscope. The proximal and distal tendons were subsequently trimmed from the muscle and the muscle was weighed. The cross-sectional area was determined by dividing the wet mass of the muscle (mg) by the product of L, (mm) and a muscle density of 1.06 mg/mm3 (1). The Vm, F,, and maximum power (P,) were also calculated in terms of &per second, newtons per centimeter squared, and watts per kilogram

0.00

0.25

0.50

Velocity Ratio, Vi/Vm 6. Power outputs of mouse EDL-1 muscle (bottom) and Sol-2 muscle (top) in isovelocity shortening tests. Open symbols are experimentally measured values. Model-predicted responses are shown for quadratic spring model by solid curves and for cubic spring model by dashed curves. Shortening was to 10% of fiber length. FIG.

Downloaded from www.physiology.org/journal/jappl by ${individualUser.givenNames} ${individualUser.surname} (163.015.154.053) on October 4, 2018. Copyright © 1991 American Physiological Society. All rights reserved.

CONTRACTION TABLE

5.

hnpk

power

in

OUtpUtS

BIOMECHANICS

809

MODEL

tests

isOvelOCity

Measured

Shortening Velocity Ratio (Vi/Vm)

Power

mW

W/kg

Predicted Power (Cubic Model)

1.50 1.67 2.26

118 138 179

1.53 1.83 2.07

Measured-toPredicted Ratio

Mouse muscles EDL-1 EDL-2 EDL-3

0.18

Sol-l Sol-2 Sol-3

0.11 0.17 0.15

0.20 0.21

Rat muscles EDL-1

0.265 0.302

0.16 0.12 0.12

EDL-2 EDL-3

28 60 53

0.369

12.0 14.2 18.0

138 111 92

0.215 0.374 0.278 12.6 16.6 17.8

0.98 0.91 1.09 1.23

0.99 1.09

0.95 0.86 1.01

Shortenings were of 10% of fiber length except for mouse Sol-l, for which it was 5%.

Comparison of Predictions and Measurements Isometric contractions. For mouse EDL and Sol muscles, measured and predicted force responses plotted against time during isometric tetanic contractions were in good agreement (Fig. 3, A and B, shows examples). Both the quadratic and the cubic spring models closely predict the measured response, and the measured response lies between the responses predicted by the two models. This was generally true for all nine muscles tested. Comparison of measured and predicted values of (dFl dt), and the times during a contraction at which the larg-

3.00

2.40

0.60

I

0.00

0.00

1

I

0.30

0.20

0.10 Velocity

I

Ratio,

0.40

1

0.50

VP/m

7. Power outputs of mouse EDL-3 muscle in isovelocity shortening tests. Open symbols are experimentally measured values. Predicted responses are shown for quadratic spring model by solid curves and for cubic spring model by dashed curves. Shortenings were to 20% (top), 10% (mi&Ue), and 5% (bottom) of fiber length. FIG.

est values occur (Table 3) show that the measured values usually lay between the quadratic and the cubic spring model predictions. The quadratic spring model predictions generally were closer to the measured values than were the cubic spring predictions. The observed times for the largest rate of force development usually lay between the times predicted by the two models. Stiffness-time curves were predicted by the model. In terms of relative stiffness, the universal response curves for an a/F* of 0.3 show (Fig. 2, A and B) that, at an FIF, of 0.5, the predicted relative stiffness-to-relative force ratio is 0.71 and 0.79 for the quadratic and cubic spring models, respectively. No independent measurements of stiffness were made. Instead, the data reported by Stein and Gordon (19) were used for comparisons. Stein and Gordon measured stiffness during continuous length oscillations of 0.1% of L, with the frequency of oscillations ranging from 700 to 1,100 Hz, which were presumed to produce the most valid estimates of stiffness. The predicted ratios of relative stiffness to relative force are in reasonable agreement with their measured values of 0.66-0.82 for whole EDL and Sol muscles of mice and rats (19). It might be noted that relative stiffness depends only on F/F, and is not dependent of alF,. In contrast, the time that an FIF, of 0.50 is reached is dependent on alF,, with larger values of a/F0 resulting in a shorter time. With regard to absolute stiffness, Table 4 reports predicted muscle stiffness at F, for each of the nine muscles. For the three mouse Sol muscles at 25OC and at F,, our mean value of 581 mN/mm predicted from the cubic spring model is somewhat less than reported values of stiffness of 700-850 mN/mm measured at F, (19). The agreement was much poorer for the 385-mN/mm value predicted by the quadratic spring model. Since stiffness is a function of force, a correction based on the ratio FJF, of 1.11 gives predicted values of 645 and 427 mN/ mm for the cubic and quadratic spring models, respectively. After-loaded isotonic contractions. Model-predicted values of TB,the time at which shortening commences in an afterloaded isotonic contraction, seldom differed from experimentally measured values by more than 5 ms (Fig. 4). Isovelocity shortening contractions. Measured and model-nredicted forces nlotted against time for tvpical

Downloaded from www.physiology.org/journal/jappl by ${individualUser.givenNames} ${individualUser.surname} (163.015.154.053) on October 4, 2018. Copyright © 1991 American Physiological Society. All rights reserved.

810

CONTRACTION

BIOMECHANICS

MODEL

isovelocity shortening contractions demon strate reasonable agreement (Fig. 5). The cubic model was generally more accurate than the quadratic model in predicting these responses, particularly in their early phases. Power outputs measured during isovelocity shortening contractions and those predicted with both the quadratic and cubic spring models were also generally in good agreement (Fig. 6). For low velocity ratios (Vi/V,), the two models predict essentially the same power output. For larger velocity ratios, the cubic spring model predicts smaller power outputs than does the quadratic spring model. Power outputs measured and predicted by the cubic spring model for each of the nine muscles are presented in Table 5. In each case, results are given for the velocity ratio (Vi/V,) which yielded the largest values for measured power outputs . For seven of the muscles the measured power outputs were within 10% of those predicted. In one muscle the measured power was 14% lower, and in another, 23% higher than that predicted. For one mouse EDL muscle, measurements and predictions of power output for isovelocity shortenings through 5, 10 and 20% of L, were made (Fig. 7). Both models predict the qualitative character of the measurements well, particularly in showing the effects of changing the total shortening allowed. Quantitatively, these data show one of the poorest of the agreements between predicted and measured values.

force development in an isometric response, both the linear and square-root spring models predict that the maximum value of dFldT during an isometric contraction will occur at the onset of the contraction. In reality, the maximum value of dWdT occurs later than this. To avoid complexity, the models incorporate some known oversimplications. The models assume that muscles obey the Hill force-velocity relation. In the experiments used to select model parameter values, the Hill relation was determined from velocity measurements at afterloads between 0.05 and 0.50 F,. Within this portion of the force-velocity relationship, the assumption is reasonably accurate. In contrast, Claflin and Faulkner (3) demonstrated that for forces 0.70 F,. Since maximum power is developed between 0.25 and 0.30 of F, , the assumption that a Hill hyperbola governs the entire force-velocity rel .ation i s a reasonable approximation of reality. Quantitative discrepancies between predicted and observed behaviors may reflect either shortcomings in the models or variability in measured behaviors. Variability in measured behaviors may arise from nonideal response of the displacement and force transducers, from physiological variability in the behavior of the muscle, or both. Until the afterload is reached and shortening starts, the response during an afterloaded isot Ionic contraction should be identica .lto that which occurs during an isometric tetanic contraction. Despite this, slightly different valDISCUSSION ues of To are needed to match predicted and observed The findings presented show that bo th the quadratic behaviors in the two different tests (Table 2). The differand cu bit spring models qualti tatively predic t well the ences in observed behaviors duri .ng the two contractions actua 1 muscle contraction responses of the kin .ds examaPPear to result from less than ideal force and motion predi ctions of significant response control in the test apparatus. Conversely, measured beined. Quantitatively, parameters were , with so me exceptions, within 15% of haviors of a muscle may vary during a test as a result of those measured. In many cases measured response v al- fatigue (5,16). Although we can correct for the effects of by monitoring force and power, faues lay between those predicted by the quadratic and fatigue intermittently cubic spring models. This suggests that a model with an tigue may affect behavior during a single contraction. To EE in which force is proportional to stretch raised to a expect model predictions made using a single set of papower between 2.0 and 3.0 might more accurately reflect rameters to match experimental observations exactly is observed behavior. Such an improvement did not seem unreasable when simplified representations of reality are warranted. More complicated formulas would be re- incorporated into the model. quired, probably with solutions by numerical integration The parameter D, was defined as the stretch in the EE at a contraction force of F,. For the muscles studied, D, rather than in closed form. A concomitant loss in insight into how parameter changes affect behavior might result. was found typically to be lo-20% of Lf. Griffith (7) meaWe therefore selected the simpler models with their sured the shortening of fibers in the medial head of gastrocnemius muscles of cats in vivo. Piezoelectric crystals closed-form expressions, despite their occasional quantitative shortcomings in representing observed behaviors. were attached to the muscle surface at the alternate ends of fibers, and shortening was measured by an ultrasound Because the two models often seem to bracket observed behaviors, both models are presented. One does not seem transit-time technique. Shortening was estimated to be ~25% of L, during maximum isometric tetanic contracto be clearly superior to the other in predicting all behaviors studied. tions. For single intact fibers the amount of sarcomere In contrast, the linear and square-root spring models shortening is - 10% (6), although intersarcomere dynamics are of considerable significance (12). Our predicwere clearly unrealistic in predicting some contraction response features. For example, with respect to stiffness, tions for shortening lay between the estimates for a the linear spring model predicts that stiffness will be in- whole in situ muscle and for a single fiber. dependent of contraction force magnitude, whereas the The comparison of data on our whole muscles studied square-root spring model predicts that stiffness will de- in vitro and in situ with data on in vitro single intact crease as contraction force magnitude increases. In realfibers and in vivo whole muscles is complicated by differity, stiffness is observed toi ncrease with increases in co n- ences in the EE component and consequently on the pretraction force magnitu .de. With respect to the rate of dictions of D,. To the extent that tendinous tissue is presDownloaded from www.physiology.org/journal/jappl by ${individualUser.givenNames} ${individualUser.surname} (163.015.154.053) on October 4, 2018. Copyright © 1991 American Physiological Society. All rights reserved.

CONTRACTION

BIOMECHANICS

ent in series with the muscle fibers, sarcomeres will shorten to compensate for the stretch in the tendinous tissue, since total length cannot change in an isometric test. For a hypothetical muscle completely devoid oftendinous tissue, no sarcomere shortening would occur during an isometric contraction. Despite that, D, would still have a physical interpretation. In the hypothetical muscle, the actin and myosin filaments and the cross-bridges would be stretched from their relaxed values by the contraction force F,. At the same time the amount of overlap of the actin and myosin filaments would increase because of displacement as cross-bridges move through their driving stroke. In the absence of tendinous tissue, these two effects would have to cancel one another in order that sarcomere length be held constant. D, would then be a measure of the filament additional overlap or, identically, of the additional filament stretching. In real muscles, then, D, probably represents a combination of tendinous tissue stretching and filament stretching with increased overlapping. Our conclusion is that a simple Hill element-nonlinear spring model provides a reasonably valid prediction of the biomechanical behaviors of whole mammalian muscles contracting in vitro. Most significant parameters associated with isometric, afterloaded isotonic, and isovelocity shortening contractions are predicted to within 15% of measured values. This simple model should be of practical value to investigators in muscle physiology in the planning and interpretation of experiments. Appropriate use of the model should permit a reduction in the number of actual experiments performed on animals and provide additional insights into how changes in individual properties affect overall mechanical behavior in a variety of situations. APPENDIX

Equations

Governing

Responses of Cubic Spring

Model

In the modelthat usesa cubic spring EE, the force-displacement relation for the EE is FIF,

= [W - WD,13

(Al)

From the same procedures that were used in deriving the quadratic spring EE relations, those for the cubic spring model can be obtained as follows. The tangent stiffness is dFldX

= 3(F,lD,)(FlFm)‘2/3’

(AZ)

The expressionfor the CE velocity in an isometric response 1s WV,

= (T,/3F,.J

(dFldt) (FIFm)(-2’3)

(A3)

The differential equation governing isometric responseis dFldt = (3FmIT,)(alFm)(F/Fm)‘2/3~(1

- FIF,)I(aIF,

+ F/F,)

(fw The solution to the above equation, accounting for the appropriate initial condition, provides the isometric responseof the model t/T,

= -r(alF,)

+ (1 + a/F,){

- R arctan [(-2r

-P/2R

- 1)/R] - 0.5 In [(? - 2r + l)/ (3 + r + 1)]}/(3alF,)

(A5)

811

MODEL

where P = 3.1416, R2 = 3, and F = (F/&J. The differential equation governing the isovelocity response

is (T,lF,)dFldt

= 3(FlFm)(2’3’((alFm)[(1 (a/F,

- F/F,)/ + FIF,)]

- Vi/V,

(A6)

The isovelocity responseis then t/T,

= -r/B X

+ (a/F, { (P/2R)

+ AIB)(K/3A) + 0.5 In [(K + r)2/(p

- Kr + ?)]

+ R arctan [(2r - K)lR] where r, P, and R are aspreviously defined and A = (alF,)( V”/Vm), B = a/F, + VilVm, and P = -AlB.

l-

SusanV. Brooks and James M. Kneebone assisted in the collection of the experimental data. This research was supported by National Institutes of Health Grants AG-06157, AG-06621, and NS-20536. Address for reprint requests: A. B. Schultz, Dept. of Mechanical Engineering, University of Michigan, Ann Arbor, MI 481092125. Received 1990.

26 September

1989;

accepted

in final

form

11 September

REFERENCES 1. BROOKS, S. V., AND J. A. FAULKNER. Contractile properties of skeletal muscles from young, adult and aged mice. J. Physiol. Lo&. 404: 71-82, 1988. 2. BROOKS, S. V., J. A. FAULKNER, AND D. A. MCCUBBFLEY. Power outputs of slow and fast skeletal muscles of mice. J. Appl. Physiol. 68: 1282-1285,1989. 3. CLAFLIN, D. R., AND J. A. FAULKNER. The force-velocity relationship at high shortening velocities in the soleus muscle of the rat. J. Physiol. Lond. 411: 627-637, 1989. 4. EDMAN, K. A. P. Double-hyperbolic force-velocity relation in frog muscle fibres. J. Physiol. Lond. 404: 301-321, 1988. 5. EDMAN, K. A. P., AND A. R. MATTIAZZI. Effects of fatigue and altered pH on isometric force and velocity of shortening at zero load in frog muscles. J. Muscle Res. CeZZ Motil. 2: 321-334, 1981. 6. EDMAN, K. A. P., AND C. REGGIANI. The sarcomere length-tension relation determined in short segments of intact muscle fibres of the frog. J. Physiol. Lond. 385: 709-731, 1987. 7. GRIFFITH, R. I. Ultrasound transit time gives direct measurement of muscle fibre length in vivo. J. Neurosci. Methods 21: 159-165, 1987. 8. HILL, A. V. The heat of shortening and the dynamic constants of muscle. Proc R. Sot. B BioZ. Sci. 126: 136-195, 1938. 9. HUXLEY, A. F. Muscle structure and theories of contraction. Prog. Biophys. 7: 255-318,1957. 10. HUXLEY, A. F., AND R. M. SIMMONS. Proposed mechanism of force generation in striated muscle. Nature Lond. 233: 533-538, 1971. 11. HUXLEY, A. F., AND R. M. SIMMONS. Mechanical transients and the origin of muscular force. Cold Spring Harbor Symp. Qwnt. BioZ. 37: 669-680,1973. 12. JULIAN, F. J., AND D. L. MORGAN. Intersarcomere dynamics during fixed-end tetanic contractions of frog muscle fibres. J. Physiol. Lond. 298: 365-378,1979. 13. JULIAN, F. J., K. R. SOLLINS, AND M. R. SOLLINS. A model for muscle contraction in which cross-bridge attachment and force generation are distinct. CoZd Spring Harbor Syrnp. Quant. BioZ. 37: 685688,1973. 14. MACMAHON, T. A. MuscZes, Reflexes, and Locomotion. Princeton, NJ: Princeton Univ. Press, 1984, p. 331. 15. METZGER, J. M., AND R. L. MOSS. Greater hydrogen ion-induced depression of tension and velocity in skinned single fibres of rat fast and slow muscles. J. Physiol. Lond. 393: 727-742, 1987. 16. MORGAN, D. L., S. MOCHON, AND F. J. JULIAN. A quantitative model of intersarcomere dynamics during fixed-end contractions of single frog muscle fibers. Biophys. J. 39: 189-196, 1982.

Downloaded from www.physiology.org/journal/jappl by ${individualUser.givenNames} ${individualUser.surname} (163.015.154.053) on October 4, 2018. Copyright © 1991 American Physiological Society. All rights reserved.

812

CONTRACTION

BIOMECHANICS

S. J. AND D. A. WINTER. Predictions of knee and ankle moments of force in walking from EMG and kinematic data. J.

MODEL

R. B., AND T. GORDON. Nonlinear stiffness-force relationships in whole mammalian skeletal muscles. Can. J. Physiol. Phur-

17. OLNEY,

19. STEIN,

Biomech. 18: 9-20, 1985. 18. SEGAL, S. S., AND J. A. FAULKNER.

nacol. 64: 1236-1244, 1986. 20. ZAJAC, F. E., R. W. WICKE,

Temperature-dependent physiological stability of rat skeletal muscle in vitro. Am. J. Physiol. 248 (Cell Physiol.

17): C265-C270,

1985.

AND W. S. LEVINE. Dependence of jumping performance on muscle properties when humans use only calf muscles for propulsion. J. Biomech. 17: 513-523, 1984.

Downloaded from www.physiology.org/journal/jappl by ${individualUser.givenNames} ${individualUser.surname} (163.015.154.053) on October 4, 2018. Copyright © 1991 American Physiological Society. All rights reserved.

A simple Hill element-nonlinear spring model of muscle contraction biomechanics.

The purpose of this study was to develop a model to predict the mechanical response of muscles during isometric tetanic, afterloaded isotonic and isov...
2MB Sizes 0 Downloads 0 Views