| HOME | HELP | FEEDBACK | SUBSCRIPTIONS | ARCHIVE | SEARCH | TABLE OF CONTENTS |
|
|
|
|||||||
|
|||||||||
1Department of Veterinary and Comparative Anatomy, Pharmacology, and Physiology, Washington State University, Pullman, Washington 99163; and 2Department of Biomedical Engineering, New Jersey Institute of Technology, Newark, New Jersey 07102
Submitted 5 November 2003 ; accepted in final form 11 December 2003
 |
ABSTRACT |
|---|
 TOP ABSTRACT METHODSRESULTSDISCUSSIONREFERENCES |
|---|
-myosin heavy chain (MHC) (mouse, rat) or 
-MHC (rabbit, ferret) were activated with solutions from pCa 6.1 to 4.3. Activated fibers were subjected to small-amplitude length perturbations [
L(t)] rich in frequency content between 0.1 and 40 Hz. In descriptive validation tests, the model was fit to the ensuing force response [F(t)] in the time domain. In fits to 118 records, the model successfully accounted for most of the measured variation in F(t) (R2 range, 0.997–0.736; median, 0.981). When some residual variations in F(t) were not accounted for by the model (as at low activation), there was very little coherence (<0.5) between these residual force variations and the applied L(t) input function, indicating that something other than L(t) was causing the measured variation in F(t). With one exception, model parameters were estimated with standard errors on the order of 1% or less. Thus parameters of the recruitment component of the model could be uniquely separated from parameters of the distortion component of the model and parameters estimated from any given fiber could be considered unique to that fiber. In explanative validation tests, we found that recruitment and distortion parameters were positively correlated with independent assessments of the physiological entity they were assumed to represent. The recruitment distortion model was judged to be valid from both descriptive and explanative perspectives and is, therefore, a useful construct for describing and explaining dynamic force-length relationships in constantly activated cardiac muscle.
muscle stiffness; cross-bridge recruitment; cross-bridge distortion; model validation; mouse; rat; ferret; rabbit
Dynamic stiffness has played a variety of roles in muscle physiology. In insect flight muscle, dynamic stiffness has been used extensively as a concept to explain and predict function of asynchronous flight muscle resulting from interaction of the muscle with its mechanical load (18, 22). In vertebrate muscle, Kawai and colleagues (14–16, 28, 33) pioneered the use of complex stiffness and its dynamic decomposition to characterize the muscle fiber contractile state, to identify specific steps in the cross-bridge (XB) cycle, and to characterize the sensitivity of these steps to ATP and its metabolic products. They extended these methods to cardiac muscle (17, 27, 34), and their pioneering work has now been applied and extended by several other groups (2, 3, 6, 13, 19, 21, 26, 29, 32).
One outcome of the work of Kawai et al. (14, 15) has been the routine practice of decomposing the measured dynamic stiffness into a sequence of first-order dynamic processes, each with increasing characteristic frequency (3, 17, 19, 21, 27, 32, 34). Thus dynamic stiffness is represented as the sum of component stiffnesses due to processes A, B, C, and so on, where process A, with characteristic frequency a, is slower than process B, with characteristic frequency b, which is slower than process C, with characteristic frequency c, and so on.
In the present study, we take a fresh look at the origins of dynamic processes that give rise to the two dominant dynamic modes normally detected in the dynamic stiffness of constantly activated cardiac muscle (equivalent to processes B and C in Kawai et al.'s treatment). We give unique interpretations for these processes based on a functional view of sarcomere length (SL)-induced XB recruitment and stretch-induced XB distortion. Furthermore, we introduce a time-domain methodology for eliciting and analyzing dynamic behavior of constantly activated cardiac muscle that has several advantages over the previously used steady-state sinusoidal method.
|
METHODS |
|---|
|
TOPABSTRACTMETHODSRESULTSDISCUSSIONREFERENCES |
|---|
Model Development
The basis for the present model was developed in earlier publications (5, 6, 23, 24) where the origin of constantly activated myocardial dynamics was considered to arise from the myofilament kinetic scheme shown in Fig. 1. This scheme included three basic kinetic processes: 1) Ca2+ binding to, and release from, troponin C; 2) thin filament regulatory protein complexes switching from non-permissive (off) to permissive (on) states; and 3) XB cycling between attached and detached states. In addition, the following kinetic-modifying mechanisms were included:
|
First, cooperative mechanisms, including 1) cooperative interactions between neighboring tropomyosin-based regulatory units (RU); 2) cooperative interactions between neighboring XBs; and 3) cooperative interactions between a strongly bound XB and its neighboring RU.
Second, length-dependent kinetic coefficients that arise from myofilament lattice structure changes when changes in SL alter the radial distance a myosin head must travel to attach to its actin binding site.
Finally, sliding filament overlap, which alters the number of available sites for XB attachment to actin.
These kinetic processes and the mechanisms that modify them were viewed relative to their net effect either on the number of parallel, attached force-generating XBs or on the average elastic distortion of these XBs (2, 6, 23). This point of view followed from the general assumption that XB-generated force is equal to the number of parallel, attached force-generating XBs times the average force in a single representative XB, which, in turn, equals the stiffness coefficient of a single XB times the average elastic stretch among all force-generating XBs. Mechanisms determining the number of attached XBs were grouped collectively into "recruitment" mechanisms, and mechanisms determining the average XB elastic stretch were grouped into "distortion" mechanisms. Kinetic features responsible for changes in the number of attached XBs (recruitment dynamics) differ substantially from the kinetic features responsible for changes in average XB elastic stretch (distortion dynamics).
When the kinetic scheme (Fig. 1) and kinetic-modifying features were used to derive a set of equations for a general dynamic force-length relation, the resulting equation set was highly nonlinear, of high dynamic order, and contained many parameters (5, 6, 24). To give such complicated equations some practical utility, we applied linearization and model-order reduction techniques to form a simplified recruitment-distortion model. This simplified linear model is valid for small perturbations (including oscillations around a steady static state). The simplified models presented earlier (6) are simplified even further here as a retrospective outcome of analysis of data presented in this report.
The time-domain equations of the simplified model, including the differential equations for recruitment and distortion, are as follows
|
| (1) |
|
| (2) |
|
| (3) |
(t) is the model-predicted variation in force in response to a measured variation in muscle length, L(t) (through recruitment dynamics, Eq. 2), and in response to the first time derivative of muscle length, 
(t) (through distortion dynamics, Eq. 3). 
(t) is the recruitment variable; it describes the incremental addition of XBs acting in parallel to produce force. x(t) is the distortion variable; it describes the average distortion of internal stretch with the elastic regions of XBs. A dot over a variable represents the first time derivative of that variable. Parameter E0 is the slope of the static FLR. Parameter E
is instantaneous stiffness, as estimated from the initial force response to a sudden stretch, i.e., the slope of the line representing change in initial tension vs. change in muscle length in quick stretch experiments. Parameter b is the rate constant governing recruitment dynamics. b is a complicated function of all the kinetic rate constants in Fig. 1 [the on and off rate constants (kon and koff), forward and backward attachment rate constants (f and f'), forward and backward power stroke rate constants (h and h'), and the XB detachment rate constant (g)] because they are determined by the effect of the nonlinear kinetic-modifying mechanisms evaluated at the static reference state around which L(t) is imposed (6). Because the analysis is restricted to small perturbations, the nonlinearities participate by contributing to the value of constant coefficients in the reduced model. The multiple kinetic constants in Fig. 1 coalesce into a single constant, b, as a result of model-order reduction procedures. Parameter c is the rate constant governing distortion dynamics. It also has a complicated functional dependence on all the rate determining factors in Fig. 1, but it is most strongly sensitive to the XB detachment rate constant g (6).
Equations 1–3 may be used to describe both transient and steadystate behavior. When the sinusoidal steady state is of interest, Fourier transformation may be applied to Eqs. 1–3 to produce a frequency-domain expression of the complex (dynamic) muscle stiffness, 
(j
), as follows
|
| (4) |
and is angular frequency. In this complex stiffness formulation, the recruitment dynamic takes on a low-pass character (passing low frequencies without distortion while attenuating high frequencies) and the distortion dynamic takes on a high-pass character (passing high frequencies without distortion while attenuating low frequencies). From Eq. 4, it can be seen that E0 and E are zero-frequency and infinite-frequency stiffness magnitudes that scale the contributions of recruitment and distortion to overall stiffness. Furthermore, the rate constants b and c represent characteristic frequencies of the low-pass recruitment and the high-pass distortion processes, respectively. If b << c, these two processes separate such that low-pass recruitment dominates at low frequencies and high-pass distortion dominates at high frequencies. Cardiac Muscle Experiments
Experiments were conducted in skinned, constantly Ca2+-activated cardiac muscle fibers to generate data for two purposes: 1) descriptive model validation by providing dynamic force-length data to which the model was fit, and 2) explanative validation by providing independent assessment of contractile kinetic processes that were putatively represented by specific model parameters. To examine a range of cardiac muscle dynamic behaviors, skinned cardiac muscle fibers were obtained from hearts of animals expressing primarily -myosin heavy chain (MHC) [adult mouse (n = 4), adult rat (n = 4)] and primarily -MHC [adult rabbit (n = 3), adult ferret (n = 4), propylthiouracil-treated (PTU) adult rat (n = 4)]. In addition, a range of dynamic behaviors was obtained by bathing with solutions of different Ca2+ concentration to activate fibers to various levels, from low to maximum Ca2+ activation.
Skinned cardiac muscle fiber bundles. Heart tissue for cardiac muscle fibers was obtained from animals according to a protocol approved by the Washington State University Institutional Animal Care and Use Committee. All animals in this study received humane care in compliance with the animal use principles of the American Physiological Society and the Principles of Laboratory Animal Care formulated by the National Society of Medical Research and the National Institutes of Health's Guide for the Care and Use of Laboratory Animals (NIH Pub. No. 85-23, Revised 1985). Procedures for the preparation of skinned cardiac muscle fiber bundles have been described elsewhere (8). Briefly, animals were anesthetized with a mixture of ketamine (100 mg/kg) and xylazine (8 mg/kg). Hearts were rapidly removed and placed in an ice-cold relaxing solution containing 20 mM MOPS (pH 7.0), l0 mM EGTA, 5.0 mM MgATP2–, 1 mM free Mg2+, 12 mM creatine phosphate, 53 mM KCl, 1 mM DTT and a cocktail of protease inhibitors (ionic strength of 150 mM). Papillary muscles from these hearts were dissected into strips 
150–200 µm in width and 2–3 mm in length. Skinning of cardiac muscle fibers was performed at 4°C overnight in relaxing solution containing 1% Triton X-100.
Skinned muscle fiber bundles were attached to a displacement motor (model 300B, Cambridge Technology) on one end and to a force transducer (AE 801, SensoNor) on the other end using aluminum T-clips. Mounted fibers were placed in a 100-µl well, where they were bathed in either activating or relaxing solutions. SL was measured by laser light diffraction methods. Visual display of the photodiode signal allowed on-line judgment as to the quality of the preparation and the accuracy of the SL assessment. SL in all preparations was set to 2.2 µm before Ca2+ activation. Fiber cross-sectional area was calculated from optical measurements of major and minor axes using an assumption of an elliptical fiber cross section. Normalized force was expressed in milliNewtons per millimeter squared.
The Ca2+ content of the activating solution was changed, according to a computer program, to produce static force-generating activations at pCa 4.3, 5.3, 5.6, 5.8, 6.0, 6.2, and 8.0. At each level of activation and its associated static force (Fs), the computer-controlled motor was commanded to first hold muscle length constant for a period while changes in UV absorbance at 340 nm by NADH were measured to determine myofilament ATPase (8, 9) and muscle length was then commanded to change according to the following two protocols.
DYNAMIC FLR. The dynamic FLR protocol was designed to provide information rich in content at all frequencies between 0.1 and 40 Hz. To do this, fiber length (LM) was commanded to change according to a constant amplitude (0.5% of LM) sinusoid of continuously varying frequency, i.e., a chirp. Two chirps were delivered over two sequential time periods: in one period of 40-s duration, chirp frequencies varied between 0.1 and 4 Hz to emphasize low-frequency behavior, and in a second time period of 5-s duration, chirp frequencies varied between 1 and 40 Hz to emphasize higher-frequency behavior. Measured force changes [F(t)] during the FLR protocol were maximally 10% of the Fs baseline. Thus a small wander of baseline Fs during signal recording represented a significant component of the variation in the recorded F(t) signal. Accordingly, baseline trends and wander were removed from the F(t) record by fitting a fourth-order polynomial in time to the F(t) signal. It can be shown that the frequency content of the fourth-order polynomial was in a range below the frequency composition of the L(t) signal. Examples of data obtained with this protocol are shown in Fig. 2.
|
DYNAMIC FORCE REDEVELOPMENT. After the dynamic FLR protocol, and while the muscle continued to maintain the constant level of force achieved during the activation, a quick release-restretch was performed to completely unload the fiber and then reestablish initial length with few XBs in the attached state, as indicated by low force after restretch (8). The rate constant of force redevelopment (ktr) was estimated by fitting the time course of force redevelopment after the release and restretch [F(t)] with the following equation: F(t) = Fss(1 – e–ktrt) + F0, where Fss is the eventual steady-state force and F0 is the initial force from which force redevelopment began.
MHC assays. To establish the relative myocardial content of - and -MHC, left ventricular tissue from freshly killed animals was rinsed in ice-cold saline and frozen immediately in liquid nitrogen. Tissue samples were extensively pulverized while still frozen. The frozen powder was resuspended in 20 volumes of 50 mM Tris (pH 6.8), 2.5% SDS, 10% glycerol, 1 mM DTT, 2 µg/ml leupeptin, 1 µg/ml pepstatin A, 1 mM PMSF, 4 mM benzamidine, 2 µM E-64, and 5 µM bestatin. Tissue suspensions were homogenized at 4°C. This was followed by sonication at 4°C for 20 min. The samples were spun, and the protein concentration in the supernatant was determined. 2-Mercaptoethanol was added to a final concentration of 5% before samples were heated. Samples were run for 6–12 h on 6.5% acrylamide gels at 4°C. Gels were stained with a modified Coomassie blue stain (Pierce Gel Code Blue) for 1–6 h and destained with water. For Western blots, proteins were transferred to polyvinylidene difluoride membranes. The membranes were blocked overnight at 4°C with 5% nonfat milk in Tris-buffered saline-Tween 20 (TBST; 10 mM Tris and 0.05% Tween 20 at pH 7.4). Blotted membranes were incubated for 1 h at room temperature with a monoclonal antibody recognizing MHC diluted 1:5,000 in TBST-milk. Blots were then incubated for 1 h at room temperature with an appropriate secondary antibody diluted 1:5,000 in TBST-milk. Proteins were visualized using Western blue stabilized color substrate (Promega). Relative amounts of either - or -MHC were quantified using a scanner and NIH Image analysis software.
Model Validation
The purpose of data analysis was to test the model for validity. Two levels of validation were sought: descriptive validation and explanative validation.
Descriptive validation: fitting the dynamic FLR data. There were three categories of descriptive validation indexes: goodness of fit indexes, indexes of systematic trends in residual errors, and indexes of precision of model parameter estimates. These indexes were obtained from a data fitting procedure as follows. Measured F(t) from the dynamic FLR protocol contained transient behavior both at the initiation of the L(t) chirp and then throughout the chirp as frequencies changed continuously during the protocol. Thus steady state at any given frequency was never achieved. Accordingly, frequency analysis was not performed, and analysis was conducted in the time domain by fitting the model to the data. Briefly, for a given set of model parameters (E0, b, E, and c) and a given recorded L(t), model-predicted force [(t)] was obtained by numerically integrating Eqs. 2 and 3 using L(t) as the input function (fourth-order Runge Kutta; integration step size equal to the data sampling interval, 0.001 s) to obtain the recruitment and distortion variables, (t) and x(t), respectively. These were then used to calculate (t) according to Eq. 1. Because (t) was to be fit to measured F(t), (t) was low-pass filtered (cutoff frequency = 200 Hz) just as measured F(t) was filtered during recording. F(t) was fit by a procedure in which model parameters (E0, b, E, and c) were adjusted (MINPACK, Argonne National Laboratories) to minimize the mean sum of square error between (t) and F(t). The time series of residual errors, 
(t) = F(t) – (t), was calculated to determine the mean sum of square error during optimization procedures and, finally, was calculated after optimized parameters were found to assess both model goodness of fit and systematic residual errors.
CATEGORY 1 OF DESCRIPTIVE VALIDATION INDEXES. Goodness of fit indexes were taken from the R2 of the nonlinear regression of model fit to data and from the linear regression of (t) on measured F(t). R2, as a measure of the fraction of total variation in measured force accounted for by the model, was used in conjunction with the coherence between L(t) input and (t) (see below) to evaluate the importance of data variation not accounted for by the model. Linear regression parameters of slope (m) and intercept (b) for (t) = f[F(t)] were also employed as indexes of fit; ideal regression parameters were m = 1 and b = 0. Confidence intervals (95%) on the m and b estimates were used to test whether either was significantly different from its respective ideal value.
CATEGORY 2 OF DESCRIPTIVE VALIDATION INDEXES. Systematic trends in residual errors are damaging to any data fitting result. Because interest in these experiments was in dynamic behavior, damaging systematic trends in (t) are those correlated with the L(t) input forcing function. To determine whether such systematic trends occurred, the frequency-dependent coherence function 
between the L(t) input and (t) was calculated as follows
 | (5) |
L(t) and (t) and GL,L and Ge,e are the autocorrelation spectra of L(t) and (t), respectively. Values of 
range between 0 and 1. Ideally, there would be no coherence between L(t) and (t) and 
. However, we were satisfied that 
indicated that there was no meaningful systematic trends in (t) due to L(t).
CATEGORY 3 OF DESCRIPTIVE VALIDATION INDEXES. As a final assessment of the ability of the model to describe the data, we examined the precision with which model parameters E0, b, E, and c were estimated, i.e., the standard errors of the parameter estimates. If the standard errors were small relative to the estimated parameter value (<5%), we concluded that the assumed distortional and recruitment variables were largely independent in their effects on force and that these independent effects were therefore clearly discernible in the measured force response.
Explanative validation: correlation with independently assessed contractile parameters. Explanative validation was based on correlation between model-estimated parameters and independently determined parameters that putatively represented the physiological entity corresponding to the particular model parameter. Strong correlation between a model parameter and independent physiological assessment was taken as evidence that the model parameter could be used to represent the underlying physiological event and, by this, the model could be used to explain underlying phenomena responsible for observed dynamic behavior. Again, variation in parameters over which the ensuing tests were conducted was obtained within fibers by variation in Ca2+ activation and between fibers by obtaining data from fibers from animals exhibiting widely different - and -MHC composition.
FS VS. E. The parameter E in the recruitment-distortion model purportedly represents the collective stiffness of all parallel attached XBs; thus E = A2, where is the stiffness of a single XB and A2 is the number of parallel XBs in the A2 attached state. Similarly, static isometric force (Fs) around which F(t) was induced and E estimated, is also proportional to the number of parallel force generators, i.e., Fs 
A2. Thus Fs and E are both proportional to A2 and both should represent the same underlying physiological quantity. For this to be true, these two quantities must be linearly correlated. Accordingly, linear regression analysis was used to relate measured Fs to model-estimated E, and the resultant correlation coefficient was calculated.
(ATPASE/FS) VS. C. The distortion rate constant of the model (c) has a strong dependence on the XB detachment rate constant g, and thus is a putative index of g. Furthermore, it can be shown that the energy cost of force development (ATPase/Fs) is also proportional to g (8, 9). Thus, as an additional test of explanative validity, we tested whether independently measured ATPase/Fs and model-estimated c were linearly correlated.
KTR VS. B. The recruitment rate constant of the model (b) putatively represents a characteristic time by which XBs accumulate in (i.e., are recruited into) the force-bearing state. Similarly, the rate constant of force redevelopment after release and restretch (ktr) also represents a characteristic time for XB accumulation in the force-bearing state. Again, because b and ktr supposedly represent the same underlying dynamic process, we tested whether independently measured ktr and model-estimated b were linearly correlated.
|
RESULTS |
|---|
|
TOPABSTRACTMETHODSRESULTSDISCUSSIONREFERENCES |
|---|
The model was fit to 118 data records: 27 data records from 4 mouse fibers (6–7 activation levels/fiber), 26 data records from 4 rat fibers (6–7 activation levels/fiber), 22 data records from 4 fibers from PTU-treated rats (4–6 activation levels/fiber), 18 data records from 3 rabbit fibers (6 activation levels/fiber), and 25 data records from 4 ferret fibers (4–7 activation levels/fiber). Stable force during constant activation was typically 50–65 mN/mm2 at pCa 4.3 (highest was 73.9 mN/mm2 in a ferret fiber) and was typically 10–15 mN/mm2 at pCa 6.1 (lowest was 9.6 mN/mm2 in a rabbit fiber). R2 from these fits ranged from 0.997 to 0.736 with a median of 0.981. More than 80% of fits achieved an R2 > 0.97 and only 5% had an R2 < 0.85. The lowest R2s (<0.85) were all obtained in fits at low activation levels (pCa 6.0 and 6.1).
Examples of two fits to data are shown in Fig. 2, where records obtained from a control rat fiber (expressing predominantly -MHC; R2 = 0.962, which was less than the median for R2) and from a PTU-treated rat fiber (expressing predominantly -MHC; R2 = 0.994, which was greater than the median for R2) are shown. Both records were obtained from partially activated fibers at pCa 5.8, which produced force equal to 59% and 67% of maximally activated force, respectively. The full response to both low- and high-frequency chirps over 45 s, as well as an expanded segment of the last 0.2 s of the response where the frequency approached 40 Hz, are shown.
When (t) was plotted versus F(t) for the two records from the two fibers in Fig. 2, the results were as shown in Fig. 3. The outcome of regressing (t) versus F(t) was as follows: control rat – m = 0.961*, b = 0.0044; and PTU-treated rat – m = 0.995*, b = 0.0562* [*regression statistic was significantly different than the ideal value (P = 0.05)]. Despite this statistically significant difference (based on 44,995 degrees of freedom in the model fit), m closely approximated 1.0 and b closely approximated 0.0 in each of these cases and in all 118 cases analyzed. The standard errors of the model parameter estimates were <0.75% for each of the four estimated model parameters (E0, E, b, and c) in each fiber, indicating that the input V(t), acting through the recruitment part of the model containing parameters E0 and b, and the input 
(t), acting through the distortion part of the model containing parameters E and c, were independent in producing their effects on F(t).
|
Features of note in these examples, as was the case in fits to all other records, included the following: 1) the model fit through the noise in the measured F(t) signal; 2) both the phase and amplitude of measured F(t) were well reproduced by the model throughout the record from the lowest frequency part of the signal (Fig. 2, far left) to the highest frequency part of the signal (Fig. 2, far right and expanded time scales); and 3) the model reliably reproduced the region of the signal at which the response amplitude became minimal (Fig. 2, arrows), demonstrating the differences in the frequency of the minimal amplitude response between the predominantly -MHC muscle fiber (1.8 Hz) versus the predominantly -MHC muscle fiber (0.8 Hz).
Even though very little of the fractional variation in the measured F(t) (equal to 1 – R2) was not accounted for by the model in the two examples shown, and in 95% of the records analyzed, there was still a need to look for systematic error in the residuals. We did this by examining the coherence [] of (t) with respect to the input L(t). Small indicate that there was no relation between (t) and L(t), demonstrating that no improvement could be made in the model regardless of the magnitude of 1 – R2. For the two fits to the data from the two fibers shown in Fig. 2, was small (<0.5) in some frequency regions and large (>0.8) in others. In both fibers, was >0.8 for frequencies above 15 Hz. In the -MHC fiber, was also >0.8 in the frequency range of 2–7 Hz. For all -MHC fibers from the rabbit and ferret, was uniformly low (<0.5) across the full range of frequencies represented in the data. For -MHC fibers from the mouse and rat, the patterns of described for the two examples shown in Figs. 2 and 3 were typical of what was found in fits at high activation levels (pCa < 5.8). However, in these mouse and rat fibers became progressively smaller and uniformly <0.5 across the full frequency range at lower levels of activation as R2 decreased and, correspondingly, as 1 – R2 increased.
We attribute these patterns of behavior to baseline instability and imperfect trend removal rather than to model imperfections. Consider two sinusoidal signals that were perfectly matched with zero difference. Now, offset one, either up or down, from the other and take the resulting difference between these signals. The difference produces a sinusoid with the same frequency as the two original signals and with an amplitude that is proportional to the offset. Furthermore, there would be a high correlation (coherence) between this difference sinusoid and any other sinusoidal signal of the same frequency. Now, consider that as a result of imperfect trend removal there was a small baseline offset in the high-frequency region of the data records. The difference between the model prediction and the recorded data would then produce a small (t) signal that would be highly correlated with the L(t) signal. For instance, (t) for the high-frequency region of the signal shown in Fig. 2 is given as the small amplitude (dashed line). Although small in value, this (t) is very systematic and highly correlated with L(t), as indicated by the high value of noted at frequencies approaching 40 Hz. However, the periodic (t) shown in Fig. 2 has an average value of 0.19, indicating that the measured F(t) was offset from the model-predicted F(t). No correction could be made in the model that would remove this offset. Thus the problem is most likely one of imperfect trend removal and not of model deficiency.
Now consider an example at low activation (pCa 6.0) where the R2 was low (0.774) and thus a significant fraction of the measured signal was not accounted for by the model (Fig. 4). Despite the low R2, was <0.5 over the full frequency range (0–40 Hz), indicating that there was no relation between (t) and V(t). That part of the variation in the measured F(t) that was not accounted for by the model was caused by something other than L(t), i.e., intrinsic baseline instability. With such marked baseline instability, our trend removal procedures were relatively ineffective. In this case, nothing could be done to improve a model that related F(t) to L(t) because some of the measured F(t) was due to factors other than L(t). By this reasoning, we concluded that low R2 in conjunction with low was an acceptable outcome from model fitting.
|
On the basis of the ability of the model to fit the data well, the recruitment-distortion model of the cardiac muscle dynamic FLR was judged to be descriptively valid.
Explanative Validation
Fs vs. E. Explanative validation required that model parameters have physiological meaning. The model parameter E putatively represents the number of parallel force-generating XBs responsible for Fs. If so, then model-estimated E should linearly correlate with measured Fs. Accordingly, the 118 estimates of E were plotted against the corresponding 118 measurements of Fs and tested for a linear correlation (Fig. 5). The correlation coefficient between E and Fs was 0.87 (P < 0.01). Because of this significant correlation, model-estimated E could be taken to represent the same underlying physiological feature as measured Fs.
|
(ATPase/Fs) vs. c. Because the model-estimated distortion rate constant c has a strong dependence on the XB detachment rate constant g (6), and because it can be shown that the energy cost of force development (ATPase/Fs), once adjusted for fiber volume and cross-sectional area, is also proportional to g (5), we tested whether (ATPase/Fs) was linearly correlated with c to establish whether model-estimated c could be used interpretively to represent g, just as (ATPase/Fs) often is. The correlation coefficient in the relation between (ATPase/Fs) and c (Fig. 6) was 0.91 (P < 0.01). Thus we concluded that model-estimated c could be used interpretively to represent g.
|
ktr vs. b. The rate constant of force redevelopment after release and restretch ktr was used as an independent assessment of a recruitment rate constant to be compared with model-estimated recruitment rate constant b (Fig. 7). The correlation between ktr and b was 0.79 (P < 0.01). Thus model parameter b represents much of the same dynamic behavior that arises from underlying XB recruitment processes as does ktr and can be used interpretively in much the same way as ktr often is.
|
|
DISCUSSION |
|---|
|
TOPABSTRACTMETHODSRESULTSDISCUSSIONREFERENCES |
|---|
The recruitment-distortion model and the methods by which it was applied are novel. We now consider these novel model features and the insights they provide.
Novel Features
Derivation. We derived the general form of the present model, consisting of a dynamic recruitment block and a dynamic distortion block, where all kinetic events within the myofilament system (Fig. 1), including thin filament processes regulating the switching from nonpermissive (off) to permissive (on) thin filament states, were taken into account, in a previous paper (6). This global view of myofilament kinetics was chosen because length-sensing mechanisms play an important role in myofilament activation (1, 10, 11, 20, 31). Thus the model derivation took into account protein interactions in addition to just those between actin and myosin. For instance, consider a sarcomere at a given length and in a state of partial Ca2+ activation where only some fraction of the total actin sites is in the on state and available for myosin binding. Length-dependent activation requires that an increase in SL causes a greater fraction of these actin sites to switch to the on state. By including thin filament regulatory dynamics, the model presented here differs from all other models designed to represent myocardial force-length dynamics; these other models have ascribed all aspects of the FLR strictly to XB interactions between actin and myosin (2, 3, 14–17, 19, 21, 25, 27–30, 32–34). As a result, different, and we believe more insightful, interpretations of the measured dynamic force-length relation (and complex stiffness) will be obtained from using this recruitment-distortion model than would be obtained from using other models.
For instance, the dynamics of XB distortion is due to XB elastic properties coupled with XB turnover and may be appreciated from what happens when a sarcomere is suddenly stretched (Fig. 8). A sudden stretch increases the average length of the elastic link within attached XBs. The consequence is that force generated by these overstretched XBs also increases. We call the overstretch a "distortion" of the XB linkage. However, as XBs continue to cycle, those with increased distortion due to the sudden stretch detach and are replaced by newly formed nondistorted XBs. Force produced by distortion returns to its original isometric level as attached XBs become progressively populated with XBs that did not experience the sudden change in length. This dynamic of reestablishing a XB population that no longer remembers the sudden stretch determines the dynamic character of the XB distortion component.
|
The dynamics of XB recruitment may be appreciated by further consideration of the sudden stretch example (Fig. 8). When the sarcomeres undergo a sudden increase in length, the length-sensing mechanisms within the myofilament system (including filament overlap, length-dependent XB attachment, and amplification of these effects by cooperativity) initiate a sequence of events that increases the number of force-bearing XBs. Force rises as more force-bearing XBs are recruited. The time course over which this occurs depends on the overall dynamic behavior of the entire myofilament system, including thin filament on-off transitions as well as XB cycling (Fig. 1). Because such recruitment embraces the entire myofilament system, its dynamics are greatly affected by cooperativity (19). We have shown that this overall recruitment dynamic is much slower than any single kinetic step in the myofilament system and is also much slower than the distortion response (4, 6).
Despite the differences in derivation and interpretation, the mathematical structure of our model is similar to the mathematical structure of Kawai et al.'s model (14, 15, 17). For instance, Kawai et al.'s model for cardiac muscle sinusoidal responsiveness, over the frequencies normally examined, takes the form of a five-parameter complex stiffness given by
 | (6) |
H, Eq. 6 reduces to a four-parameter version
 | (7) |
In a previous theoretical publication where no data were analyzed (6), we considered that reduction of the linearized recruitment-distortion model was possible and left open the question of the number of parameters needed in such a model to fit the data. Our present results now demonstrate that reduction of each of the recruitment and distortion blocks to first-order dynamics, resulting in a model with only four parameters, is adequate to fit force-response data. Because virtually all the variation in F(t) due to L(t) is accounted for by the four-parameter model, the addition of a fifth parameter would do nothing to improve the fit to the data and would result in overparameterization of the model with degradation of the precision and accuracy of parameter estimates.
Model application method. An additional novel feature of this study is our experimental approach, which relies on model fitting in the time domain. Because we opted to work in the time domain, our model was applied to data as the solutions of its differential equations. This differs from previous work where the chosen model is expressed in a frequency-domain formulation of complex stiffness and then applied to sinusoidal steady-state data. The advantage of the time-domain approach is that there is no need to achieve sinusoidal steady state and there is no need to deliver perfect sinusoids as forcing functions. The time-domain approach accommodates transient behaviors as well as forcing functions of any shape providing there is sufficient information in these forcing functions to elicit behavior in all time scales of importance to the model identification problem. Because of the need to elicit responses over a broad range of time scales, we opted to use a two-part chirp as a forcing function. One part of the chirp, of 40-s duration, swept through low frequencies between 0.1 and 4 Hz and elicited force responses in time scales where first recruitment dynamics dominated and then distortion dynamics began to appear (see below). The second part of the chirp, of 5-s duration, swept through frequencies between 1 and 40 Hz where distortion was the dominant dynamic process. Thus the information needed to characterize the dynamic features of the system can be collected in a short period of time during a single activation.
Advantages
Whereas we take the model dependency of our time-domain approach to be an advantage, others will perceive it to be a disadvantage. For example, results from steady-state sinusoidal forcings or from approximation of the steady state with continuous changing frequency (13, 26) or from frequency analysis of responses to pseudorandom binary sequence perturbations (13, 25, 30) can be presented in model-independent format as complex stiffness frequency spectra. However, interpretations of such results must eventually rely on some understanding of the system, i.e., a model. Although advocates of the sinusoidal method may view model-independent results as a benefit because the user is not locked into a structured, and necessarily, simplified representation of the system, there are substantial dangers to a model-independent approach when uncritical and nonformal model interpretations are applied. To wit, we point to the common, but erroneous, interpretation of the frequency of minimum complex stiffness (fmin) as a measure of the XB cycling rate.
To make our point, we show in Fig. 9 the model-predicted stiffness spectra for the - and -MHC fibers whose force responses are shown in Fig. 2. Stiffness magnitude and phase, found by us using our time-domain experimental and analytic techniques, exhibits all the essential features of cardiac muscle stiffness reported by others using model-independent frequency-domain approaches (2, 3, 7, 13, 19, 21, 26, 29, 32). These features include the following. First, there is a quasiplateau (at least in a log-frequency scale) in (j) magnitude at the lowest frequencies. Second, as frequencies increase from zero, (j) magnitude first decreases to a minimum at fmin and then rises to a high-frequency plateau that is approximately three times the low-frequency plateau. Finally, at frequencies less than fmin, (j) exhibits a negative phase; at frequencies greater than fmin, (j) exhibits a positive phase (phase data not shown).
|
We decomposed (j) of the -MHC fiber into its recruitment and distortion components (Fig. 9, dashed lines). It can be seen that recruitment is the dominant process responsible for (j) at frequencies less than fmin, giving rise to the negative phase and falling (j) magnitude with increasing frequency in the low-frequency range. Furthermore, distortion is the dominant process responsible for (j) at frequencies greater than fmin, giving rise to the positive phase and rising (j) magnitude with increasing frequency at high frequencies. Minimum (j) at fmin arises because the contribution by recruitment becomes progressively less effective as frequency increases, and this loss of effectiveness occurs at frequencies lower than those at which distortion begins to effectively contribute to (j). Thus fmin represent the frequency where recruitment dominance gives way to distortion dominance in determining (j).
The appearance of a minimal (j) (and, consequently, fmin) depends on the characteristic frequency of the recruitment process, b, being much less than the characteristic frequency of the distortion process, c. Note that b has a functional dependence on all the kinetic coefficients in Fig. 1, which also includes the effect of length dependencies and cooperativities in modifying these coefficients (6, 24). Thus b obviously includes much more than XB cycling. On the other hand, c depends most strongly on the rate constant governing XB detachment (g; Fig. 1). Thus if one takes detachment to be a limiting step in the XB cycle, c may be used as a measure of XB cycling speed. Therefore, fmin is a marker of transition of dynamic dominance from processes responsible for global myofilament dynamics to dominance by processes within the XB cycle, but it is not an explicit statement of XB cycling speed.
Limitations
A limitation in our methodology occurs when b and c are not different in value. When the c-to-b ratio is <2, there becomes a less obvious minimum in (j). Consequently, it becomes difficult or impossible to identify fmin. An anecdotal observation from this experiment, which is consistent with the findings of Shibata et al. (27), is that skinned fiber preparations that otherwise do not meet the criteria for acceptable viability also do not exhibit a minimum (j). Essentially, this means that b and c have come close together in these fibers of poor viability. We encountered this condition in one of the ferret fibers included in this data set (data from this particular fiber were not discarded because the fiber met other criteria for preparational viability). The force response of this fiber exhibited no discernible dip, and the model estimates of b and c were 11.5 s–1 and 11.2 s–1, respectively. The standard error of the b estimate was 12.6% of the estimated value, and this standard error was more than 10 times the standard error of the b estimates in all 112 records obtained from the other 18 fibers!
Compared with the other ferret fibers, where the average values of b and c were, at comparable pCa, 4.84 s–1 and 11.23 s–1, respectively, this fiber with no dip had an abnormally high b but a normal c. We suspect that the abnormally elevated b in this fiber was due to fiber deterioration during preparation of the specimen and further suspect that such deterioration is the primary reason for the loss of discernible (j) minimum in most fibers where no minimum is observed.
Why would b be higher than normal in a deteriorating fiber? Recall that the internal myofilament regulatory processes, such as length sensing and cooperativity, participate importantly in determining the magnitude and speed of recruitment, E0 and b, respectively. The effect of these internal regulations is to enhance and slow the XB recruitment response of the system to changes in length (7). If these internal regulatory processes are lost, E0 will decrease and the recruitment process will speed up, i.e., b will increase. Indeed, E0 at comparable pCa in the three ferret fibers with discernible dip in the force response was 300 mN/mm2, whereas it was only 235 mN/mm2 in the fiber without a discernible dip.
Thus we suspect that one of the first mechanisms to be lost in a deteriorating fiber is the internal myofilament regulatory processes that control length-dependent activation. The interactions among proteins responsible for this regulation may be more sensitive to factors that sustain viability than is the interaction between actin and myosin in the formation of a XB. The problem of fiber viability is not related to the model or the analysis that we used here. Rather, it is a problem general to all experiments that use skinned fiber preparations. The experimental method that we employed gives the experimenter a quick way to evaluate the preparation while data are being collected. If no dip is visible in the force record, then the fiber may be discarded from the experiment.
Applications to Physiology
There were many physiologically interesting items in the data collected for this methodological study. However, in keeping with the methodological theme, we mention only one physiological item to illustrate the power and utility of the method for physiological studies. We obtained dynamic FLR records from control rats expressing 80% -MHC and 20% -MHC and from PTU-treated rats expressing virtually 100% -MHC. As expected, the value of the parameter most closely associated with MHC kinetics, c, differed markedly between these two groups of animals: mean of 30.3 s–1 in the control group versus 13.2 s–1 in the PTU-treated group at pCa 4.3. Also, however, parameter b varied between these groups: mean of 5.35 s–1 for the control group versus 3.75 s–1 for the PTU-treated group at pCa 4.3. There are three possible explanations as to why b (which is a global myofilament dynamic constant) is affected with changes in MHC isoform expression. First, even though b is a global constant, myosin kinetics enter into the mix that determines its value; therefore, a slowing of myosin kinetics is expected to decrease b. Second, the expression of genes coding for myofilament proteins may not be independent of each other. Thus a change in the regulation of MHC gene expression, such that the slow -isoform of MHC is expressed, may be linked to a change in the regulation of other myofilament proteins or to factors controlling protein interaction such that overall XB recruitment dynamics are slowed along with slowed XB detachment. Finally, cooperativity within the myofilament system may be such that slowed XB detachment causes slowed recruitment through cooperative mechanisms that act on internal myofilament processes. For instance, XB-based cooperativity would be enhanced with longer dwell time in the attached state, as would occur with decreased c. This enhanced cooperativity could then act to slow b. Indeed, recruitment magnitude, as measured by E0, increased when c decreased; average E0 at pCa 4.3 of the fibers from the PTU-treated rat was 261 mN/mm2 compared with 162 mN/mm2 for fibers from the control rat. This is consistent with cooperativity being greater in the fiber with -MHC. These three possibilities are all potential mechanisms by which the dynamics of recruitment and distortion are coordinated in the myofilament system to preserve important functional features of the fiber such as a minimum in (j) magnitude and fmin.
As this example demonstrates, outcomes from the application of the recruitment-distortion model to physiological studies can be 1) the illustration of specific differences induced by treatment to test an a priori hypothesis and 2) unanticipated observations for which alternative hypotheses may be posited and future experimental studies proposed.
In summary, we proposed a recruitment-distortion model for describing and explaining dynamic FLRs in constantly activated cardiac muscle. We first thoroughly evaluated this model with respective to descriptive validity criteria, asking whether the model could account for most of the observed variation in force induced by changes in length. It did. We further asked whether the unaccounted-for variation in force was related to length changes. It was not. Finally, we asked whether model parameters could be estimated with sufficient precision that distortion parameters could be uniquely separated from recruitment parameters and also that the parameters could be considered unique to the preparation from which they were estimated. They could. We then evaluated the model with respect to explanative validity criteria, asking whether specific parameters were positively correlated with independent assessments of the physiological entity they supposedly represented. They were. We compared the recruitment-distortion model with competing models now in current use and found similarities in mathematical structure but differences in interpretations. Among those models that we reviewed here, only the recruitment-distortion model incorporates the full spectrum of myofilament kinetic processes responsible for both activation and XB cycling. Therefore, only the recruitment-distortion model takes into account the effect of SL on myofilament activation. The recruitment-distortion model avoids potentially erroneous interpretations of complex stiffness and is a useful construct for describing and explaining dynamic FLRs in constantly activated cardiac muscle.
|
FOOTNOTES |
|---|
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. Section 1734 solely to indicate this fact.
|
REFERENCES |
|---|
|
TOPABSTRACTMETHODSRESULTSDISCUSSIONREFERENCES |
|---|
