# Distribution of Myosin Attachment Times Predicted from Viscoelastic Mechanics of Striated Muscle

Authors: Bradley M. Palmer; Yuan Wang; Mark S. Miller
Journal: Journal of Biomedicine and Biotechnology
 (2011)

---

## Abstract

We demonstrate that viscoelastic mechanics of striated muscle, measured as elastic and viscous moduli, emerge directly from the myosin crossbridge attachment time,tatt, also called time-on. The distribution of tatt was modeled using a gamma distribution with shape parameter, p, and scale parameter, β. At 5 mM MgATP, β was similar between mouse α-MyHC (16.0±3.7 ms) and β-MyHC (17.9±2.0 ms), and p was higher (P<0.05) for β-MyHC (5.6±0.4 no units) compared to α-MyHC (3.2±0.9). At 1 mM MgATP, p approached a value of 10 in both isoforms, but β rose only in the β-MyHC (34.8±5.8 ms). The estimated mean tatt (i.e., pβ product) was longer in the β-MyHC compared to α-MyHC, and became prolonged in both isoforms as MgATP was reduced as expected. The application of our viscoelastic model to these isoforms and varying MgATP conditions suggest that tatt is better modeled as a gamma distribution due to its representing multiple temporal events occurring within tatt compared to a single exponential distribution which assumes only one temporal event within tatt.


---

## Body

## 1. Introduction

Advances in optical techniques have allowed detailed analysis of the isolated single myosin crossbridge. Upon formation of a myosin crossbridge, isomerization of the myosin molecule provides a unitary force proportional to the length of the lever arm and the extent of the lever arm swing, also called the power stroke [1–4]. The time duration of crossbridge attachment, sometimes called time-on, is also measured by optical trapping. The distribution of attachment times has been found to depend upon the myosin isoform, nucleotide availability, point mutations, and several other factors [5–12]. This time duration of crossbridge attachment plays a significant role in determining muscle performance observed as force and velocity of contraction at the single fiber level.There remains, however, a considerable challenge to understand and describe how the unitary force and the temporary attachment times of myosin molecules are manifested at the level of muscle tissue, which possesses a three-dimensional lattice structure and can be studied intact or submerged in physiological solutions after removal of the plasma membrane. Viscoelastic mechanics of muscle tissue represents one such macroscopic consequence of these molecular phenomena. We describe in this paper a quantitative justification and methodology for estimating the distribution of myosin crossbridge attachment times based on viscoelastic mechanics measured in striated muscle fibers.Length perturbation analysis of muscle provides a means to quantify viscoelastic mechanics and entails applying a small length change at one end of a muscle and recording the force response on the other end. The length perturbation and force response are then subjected to Fourier transformation, and the dynamic mechanical properties of the muscle are characterized by the complex ratio of the transformed force response normalized to cross-sectional area (orσ(ω)= Fourier transformed tensile stress) divided by the transformed length perturbation normalized to initial length (or ε(ω)= Fourier transformed strain). This complex ratio is also called the mechanical transfer function or complex modulus of the muscle, Y(ω)=σ(ω)/ε(ω). The real and imaginary parts are, respectively, termed the elastic modulus and viscous modulus of the muscle [13, 14].The frequency characteristic dips (low values) and shoulders (high values) observed in the elastic and viscous moduli of activated muscle are sensitive to species of origin, myosin isoform, and varying concentrations of MgATP and therefore reflect the macroscopic consequence of myosin crossbridge kinetics [4, 11, 15–18]. The dips in the moduli, most notably the negatively valued viscous modulus, occur near the frequencies at which the muscle operates in vivo[14]. The shoulders of the moduli, however, are more prominent in magnitude. According to our previous modeling endeavor [19], the shoulders appearing at these higher frequencies of the moduli reflect the mechanical consequences of intermittent myosin crossbridge formation. In that work, we proposed that a two-state model of the acto-myosin crossbridge governed by first-order kinetics gives rise to a viscoelastic work-absorbing property, termed the C-process by Kawai and colleagues [13, 16, 20], which is characterized by an exponential rate constant equivalent to the myosin crossbridge off rate termed g by Huxley [21]. The mean myosin attachment time, based on a single exponential distribution of attachment times, could be estimated as the reciprocal of this exponential rate constant, 2πc, after fitting (1) to a measured complex modulus [14, 19]:(1)Y(ω)=A(iω)k-B(iω2πb+iω)+C(iω2πc+iω).While the single exponential representation of C-process and its interpretation have been valuable for examining a variety of muscle types under a number of conditions [7, 11, 13, 14, 19, 20, 22], the assumptions of first-order kinetics and a single exponential distribution of myosin attachment times are limiting. It is known, for example, that multiple time periods associated with multiple biochemical states make up the myosin crossbridge cycle. These multiple states constitute the entire myosin crossbridge attachment time. Such an addition of multiple smaller time periods would result in a multiple exponential distribution or a gamma distribution of the crossbridge attachment times rather than a single exponential distribution assumed by first-order kinetics.In the present study, we examine the mechanical consequences of the distribution of myosin attachment times represented by a gamma distribution, which allows the description of more discrete and longer lived time-on and yet also allows for the possibility of a single exponential function. The consideration of the gamma distribution effectively poses the hypothesis that multiple temporal events which occur within the time of attachment can be discerned in the viscoelastic response of striated muscle at the macroscopic level. We provide an analytical solution to the mechanical consequences that emerge as the shoulders of the elastic and viscous moduli. We also demonstrate the validity of the analytical solution with computer simulations. Finally, we demonstrate the application of the gamma distribution representation of the C-process in a comparison of mouse cardiac myosin heavy chain isoforms,α-MyHC and β-MyHC, subjected to varying concentrations of MgATP.


## 2. Mathematical Modeling

Our goal in this section is to derive a mathematical expression that describes the viscoelastic component of a force response to a length perturbation that has been applied to an ensemble of myosin molecules intermittently attaching and detaching from actin within a half sarcomere of striated muscle. The force response and length perturbation will form the bases for a mechanical transfer function, which is equivalent to a viscoelastic complex modulus measured in frequency space. We use a two-state model of myosin attachment and detachment as our starting point in modeling the mechanical consequences of intermittently attached actomyosin crossbridges. It is important to point out that in the following development we do not assume first-order kinetics as governing myosin attachment and detachment; in doing so we will be able to provide a more generalized model of the viscoelastic mechanics due to the myosin crossbridges.The time periods of attachment,tatt, and detachment, tdet, are considered here to be random variables whose values are governed by stochastic processes independent of force, stress, strain, length, velocity, and each other. The total time period for a myosin crossbridge cycle is also a random variable, tcycle=tatt+tdet. All myosin heads are assumed independent of each other. During detachment, force due to an individual myosin head is zero. During attachment, force is given as the product of the unitary force due to the power stroke, the length displacement of an elastic element in series with the crossbridge, and the stiffness of that elastic element (Figure 1). We assume that length perturbations are sufficiently small to elicit a linear force response of the elastic element in series with the crossbridge.Force due to myosin crossbridges. (a) A unitary force,Funi, is generated over a finite time of attachment, tatt (also called time-on), when the myosin power stroke swings the lever arm a net unitary displacement, duni. The time of detachment, tdet, refers to the time period between a crossbridge detaching from actin and before reattaching. (b) An additional crossbridge-dependent force is transmitted through the half sarcomere when an elastic element in series with the crossbridge and with stiffness kstiff undergoes a length change due to an externally applied perturbation of the half sarcomere, Lhs. The force recorded at any time t, f(t), due to the perturbation reflects the displacement of the elastic element since the initial time of crossbridge attachment, tini.
(a)(b)The elastic element represents the most compliant structures between the M- and Z-lines of the sarcomere in series with the attached crossbridge. Previous work by others has identified the most compliant structures as the head and neck regions of the myosin S1 segment including the lever arm bound with the essential and regulatory light chains and that portion of myosin S2 segment not incorporated into the thick filament backbone [4, 23, 24]. We also assume that when the myosin is in a postpower stroke state the most compliant elements are taut, that is, not slack as might occur if the myosin bound to actin at a point too close to the M-line to stretch the S2 segment [4, 23, 24].


### 2.1. Viscoelastic Mechanical Response of Half Sarcomere

Here we provide a theoretical description of the force response of a half sarcomere to an externally applied length perturbation. We assume here that the total force that might be recorded from a half sarcomere,Ftotal(t), is given as the sum of the forces produced by each attached myosin crossbridge due to its unitary force, Funi, and its dynamic response due to the length perturbation, f(t) (Figure 1). The total force produced by a summation of attached crossbridges will be denoted as(2)Ftotal(t)=∑i=1MFuni(i)+f(i)(t),
where M= the number of attached crossbridges in the half sarcomere at time t. When an external length perturbation is applied, the dynamic response of the ith attached crossbridge is given as the stiffness-displacement product, which is a linear approximation suitable for very small length perturbations:(3)f(i)(t)=kstiff(i)(Lhs(t)-Lhs(tini(i))),
where kstiff = the stiffness coefficient of the crossbridge elastic element, Lhs(t) = displacement imposed upon of the half sarcomere by an externally driven perturbation, and tini = the instant of initial crossbridge attachment. Under isometric conditions, the total force is simply equal to the total number of attached crossbridges, M, multiplied by the mean unitary force, F̅uni. The dynamic force recorded over the half sarcomere, Fhs(t), can then be defined as the total force minus the isometric force and written as follows:(4)Fhs(t)=Ftotal(t)-MF¯uni=∑i=1Mkstiff(i)(Lhs(t)-Lhs(tini(i))).We now use the Inverse Fourier Transform definitions forFhs(t) and Lhs(t), namely,(5)Fhs(t)=12π∫-∞∞F̃hs(ω)eiωtdω,Lhs(t)=12π∫-∞∞L̃hs(ω)eiωtdω,
where F̃hs(ω) and L̃hs(ω) are the Fourier Transform representations of Fhs(t) and Lhs(t), respectively; then (4) can be written as follows:(6)12π∫-∞∞F̃hs(ω)eiωtdω=12π∫-∞∞∑i=1Mkstiff(i)L̃hs(ω)(eiωt-eiωtini(i))dω.
We now define the variable, τ=t-tini, as the time period between time, t, and the instant of the most recent attachment of a crossbridge, tini. For the ith crossbridge, tini(i)=t-τ(i). The term (eiωt-eiωtlni(i)) in (6) then becomes (1-e-iωτ(i))eiωt and we can remove the Inverse Fourier Integrals of (6):(7)F̃hs(ω)=∑i=1Mkstiff(i)L̃hs(ω)(1-e-iωτ(i)).The summation over a large number of attached crossbridges can be replaced as the product ofM and the expected value of the summed terms. To do so would require our providing the probability density functions (PDFs) for those random variables represented as bearing the i-subscript in (7). We assume here that the random variables kstiff and τ are independent of each other. Thus, the result for the expected value of kstiff would be independent of its probability density function and equal to the mean value for kstiff. Equation (7) can be written as(8)F̃hs(ω)=Mk¯stiffL̃hs(ω)∫0∞PDFτ(t)(1-e-iωt)dt,
where k̅stiff = mean stiffness of the elastic element and PDFτ(t)= the probability density function for the random variable τ as a function of time, t. Equation (8) represents the force response in frequency space due to an externally applied length perturbation on striated muscle with no strain dependence on myosin crossbridge kinetics.The number of crossbridges attached at any time,M, can be replaced by the total number of myosin heads available to form crossbridges, N, multiplied by the probability that any one myosin head has formed a crossbridge, that is, the duty ratio t̅att/t̅cycle. The mechanical transfer function, which represents a measured viscoelastic complex modulus, can then be written as follows:(9)F̃hs(ω)L̃hs(ω)=N(t̅attt̅cycle⁡)k¯stiff∫0∞PDFτ(t)(1-e-iωt)dt,
where t̅att = mean time duration of crossbridge attachment and t̅cycle = mean time of a complete crossbridge cycle. Equation (9) represents a generalized mathematical description of the visco-elastic response of intermittently attaching myosin crossbridges. The form of (9) suggests that the force response is dictated by PDFτ(t), so long as the other random variables are sufficiently represented by their respective means. A physical interpretation of PDFτ(t) and its relation to the probability density function for tatt, PDFtatt(t), are presented below.


### 2.2. Survival Function

As defined above, the random variableτ represents the time period between any given time t and the instant any bound crossbridge was formed (Figure 1(b)). In other words, τ represents the time period any bound crossbridge has survived up to time t. The probability (Pr) that a crossbridge has survived for a time period τ, we will call it Sτ(t), is defined as the probability that τ is less than the time of crossbridge attachment, tatt [25]. The probability density function of crossbridge attachment time PDFtatt(t) and the survival function Sτ(t) are related as follows:(10a)Sτ(t)=Pr(t<tatt)=∫t∞PDFtatt(t)dt,(10b)-dSτ(t)dt=PDFtatt(t).The PDFτ(t) needed in (8) is simply a normalization of Sτ(t) and is defined as the survival function of (10a) divided by its integral:(11)PDFτ(t)=Sτ(t)∫0∞Sτ(t)dt.With (9) we have provided a mathematical representation of the energy-absorbing viscoelastic complex modulus that arises from intermittently attached myosin crossbridges, and we have not assumed a specific scheme for the biochemical steps that govern the distributions of the random variable tatt. Thus, any measure of the viscoelastic complex modulus, given by the left hand side of (9), can be used to calculate an estimate of PDFτ(t) and by extension through (10b) also be used to calculate an estimate of PDFtatt(t).


### 2.3. Gamma Distribution Representation ofPDFtatt(t)

In our previous work we chosePDFtatt(t) to be represented by a single exponential distribution [19], which by extension of (10a), (10b), and (11) also defined PDFτ(t) as the same single exponential distribution. That choice reflected an assumption of first-order kinetics governing the lifetime of the myosin crossbridge in the attached state. The result was a mathematical representation of the mechanical consequences of length perturbation that was equivalent to the historically represented C-process shown in (1).Here we do not assume first-order kinetics. Instead, we considerPDFtatt(t) to be represented by a gamma distribution, which is used to describe random variables that represent a time period made up of the sum of several smaller time periods [26, 27](12)PDFtatt(t)=1βΓ(p)(tβ)p-1e-t/β,
where p= a shape parameter reflecting the number of smaller time periods summed together to produce total time attached, tatt, β= a scale parameter reflecting the average duration of the smaller time periods. The mean tatt is calculated as pβ, and Γ(p)= the gamma function evaluated at p, a normalization factor.The two parametersp and β are sufficient to provide a gamma distribution description of PDFtatt(t). The case of p=1 is a special case that is equivalent to a single exponential function with mean β. Several examples of gamma distributions, including a case of p=1, are illustrated in Figures 2(a) and 2(b), which demonstrate the effectiveness of the gamma distribution representation in accommodating a wide variety of probability density functions for tatt. Upon applying (10a), we attain the associated survival function Sτ(t):(13a)Sτ(t)=Γ(p,t/β)Γ(p),
where Γ(p,t/β)= the upper incomplete gamma function [26]. Example survival functions are illustrated in Figures 2(c) and 2(d) and correspond to the PDFtatt(t) examples in Figures 2(a) and 2(b). The integral of Sτ(t) over the abscissa, which is useful for normalization, is(13b)∫0∞Sτ(t)dt=pβ.Figure 2
Examples of gamma distributions and associated survival functions. (a, b) The two parameters of the gamma distribution (p and β) permit description of a wide variety of distributions useful for representing the probability density function for myosin crossbridge attachment time, PDFtatt(t). Notably, a single exponential results when p=1 (a), and a Gaussian distribution is approximated as p increases. (c, d) The survival function, Sτ(t), refers to the probability that an attached crossbridge will survive to time t. The analytical relationship between the PDFtatt(t) and Sτ(t) is provided by (10a) and (10b).Upon applying (13a) and (13b) to (11) we have(14)PDFτ(t)=Γ(p,t/β)pβΓ(p).
We can now state that the complex modulus that would arise from an ensemble of myosin crossbridges, whose distribution of tatt can be represented by a gamma distribution presented in (12), is given by the insertion of PDFτ(t) of (14) into the expression for the mechanical transfer function of (9). Details of the evaluation of the integral in (9) after such a substitution are provided in the Appendix. The result is(15)F̃hs(ω)L̃hs(ω)=ξ{1-1piωβ((iωβ+1)p-1(iωβ+1)p)},
whereξ=N(t¯att/t¯cycle)k¯stiff.Equation (15) describes the predicted complex modulus that would arise from length perturbation analysis of an ensemble of myosin crossbridges whose distribution of attachment times is represented by a gamma distribution. From a measured viscoelastic complex modulus represented on the left-hand side of (15) we could use nonlinear least-squares methods to estimate the values of the three independent parameters ξ,p, and β, on the right-hand side of (15).


## 2.1. Viscoelastic Mechanical Response of Half Sarcomere

Here we provide a theoretical description of the force response of a half sarcomere to an externally applied length perturbation. We assume here that the total force that might be recorded from a half sarcomere,Ftotal(t), is given as the sum of the forces produced by each attached myosin crossbridge due to its unitary force, Funi, and its dynamic response due to the length perturbation, f(t) (Figure 1). The total force produced by a summation of attached crossbridges will be denoted as(2)Ftotal(t)=∑i=1MFuni(i)+f(i)(t),
where M= the number of attached crossbridges in the half sarcomere at time t. When an external length perturbation is applied, the dynamic response of the ith attached crossbridge is given as the stiffness-displacement product, which is a linear approximation suitable for very small length perturbations:(3)f(i)(t)=kstiff(i)(Lhs(t)-Lhs(tini(i))),
where kstiff = the stiffness coefficient of the crossbridge elastic element, Lhs(t) = displacement imposed upon of the half sarcomere by an externally driven perturbation, and tini = the instant of initial crossbridge attachment. Under isometric conditions, the total force is simply equal to the total number of attached crossbridges, M, multiplied by the mean unitary force, F̅uni. The dynamic force recorded over the half sarcomere, Fhs(t), can then be defined as the total force minus the isometric force and written as follows:(4)Fhs(t)=Ftotal(t)-MF¯uni=∑i=1Mkstiff(i)(Lhs(t)-Lhs(tini(i))).We now use the Inverse Fourier Transform definitions forFhs(t) and Lhs(t), namely,(5)Fhs(t)=12π∫-∞∞F̃hs(ω)eiωtdω,Lhs(t)=12π∫-∞∞L̃hs(ω)eiωtdω,
where F̃hs(ω) and L̃hs(ω) are the Fourier Transform representations of Fhs(t) and Lhs(t), respectively; then (4) can be written as follows:(6)12π∫-∞∞F̃hs(ω)eiωtdω=12π∫-∞∞∑i=1Mkstiff(i)L̃hs(ω)(eiωt-eiωtini(i))dω.
We now define the variable, τ=t-tini, as the time period between time, t, and the instant of the most recent attachment of a crossbridge, tini. For the ith crossbridge, tini(i)=t-τ(i). The term (eiωt-eiωtlni(i)) in (6) then becomes (1-e-iωτ(i))eiωt and we can remove the Inverse Fourier Integrals of (6):(7)F̃hs(ω)=∑i=1Mkstiff(i)L̃hs(ω)(1-e-iωτ(i)).The summation over a large number of attached crossbridges can be replaced as the product ofM and the expected value of the summed terms. To do so would require our providing the probability density functions (PDFs) for those random variables represented as bearing the i-subscript in (7). We assume here that the random variables kstiff and τ are independent of each other. Thus, the result for the expected value of kstiff would be independent of its probability density function and equal to the mean value for kstiff. Equation (7) can be written as(8)F̃hs(ω)=Mk¯stiffL̃hs(ω)∫0∞PDFτ(t)(1-e-iωt)dt,
where k̅stiff = mean stiffness of the elastic element and PDFτ(t)= the probability density function for the random variable τ as a function of time, t. Equation (8) represents the force response in frequency space due to an externally applied length perturbation on striated muscle with no strain dependence on myosin crossbridge kinetics.The number of crossbridges attached at any time,M, can be replaced by the total number of myosin heads available to form crossbridges, N, multiplied by the probability that any one myosin head has formed a crossbridge, that is, the duty ratio t̅att/t̅cycle. The mechanical transfer function, which represents a measured viscoelastic complex modulus, can then be written as follows:(9)F̃hs(ω)L̃hs(ω)=N(t̅attt̅cycle⁡)k¯stiff∫0∞PDFτ(t)(1-e-iωt)dt,
where t̅att = mean time duration of crossbridge attachment and t̅cycle = mean time of a complete crossbridge cycle. Equation (9) represents a generalized mathematical description of the visco-elastic response of intermittently attaching myosin crossbridges. The form of (9) suggests that the force response is dictated by PDFτ(t), so long as the other random variables are sufficiently represented by their respective means. A physical interpretation of PDFτ(t) and its relation to the probability density function for tatt, PDFtatt(t), are presented below.


## 2.2. Survival Function

As defined above, the random variableτ represents the time period between any given time t and the instant any bound crossbridge was formed (Figure 1(b)). In other words, τ represents the time period any bound crossbridge has survived up to time t. The probability (Pr) that a crossbridge has survived for a time period τ, we will call it Sτ(t), is defined as the probability that τ is less than the time of crossbridge attachment, tatt [25]. The probability density function of crossbridge attachment time PDFtatt(t) and the survival function Sτ(t) are related as follows:(10a)Sτ(t)=Pr(t<tatt)=∫t∞PDFtatt(t)dt,(10b)-dSτ(t)dt=PDFtatt(t).The PDFτ(t) needed in (8) is simply a normalization of Sτ(t) and is defined as the survival function of (10a) divided by its integral:(11)PDFτ(t)=Sτ(t)∫0∞Sτ(t)dt.With (9) we have provided a mathematical representation of the energy-absorbing viscoelastic complex modulus that arises from intermittently attached myosin crossbridges, and we have not assumed a specific scheme for the biochemical steps that govern the distributions of the random variable tatt. Thus, any measure of the viscoelastic complex modulus, given by the left hand side of (9), can be used to calculate an estimate of PDFτ(t) and by extension through (10b) also be used to calculate an estimate of PDFtatt(t).


## 2.3. Gamma Distribution Representation ofPDFtatt(t)

In our previous work we chosePDFtatt(t) to be represented by a single exponential distribution [19], which by extension of (10a), (10b), and (11) also defined PDFτ(t) as the same single exponential distribution. That choice reflected an assumption of first-order kinetics governing the lifetime of the myosin crossbridge in the attached state. The result was a mathematical representation of the mechanical consequences of length perturbation that was equivalent to the historically represented C-process shown in (1).Here we do not assume first-order kinetics. Instead, we considerPDFtatt(t) to be represented by a gamma distribution, which is used to describe random variables that represent a time period made up of the sum of several smaller time periods [26, 27](12)PDFtatt(t)=1βΓ(p)(tβ)p-1e-t/β,
where p= a shape parameter reflecting the number of smaller time periods summed together to produce total time attached, tatt, β= a scale parameter reflecting the average duration of the smaller time periods. The mean tatt is calculated as pβ, and Γ(p)= the gamma function evaluated at p, a normalization factor.The two parametersp and β are sufficient to provide a gamma distribution description of PDFtatt(t). The case of p=1 is a special case that is equivalent to a single exponential function with mean β. Several examples of gamma distributions, including a case of p=1, are illustrated in Figures 2(a) and 2(b), which demonstrate the effectiveness of the gamma distribution representation in accommodating a wide variety of probability density functions for tatt. Upon applying (10a), we attain the associated survival function Sτ(t):(13a)Sτ(t)=Γ(p,t/β)Γ(p),
where Γ(p,t/β)= the upper incomplete gamma function [26]. Example survival functions are illustrated in Figures 2(c) and 2(d) and correspond to the PDFtatt(t) examples in Figures 2(a) and 2(b). The integral of Sτ(t) over the abscissa, which is useful for normalization, is(13b)∫0∞Sτ(t)dt=pβ.Figure 2
Examples of gamma distributions and associated survival functions. (a, b) The two parameters of the gamma distribution (p and β) permit description of a wide variety of distributions useful for representing the probability density function for myosin crossbridge attachment time, PDFtatt(t). Notably, a single exponential results when p=1 (a), and a Gaussian distribution is approximated as p increases. (c, d) The survival function, Sτ(t), refers to the probability that an attached crossbridge will survive to time t. The analytical relationship between the PDFtatt(t) and Sτ(t) is provided by (10a) and (10b).Upon applying (13a) and (13b) to (11) we have(14)PDFτ(t)=Γ(p,t/β)pβΓ(p).
We can now state that the complex modulus that would arise from an ensemble of myosin crossbridges, whose distribution of tatt can be represented by a gamma distribution presented in (12), is given by the insertion of PDFτ(t) of (14) into the expression for the mechanical transfer function of (9). Details of the evaluation of the integral in (9) after such a substitution are provided in the Appendix. The result is(15)F̃hs(ω)L̃hs(ω)=ξ{1-1piωβ((iωβ+1)p-1(iωβ+1)p)},
whereξ=N(t¯att/t¯cycle)k¯stiff.Equation (15) describes the predicted complex modulus that would arise from length perturbation analysis of an ensemble of myosin crossbridges whose distribution of attachment times is represented by a gamma distribution. From a measured viscoelastic complex modulus represented on the left-hand side of (15) we could use nonlinear least-squares methods to estimate the values of the three independent parameters ξ,p, and β, on the right-hand side of (15).


## 3. Methods

### 3.1. Computer Simulations

The force response of a virtual half sarcomere was simulated for a time period of two seconds for many (20,000) independent myosin heads alternately attaching and detaching to actin according to independent stochastic processes governing the random variablestatt and tdet. For any one computer simulation, random numbers representing tatt were generated to conform with a gamma distribution with values for p ranging from 1 to 48 and for β ranging from 2 to 24 such that the mean time attached was either 24 ms or 96 ms. Random numbers representing tdet were generated with values for p=2 and β=48 ms such that the mean time detached was always 96 ms. Figure 3(a) illustrates the PDFtatt(t) when p=3 and β=8. Each random number generated from this PDFtatt(t) was used to dictate the time period over which a virtual crossbridge was attached for each occasion, as indicated in Figure 3(b). Figure 3(c) illustrates the PDF for tdet, PDFtdet(t), with p=2 and β=48. Each random number generated from this PDFtdet(t) was used to dictate the time over which a virtual crossbridge was detached for each occasion. Figure 3(b) illustrates an example of successive attachment and detachment time periods.Figure 3
Schematic representation of computer-generated elastic and viscous moduli. (a)–(c) Random numbers generated from variousPDFtatt(t) (one example shown in (a)) were used to dictate the attachment times, tatt, of a virtual crossbridge. (b) shows an example attachment and detachment sequence. (c) illustrates the PDFtdet(t), with p=2 and β=48, used to dictate the intermittent times of detachment. (d) A sinusoidal length perturbation with amplitude 1 nm was used to simulate an externally applied length perturbation. This example is at 1 Hz. (e) An example series of successive crossbridge attachments and detachments would be subject to the length perturbation. (f) The length change that occurs during the time of attachment results in a force deflection. (g) The sum of these force deflections from 20,000 crossbridges is then fit to a sinusoidal (solid smooth line). (h) The magnitude of the component of the fitted sinusoid in phase with the length perturbation represents the elastic modulus at this frequency, the magnitude of the component out of phase represents the viscous modulus.Sinusoidal length perturbations of the virtual half sarcomere were simulated as having amplitude 1 nm and frequencies over the range 1–250 Hz. Figure3(d) illustrates an example perturbation at 1 Hz over a one-second time interval. Each crossbridge was assigned a stiffness constant of 1 pN/nm [28]. The change in length of the crossbridge relative to the time of initial attachment was multiplied by the stiffness constant to simulate the force resulting from the strain on the elastic element of the crossbridge. Figure 3(e) illustrates an example series of successive attachments and detachments over one second. Figure 3(f) illustrates the force deflection that would occur for a crossbridge during each time of attachment. The resulting force deflections of 20,000, independent crossbridges were summed to provide an equivalent two-second period for the ensemble in the virtual half sarcomere. An example of the resultant force from 20,000 crossbridges is given in Figure 3(g). This simulated force of the half sarcomere, Fhs(t), was fit using a simplex method to a sine function, whose amplitude (Amp) and phase (ϕ) permitted the calculation of two components that were in phase and out of phase with respect to the length perturbation. The value of each moduli was then calculated from the amplitude and phase as elastic modulus = Amp cos(ϕ) and viscous modulus = Amp sin(ϕ). The smooth line in Figure 3(g) presents an example fit of a sine function to the resultant force, and Figure 3(h) presents the in-phase and out-of-phase components provided by the Simplex method. The resulting elastic and viscous moduli, given in units of pN/nm due to the lack of normalization factors used in the simulation, were fit simultaneously to (15) using a nonweighted Levenberg-Marquardt nonlinear least-squares routine. The diagonal of the resulting covariance matrix indicated the variances of parameter estimates. The square root of the variance thus provided the standard error for each parameter estimate, which was then multiplied 1.96 to provide the 95% confidence range. All computer simulations, random number generation, curve-fittings, and parameter estimation routines were performed using IDL Version 7.0 (ITT, Boulder, CO).


### 3.2. Solutions

All reagents were purchased from Sigma (St. Louis, Mo). Solutions were formulated by solving equations describing ionic equilibria [29]. Concentrations are expressed in mmol/L unless otherwise noted. Relaxing solution: pCa 8.0, 5.0 ethylene glycol-tetra-acetic acid (EGTA), 5.0 ATP, 1.0 Mg2+, 20 hydroxyethyl-aminoethanesulfonic acid (BES), 35 phosphocreatine, 300 U/mL creatine kinase, ionic strength 200, pH 7.0. Activating solution: same as relaxing solution with pCa  4.0. Rigor solution: same as activating without added ATP, creatine kinase, or phosphocreatine. Storage solution: same as relaxing with 10 μg/mL leupeptin and 50% wt/vol glycerol. Skinning solution: same as storage with 30 mM 2,3-butanedione monoxime (BDM) and 1% wt/vol Triton X-100.


### 3.3. Viscoelastic Mechanics

All procedures were reviewed and approved by the Institutional Animal Care and Use Committees of The University of Vermont. Male wild-type mice were fed either a normal mouse diet (WT) or an iodine deficient, 0.15% 6-N-propyl-2-thiouracil (PTU) diet for at least nine weeks prior to their being killed by rapid cervical dislocation. The PTU mice therefore became hypothyroid resulting in the expression ofβ-MyHC in the myocardium in contrast to the α-MyHC expressed in the WT mouse heart.Mouse left ventricular skinned myocardial strips were prepared using methods similar to those described previously to yield thin strips (~140μm diameter, ~800 μm length) with longitudinally oriented parallel fibers [7]. These strips were chemically skinned for 2 hr at 22°C and stored at −20°C for no more than 5 days. At the time of study, aluminum T-clips were attached to the ends of a strip ~150 μm apart. The strip was mounted between a piezoelectric motor (Physik Instrumente, Auburn, MA) and a strain gauge (SensoNor, Horten, Norway), lowered into a 30 μL droplet of relaxing solution maintained at 37°C and incrementally stretched to and maintained at 2.2 μm sarcomere length detected by videography and digital Fourier transform techniques (IonOptix, Milton, MA).Strips were calcium activated at pCa 4.5 and subjected to decreasing concentrations of 5, 2, and 1 mM MgATP by exchanging equal volumes of rigor solution. Sinusoidal perturbations of amplitude 0.125% strip length were applied over the frequency range 0.125–250 Hz. The elastic and viscous moduli were calculated from the recorded tension transient as the relative magnitudes of the in-phase and out-of-phase components with respect to the imposed sinusoidal length perturbations [13, 20, 30]. The measured complex modulus was fit to (1), representing the single exponential distribution model of tatt, and to (16) below, representing the gamma distribution model of tatt, using a non weighted Levenburg-Marquardt non linear least-squares routine running within IDL Version 7.0 (ITT, Boulder, CO)(16)Y(ω)=A(iω)k-B(iω2πb+iω)+ξ{1-1piωβ((iωβ+1)p-1(iωβ+1)p)}.Mean myosin time attached was calculated as(2πc)-1 when the fit was performed with (1) and as pβ for fits using (16). The correlation coefficients between the recorded data and the fitted models were calculated as a Pearson correlation coefficient. The two models used to calculate the mean myosin time attached were compared using linear correlation. Data points are presented as mean ± sem.


## 3.1. Computer Simulations

The force response of a virtual half sarcomere was simulated for a time period of two seconds for many (20,000) independent myosin heads alternately attaching and detaching to actin according to independent stochastic processes governing the random variablestatt and tdet. For any one computer simulation, random numbers representing tatt were generated to conform with a gamma distribution with values for p ranging from 1 to 48 and for β ranging from 2 to 24 such that the mean time attached was either 24 ms or 96 ms. Random numbers representing tdet were generated with values for p=2 and β=48 ms such that the mean time detached was always 96 ms. Figure 3(a) illustrates the PDFtatt(t) when p=3 and β=8. Each random number generated from this PDFtatt(t) was used to dictate the time period over which a virtual crossbridge was attached for each occasion, as indicated in Figure 3(b). Figure 3(c) illustrates the PDF for tdet, PDFtdet(t), with p=2 and β=48. Each random number generated from this PDFtdet(t) was used to dictate the time over which a virtual crossbridge was detached for each occasion. Figure 3(b) illustrates an example of successive attachment and detachment time periods.Figure 3
Schematic representation of computer-generated elastic and viscous moduli. (a)–(c) Random numbers generated from variousPDFtatt(t) (one example shown in (a)) were used to dictate the attachment times, tatt, of a virtual crossbridge. (b) shows an example attachment and detachment sequence. (c) illustrates the PDFtdet(t), with p=2 and β=48, used to dictate the intermittent times of detachment. (d) A sinusoidal length perturbation with amplitude 1 nm was used to simulate an externally applied length perturbation. This example is at 1 Hz. (e) An example series of successive crossbridge attachments and detachments would be subject to the length perturbation. (f) The length change that occurs during the time of attachment results in a force deflection. (g) The sum of these force deflections from 20,000 crossbridges is then fit to a sinusoidal (solid smooth line). (h) The magnitude of the component of the fitted sinusoid in phase with the length perturbation represents the elastic modulus at this frequency, the magnitude of the component out of phase represents the viscous modulus.Sinusoidal length perturbations of the virtual half sarcomere were simulated as having amplitude 1 nm and frequencies over the range 1–250 Hz. Figure3(d) illustrates an example perturbation at 1 Hz over a one-second time interval. Each crossbridge was assigned a stiffness constant of 1 pN/nm [28]. The change in length of the crossbridge relative to the time of initial attachment was multiplied by the stiffness constant to simulate the force resulting from the strain on the elastic element of the crossbridge. Figure 3(e) illustrates an example series of successive attachments and detachments over one second. Figure 3(f) illustrates the force deflection that would occur for a crossbridge during each time of attachment. The resulting force deflections of 20,000, independent crossbridges were summed to provide an equivalent two-second period for the ensemble in the virtual half sarcomere. An example of the resultant force from 20,000 crossbridges is given in Figure 3(g). This simulated force of the half sarcomere, Fhs(t), was fit using a simplex method to a sine function, whose amplitude (Amp) and phase (ϕ) permitted the calculation of two components that were in phase and out of phase with respect to the length perturbation. The value of each moduli was then calculated from the amplitude and phase as elastic modulus = Amp cos(ϕ) and viscous modulus = Amp sin(ϕ). The smooth line in Figure 3(g) presents an example fit of a sine function to the resultant force, and Figure 3(h) presents the in-phase and out-of-phase components provided by the Simplex method. The resulting elastic and viscous moduli, given in units of pN/nm due to the lack of normalization factors used in the simulation, were fit simultaneously to (15) using a nonweighted Levenberg-Marquardt nonlinear least-squares routine. The diagonal of the resulting covariance matrix indicated the variances of parameter estimates. The square root of the variance thus provided the standard error for each parameter estimate, which was then multiplied 1.96 to provide the 95% confidence range. All computer simulations, random number generation, curve-fittings, and parameter estimation routines were performed using IDL Version 7.0 (ITT, Boulder, CO).


## 3.2. Solutions

All reagents were purchased from Sigma (St. Louis, Mo). Solutions were formulated by solving equations describing ionic equilibria [29]. Concentrations are expressed in mmol/L unless otherwise noted. Relaxing solution: pCa 8.0, 5.0 ethylene glycol-tetra-acetic acid (EGTA), 5.0 ATP, 1.0 Mg2+, 20 hydroxyethyl-aminoethanesulfonic acid (BES), 35 phosphocreatine, 300 U/mL creatine kinase, ionic strength 200, pH 7.0. Activating solution: same as relaxing solution with pCa  4.0. Rigor solution: same as activating without added ATP, creatine kinase, or phosphocreatine. Storage solution: same as relaxing with 10 μg/mL leupeptin and 50% wt/vol glycerol. Skinning solution: same as storage with 30 mM 2,3-butanedione monoxime (BDM) and 1% wt/vol Triton X-100.


## 3.3. Viscoelastic Mechanics

All procedures were reviewed and approved by the Institutional Animal Care and Use Committees of The University of Vermont. Male wild-type mice were fed either a normal mouse diet (WT) or an iodine deficient, 0.15% 6-N-propyl-2-thiouracil (PTU) diet for at least nine weeks prior to their being killed by rapid cervical dislocation. The PTU mice therefore became hypothyroid resulting in the expression ofβ-MyHC in the myocardium in contrast to the α-MyHC expressed in the WT mouse heart.Mouse left ventricular skinned myocardial strips were prepared using methods similar to those described previously to yield thin strips (~140μm diameter, ~800 μm length) with longitudinally oriented parallel fibers [7]. These strips were chemically skinned for 2 hr at 22°C and stored at −20°C for no more than 5 days. At the time of study, aluminum T-clips were attached to the ends of a strip ~150 μm apart. The strip was mounted between a piezoelectric motor (Physik Instrumente, Auburn, MA) and a strain gauge (SensoNor, Horten, Norway), lowered into a 30 μL droplet of relaxing solution maintained at 37°C and incrementally stretched to and maintained at 2.2 μm sarcomere length detected by videography and digital Fourier transform techniques (IonOptix, Milton, MA).Strips were calcium activated at pCa 4.5 and subjected to decreasing concentrations of 5, 2, and 1 mM MgATP by exchanging equal volumes of rigor solution. Sinusoidal perturbations of amplitude 0.125% strip length were applied over the frequency range 0.125–250 Hz. The elastic and viscous moduli were calculated from the recorded tension transient as the relative magnitudes of the in-phase and out-of-phase components with respect to the imposed sinusoidal length perturbations [13, 20, 30]. The measured complex modulus was fit to (1), representing the single exponential distribution model of tatt, and to (16) below, representing the gamma distribution model of tatt, using a non weighted Levenburg-Marquardt non linear least-squares routine running within IDL Version 7.0 (ITT, Boulder, CO)(16)Y(ω)=A(iω)k-B(iω2πb+iω)+ξ{1-1piωβ((iωβ+1)p-1(iωβ+1)p)}.Mean myosin time attached was calculated as(2πc)-1 when the fit was performed with (1) and as pβ for fits using (16). The correlation coefficients between the recorded data and the fitted models were calculated as a Pearson correlation coefficient. The two models used to calculate the mean myosin time attached were compared using linear correlation. Data points are presented as mean ± sem.


## 4. Results

### 4.1. Computer Simulations

The computer simulations presented here served two purposes: first, to demonstrate the accuracy of the analytical derivation resulting in (15) and second, to demonstrate that estimating the parameters of (15) by the nonlinear least-squares methods provided reasonable estimates of the parameters known a priori to underlie the computer generated elastic and viscous moduli.Figures4(a) and 4(b) illustrate example elastic and viscous moduli, which resulted from three separate computer simulations using the parameter pairs p=1 and β=24, p=8 and β=3, and p=24 and β=4. The thin solid line represents the best fit of (15) to the simulated data. It is interesting to note that the first example, p=1 and β=24, represents the results from a single exponential model of PDFtatt(t). The elastic modulus monotonically rises to a maximum, and the viscous modulus monotonically rises to a max at 6.6 Hz and then monotonically falls. It is important to note that the shape of the viscous modulus when p=1 is symmetrical about the peak when shown on the log-scaled axis. The other examples with p>1 shown in Figures 4(a) and 4(b) do not demonstrate similarly smooth monotonic behavior. The elastic and viscous moduli instead show fluctuations, which result from the narrower distributions of tatt when p>1 compared to the single exponential distribution. Also, the viscous modulus is no longer symmetrical about its peak.Figure 4
Examples of results from computer simulations. (a, b) The elastic and viscous moduli are plotted against frequency for threep and β pairs. Nonlinear least-squares fitting of (15) to the data produced a good fit shown as the line for each pair of p and β. (c) The parameters estimated from fitting (15) to the simulated data provided a comparison between the PDFtatt(t) used to generate the data and the gamma distribution predicted from the fits. In this case, the predicted gamma distributions are nearly indistinguishable from the original PDFtatt(t). (d) The survival function, Sτ(t), was also predicted from the parameter estimates and again very closely resembles the original Sτ(t).Table1 provides a comparison of parameter values used to generate the simulations and the parameter estimates calculated by the nonlinear least-squares routine fitting (15). For the three fitted parameters p,β, and ξ, the 95% confidence intervals for the parameter estimates were narrow compared to the magnitude of the values. These confidence intervals demonstrate the high precision and uniqueness of the parameter estimates. In addition, the correlation coefficients between values used and values estimated were consistently greater than 0.997. Based on the narrow confidence intervals for parameter estimates and the high correlation between parameter estimates and actual parameter values, we are confident that the analytical expression provided by (15) accurately represents the mechanical consequences of temporarily attached myosin crossbridges, whose times of attachment are gamma distributed.Table 1
Parameter values used to generate random numbers in computer simulations of a force response due to myosin crossbridges subjected to sinusoidal length perturbations, and estimated parameter values from that force response fit to (15). For all simulations t̅det=96ms, N=20,000 and kstiff = 1 pN/nm. Est[] = estimated value of bracketed parameter, () indicates the ±95% confidence intervals for estimated value. The correlation between actual and estimated values resulted in a correlation coefficient greater than 0.997 for each parameter.
Sim numberpβ (ms)t̅att (ms)ξ (pN/nm)Est[p]Est[β] (ms)Est[t̅att] (ms)Est[ξ] (pN/nm)(1)1242440000.91 (0.01)25.13 (0.06)22.894005.0 (0.5)(2)2122440002.01 (0.01)11.98 (0.03)24.094019.8 (0.4)(3)382440003.22 (0.01)7.48 (0.02)24.083981.2 (0.4)(4)642440006.11 (0.02)3.89 (0.01)23.793986.4 (0.4)(5)832440008.70 (0.03)2.75 (0.01)23.953997.4 (0.4)(6)42496100003.99 (0.01)24.03 (0.03)95.8810007.0 (0.3)(7)81296100006.83 (0.01)13.93 (0.02)95.1210003.8 (0.3)(8)128961000012.00 (0.02)7.98 (0.01)95.7110012.8 (0.3)(9)166961000016.04 (0.03)6.05 (0.01)97.0310022.6 (0.3)(10)244961000024.10 (0.06)3.99 (0.01)96.0710017.2 (0.3)(11)482961000049.61 (0.15)1.92 (0.01)95.1310004.0 (0.3)Note:ξ=N(t̅att/t̅cycle)k̅stiff.For the three examples given in Figures4(a) and 4(b), the estimated parameter values were used to calculate the PDFtatt(t) and Sτ(t) using (12) and (13a), respectively.  Figure 4(c) illustrates a comparison between the PDFtatt(t) used to generate the simulation and that calculated from the parameter estimates. Figure 4(d) provides a similar comparison for S(τ). For both PDFtatt(t) and Sτ(t), the estimated parameters reasonably reproduced the original functions used to generate the simulations. These comparisons illustrate the robustness with which parameters of the gamma distribution can be estimated using the nonlinear least-squares methods, which fit (15) to the elastic and viscous moduli.


### 4.2. Viscoelastic Mechanics

Muscle strips from both groups activated to a similar maximum developed tension of 20.7 ± 1.8 mN·mm−2 for WT (n=6) and 22.1 ± 1.7 mN·mm−2 for PTU (n=6). Figure 5 illustrates the elastic and viscous moduli recorded for maximally activated muscle strips isolated from WT and PTU mice, representing, respectively, the α-MyHC (Figures 5(a) and 5(b)) and β-MyHC (Figures 5(c) and 5(d)) isoforms. The myosin crossbridge kinetics in the α-MyHC are faster than those in β-MyHC and are reflected in the higher frequencies ranges for the dips and shoulders in the α-MyHC. For example, the most prominent dip and shoulder of the viscous modulus at 5 mM MgATP occur at ~9 and 90 Hz, respectively, in the α-MyHC and occur at lower frequencies 2 and 30 Hz in the β-MyHC. These observations are consistent with a longer-lived tatt in the β-MyHC compared to α-MyHC [5, 9].Figure 5
Example elastic and viscous moduli from mouse cardiac muscle strips. (a, b) The elastic and viscous moduli recorded at frequencies ranging 0.125–250 Hz for wild-type mouse cardiac muscle expressingα-MyHC at maximum calcium activation and 37°C. The frequency characteristics of the dips and shoulders are reflective of the underlying enzymatic activity and shift to lower frequencies as MgATP is reduced from 5 to 1 mM. (c, d) The elastic and viscous moduli recorded for PTU-fed mouse cardiac muscle expressing β-MyHC. Frequency characteristics appear at lower frequencies compared to α-MyHC. Again, characteristics shift to lower frequencies as MgATP is reduced.As MgATP is reduced the dips and shoulders of the moduli shift to lower frequencies, as best seen in the viscous modulus. The muscle also becomes stiffer, which is reflected in higher values of the elastic modulus, due to a greater fraction of crossbridges formed at any one time. These observations with lowering MgATP concentrations are consistent with a prolongedtatt due to a prolonged rigor state prior to MgATP binding to myosin and subsequent myosin detachment from actin.


### 4.3. Parameters of Gamma Distribution Model oftatt

Fits of (16) to the elastic and viscous moduli provided estimates of the gamma distribution parameters, p and β, and a calculation of mean tatt as the product pβ. At 5 mM MgATP the shape parameter p was lower in the α-MyHC compared to β-MyHC (Figure 6(a)). At lower MgATP (1 and 2 mM), the shape parameter rises in both myosin isoforms and is not different between them. At all MgATP conditions and in both isoforms, the estimate for p was consistently near or greater than 2. This finding suggests that a single exponential distribution, which would have been represented by the gamma distribution with p=1, did not fit the data as well as a gamma distribution with p>2. Our finding p>2 suggests two or more intermediate time periods within tatt. These results imply that the PDFtatt(t) reflected in the recorded elastic and viscous moduli is better represented by a distribution that describes multiple time intervals within tatt. In the strictest interpretation of the gamma distribution shape parameter, the rise in p with decreasing MgATP implies additional elemental time periods summed together to produce the total tatt. We would, however, caution that p approaching 10 as MgATP is lowered may be more reflective of an increasing tatt and not suggestive of additional biochemical states, as the lower MgATP would prolong only one biochemical state, the rigor state [4].Figure 6
Results of fitting elastic and viscous modulus to gamma distributed and single exponential distributed models ofPDFtatt(t). (a, b) The shape parameter, p, and the scale parameter, β, were estimated from fitting the gamma distribution model of (16) to the measured elastic and viscous moduli. P*<0.05 by t-test between MyHC isoforms at same MgATP. (c) The elastic modulus recorded from one PTU β-MyHC muscle preparation at 1 mM MgATP is shown with its fit to (16), that is, the gamma distribution (GD) model. The goodness of fit is representative of the fits for all samples in this study. The A- and B-processes of (16) were then subtracted from the recorded elastic modulus to demonstrate fluctuations in the remaining C-process, which are more apparent when compared against the single exponential model (p=1) of the tatt distribution (dotted line). These fluctuations are mechanical consequences of the gamma distributed tatt. (d) The corresponding viscous modulus was also subjected to fit and subtraction of A- and B-processes. The resulting C-process is compared against the single exponential model of the tatt distribution (dotted line). The asymmetry in the peak of the viscous modulus is subtle, but apparent at frequencies higher than 11 Hz.The parameterβ effectively represents the duration of the intermediate time period that is summed p times to make up the total tatt. As illustrated in Figure 6(b), parameter β was found to be on the order of 1.5–2 ms for both α-MyHC and β-MyHC, except at 5 mM MgATP where β was 3.5 ms in the β-MyHC. We expected this parameter to have increased with decreasing MgATP, as lowering MgATP prolongs the rigor state. This parameter, however, was not as sensitive as p in describing the prolongation of tatt with decreasing MgATP.Fluctuations and asymmetries in the elastic and viscous moduli, which are characteristic of the gamma distribution forPDFtatt(t), may not be easily seen in the averaged values for the moduli presented in Figure 5. To better demonstrate these fluctuations, an example pair of elastic and viscous moduli recorded from one β-MyHC muscle preparation at 1 mM MgATP is presented in Figures 6(c) and 6(d). The A- and B-processes have also been subtracted from the measured elastic and viscous moduli resulting in the C-process, which reflects the viscoelastic mechanics modeled here. The dotted lines shown in Figures 6(c) and 6(d) represent the expected C-process when p=1, the single exponential distribution. The recorded C-process, however, demonstrates fluctuations in the elastic modulus and an asymmetry in the viscous modulus different from the p=1 case and consistent with the fluctuations and asymmetries for p>1 illustrated in Figure 4.The meantatt was prolonged in the β-MyHC compared to α-MyHC at each MgATP examined and was also prolonged as MgATP was reduced (Figure 7(a)). These findings are consistent with known consequences of myosin isoform and MgATP availability [5, 9, 16]. The calculation of the mean tatt based on the gamma distributed tatt was generally higher than that calculated based on the exponential distributed tatt, showing that the shape of the PDFtatt(t) affects the calculated value of mean tatt. Notably, the correlation between tatt based on the gamma distribution and that based on the exponential distribution model was very strong, r2=0.932 (Figure 7(b)), which indicates that both methods are capable of making comparisons of tatt among multiple groups or conditions. The Pearson correlation coefficient between the recorded data and the fitted models was very high, but was higher in the gamma distribution model (r2=0.986±0.002) compared to single exponential model (r2=0.980±0.002). We would expect the estimate of mean tatt based on the gamma distribution to be somewhat more accurate than those based on the single exponential distribution. Unfortunately, at this time and without some other independent method to verify the actual mean tatt in the muscle strip, we cannot say with certainty what bias is introduced into our estimate of tatt using either model.Mean time of crossbridge attachment. (a) The mean time of crossbridge attachment using both models ofPDFtatt(t) demonstrated a longer tatt for β-MyHC and a prolongation in tatt as MgATP was reduced. The two models agreed qualitatively in these respects, but did not agree quantitatively as estimates using the gamma distribution (GD) model were ~40% higher than those using the single exponential distribution model (SED). (b) Values for tatt using the two models correlated strongly. Despite the relative differences in the estimates of tatt, either model produces similar results for comparing the effects of isoform and MgATP. P*<0.05 by t-test between MyHC isoforms at same MgATP.
(a)(b)The parameter estimates resulting from the fits to (16) were used to predict the distributions of tatt, that is, PDFtatt(t), for each MyHC isoform and at each MgATP condition examined (Figure 8). The distributions for the α-MyHC (Figure 8(a)) appear at shorter time periods compared to those of β-MyHC (Figure 8(b)). The distributions in both isoforms shift to longer time periods as MgATP is reduced. According to results from optical trapping [2, 4], in both isoforms inorganic phosphate Pi is released either prior to force production or very quickly after force production. Thus we depict Pi release in Figures 8(c)–8(f) to occur quickly. At saturating MgATP concentrations (i.e., 5 mM), the rigor state would be very short and the longest-lived state and the most significant contributor to the total tatt would be the ADP state (Figures 8(c) and 8(e)). Because lowering MgATP concentration would lengthen only the rigor state (Figures 8(d) and 8(f)), the mean tatt in both isoforms would have been prolonged due to longer rigor states. Particularly in the cases of lower MgATP, the value for tatt would be made up of the sum of multiple times periods, as indicated by the increased shape parameter, and the distribution of tatt is better represented by the gamma distribution.Figure 8
Predicted gamma distributedPDFtatt(t) for the two myosin isoforms under varying MgATP conditions at 37°C. (a) The parameter values estimated for the α-MyHC demonstrate an increasing shape parameter and prolonged tatt as MgATP was reduced. The resulting PDFtatt(t) looks more Gaussian in shape rather than exponential as MgATP was reduced. This would be consistent with a longer rigor state prior to MgATP binding and detachment of the myosin crossbridge. (b) The distributions for β-MyHC represent longer tatt compared to α-MyHC. Again, a longer rigor state would arise with lower MgATP. (c, d) In the α-MyHC at 5 mM MgATP, the release of Pi, release of ADP, and binding of ATP occur in relatively quick succession. The time periods that represent the Pi and ADP states, however, would not be affected by MgATP availability. With the reduction in MgATP concentration, the rigor state would be prolonged as detected in the viscoelastic mechanics. (e, f). In the β-MyHC, the ADP state is longer than that in the α-MyHC and the rigor state is more dramatically prolonged with decreasing MgATP compared to α-MyHC.


## 4.1. Computer Simulations

The computer simulations presented here served two purposes: first, to demonstrate the accuracy of the analytical derivation resulting in (15) and second, to demonstrate that estimating the parameters of (15) by the nonlinear least-squares methods provided reasonable estimates of the parameters known a priori to underlie the computer generated elastic and viscous moduli.Figures4(a) and 4(b) illustrate example elastic and viscous moduli, which resulted from three separate computer simulations using the parameter pairs p=1 and β=24, p=8 and β=3, and p=24 and β=4. The thin solid line represents the best fit of (15) to the simulated data. It is interesting to note that the first example, p=1 and β=24, represents the results from a single exponential model of PDFtatt(t). The elastic modulus monotonically rises to a maximum, and the viscous modulus monotonically rises to a max at 6.6 Hz and then monotonically falls. It is important to note that the shape of the viscous modulus when p=1 is symmetrical about the peak when shown on the log-scaled axis. The other examples with p>1 shown in Figures 4(a) and 4(b) do not demonstrate similarly smooth monotonic behavior. The elastic and viscous moduli instead show fluctuations, which result from the narrower distributions of tatt when p>1 compared to the single exponential distribution. Also, the viscous modulus is no longer symmetrical about its peak.Figure 4
Examples of results from computer simulations. (a, b) The elastic and viscous moduli are plotted against frequency for threep and β pairs. Nonlinear least-squares fitting of (15) to the data produced a good fit shown as the line for each pair of p and β. (c) The parameters estimated from fitting (15) to the simulated data provided a comparison between the PDFtatt(t) used to generate the data and the gamma distribution predicted from the fits. In this case, the predicted gamma distributions are nearly indistinguishable from the original PDFtatt(t). (d) The survival function, Sτ(t), was also predicted from the parameter estimates and again very closely resembles the original Sτ(t).Table1 provides a comparison of parameter values used to generate the simulations and the parameter estimates calculated by the nonlinear least-squares routine fitting (15). For the three fitted parameters p,β, and ξ, the 95% confidence intervals for the parameter estimates were narrow compared to the magnitude of the values. These confidence intervals demonstrate the high precision and uniqueness of the parameter estimates. In addition, the correlation coefficients between values used and values estimated were consistently greater than 0.997. Based on the narrow confidence intervals for parameter estimates and the high correlation between parameter estimates and actual parameter values, we are confident that the analytical expression provided by (15) accurately represents the mechanical consequences of temporarily attached myosin crossbridges, whose times of attachment are gamma distributed.Table 1
Parameter values used to generate random numbers in computer simulations of a force response due to myosin crossbridges subjected to sinusoidal length perturbations, and estimated parameter values from that force response fit to (15). For all simulations t̅det=96ms, N=20,000 and kstiff = 1 pN/nm. Est[] = estimated value of bracketed parameter, () indicates the ±95% confidence intervals for estimated value. The correlation between actual and estimated values resulted in a correlation coefficient greater than 0.997 for each parameter.
Sim numberpβ (ms)t̅att (ms)ξ (pN/nm)Est[p]Est[β] (ms)Est[t̅att] (ms)Est[ξ] (pN/nm)(1)1242440000.91 (0.01)25.13 (0.06)22.894005.0 (0.5)(2)2122440002.01 (0.01)11.98 (0.03)24.094019.8 (0.4)(3)382440003.22 (0.01)7.48 (0.02)24.083981.2 (0.4)(4)642440006.11 (0.02)3.89 (0.01)23.793986.4 (0.4)(5)832440008.70 (0.03)2.75 (0.01)23.953997.4 (0.4)(6)42496100003.99 (0.01)24.03 (0.03)95.8810007.0 (0.3)(7)81296100006.83 (0.01)13.93 (0.02)95.1210003.8 (0.3)(8)128961000012.00 (0.02)7.98 (0.01)95.7110012.8 (0.3)(9)166961000016.04 (0.03)6.05 (0.01)97.0310022.6 (0.3)(10)244961000024.10 (0.06)3.99 (0.01)96.0710017.2 (0.3)(11)482961000049.61 (0.15)1.92 (0.01)95.1310004.0 (0.3)Note:ξ=N(t̅att/t̅cycle)k̅stiff.For the three examples given in Figures4(a) and 4(b), the estimated parameter values were used to calculate the PDFtatt(t) and Sτ(t) using (12) and (13a), respectively.  Figure 4(c) illustrates a comparison between the PDFtatt(t) used to generate the simulation and that calculated from the parameter estimates. Figure 4(d) provides a similar comparison for S(τ). For both PDFtatt(t) and Sτ(t), the estimated parameters reasonably reproduced the original functions used to generate the simulations. These comparisons illustrate the robustness with which parameters of the gamma distribution can be estimated using the nonlinear least-squares methods, which fit (15) to the elastic and viscous moduli.


## 4.2. Viscoelastic Mechanics

Muscle strips from both groups activated to a similar maximum developed tension of 20.7 ± 1.8 mN·mm−2 for WT (n=6) and 22.1 ± 1.7 mN·mm−2 for PTU (n=6). Figure 5 illustrates the elastic and viscous moduli recorded for maximally activated muscle strips isolated from WT and PTU mice, representing, respectively, the α-MyHC (Figures 5(a) and 5(b)) and β-MyHC (Figures 5(c) and 5(d)) isoforms. The myosin crossbridge kinetics in the α-MyHC are faster than those in β-MyHC and are reflected in the higher frequencies ranges for the dips and shoulders in the α-MyHC. For example, the most prominent dip and shoulder of the viscous modulus at 5 mM MgATP occur at ~9 and 90 Hz, respectively, in the α-MyHC and occur at lower frequencies 2 and 30 Hz in the β-MyHC. These observations are consistent with a longer-lived tatt in the β-MyHC compared to α-MyHC [5, 9].Figure 5
Example elastic and viscous moduli from mouse cardiac muscle strips. (a, b) The elastic and viscous moduli recorded at frequencies ranging 0.125–250 Hz for wild-type mouse cardiac muscle expressingα-MyHC at maximum calcium activation and 37°C. The frequency characteristics of the dips and shoulders are reflective of the underlying enzymatic activity and shift to lower frequencies as MgATP is reduced from 5 to 1 mM. (c, d) The elastic and viscous moduli recorded for PTU-fed mouse cardiac muscle expressing β-MyHC. Frequency characteristics appear at lower frequencies compared to α-MyHC. Again, characteristics shift to lower frequencies as MgATP is reduced.As MgATP is reduced the dips and shoulders of the moduli shift to lower frequencies, as best seen in the viscous modulus. The muscle also becomes stiffer, which is reflected in higher values of the elastic modulus, due to a greater fraction of crossbridges formed at any one time. These observations with lowering MgATP concentrations are consistent with a prolongedtatt due to a prolonged rigor state prior to MgATP binding to myosin and subsequent myosin detachment from actin.


## 4.3. Parameters of Gamma Distribution Model oftatt

Fits of (16) to the elastic and viscous moduli provided estimates of the gamma distribution parameters, p and β, and a calculation of mean tatt as the product pβ. At 5 mM MgATP the shape parameter p was lower in the α-MyHC compared to β-MyHC (Figure 6(a)). At lower MgATP (1 and 2 mM), the shape parameter rises in both myosin isoforms and is not different between them. At all MgATP conditions and in both isoforms, the estimate for p was consistently near or greater than 2. This finding suggests that a single exponential distribution, which would have been represented by the gamma distribution with p=1, did not fit the data as well as a gamma distribution with p>2. Our finding p>2 suggests two or more intermediate time periods within tatt. These results imply that the PDFtatt(t) reflected in the recorded elastic and viscous moduli is better represented by a distribution that describes multiple time intervals within tatt. In the strictest interpretation of the gamma distribution shape parameter, the rise in p with decreasing MgATP implies additional elemental time periods summed together to produce the total tatt. We would, however, caution that p approaching 10 as MgATP is lowered may be more reflective of an increasing tatt and not suggestive of additional biochemical states, as the lower MgATP would prolong only one biochemical state, the rigor state [4].Figure 6
Results of fitting elastic and viscous modulus to gamma distributed and single exponential distributed models ofPDFtatt(t). (a, b) The shape parameter, p, and the scale parameter, β, were estimated from fitting the gamma distribution model of (16) to the measured elastic and viscous moduli. P*<0.05 by t-test between MyHC isoforms at same MgATP. (c) The elastic modulus recorded from one PTU β-MyHC muscle preparation at 1 mM MgATP is shown with its fit to (16), that is, the gamma distribution (GD) model. The goodness of fit is representative of the fits for all samples in this study. The A- and B-processes of (16) were then subtracted from the recorded elastic modulus to demonstrate fluctuations in the remaining C-process, which are more apparent when compared against the single exponential model (p=1) of the tatt distribution (dotted line). These fluctuations are mechanical consequences of the gamma distributed tatt. (d) The corresponding viscous modulus was also subjected to fit and subtraction of A- and B-processes. The resulting C-process is compared against the single exponential model of the tatt distribution (dotted line). The asymmetry in the peak of the viscous modulus is subtle, but apparent at frequencies higher than 11 Hz.The parameterβ effectively represents the duration of the intermediate time period that is summed p times to make up the total tatt. As illustrated in Figure 6(b), parameter β was found to be on the order of 1.5–2 ms for both α-MyHC and β-MyHC, except at 5 mM MgATP where β was 3.5 ms in the β-MyHC. We expected this parameter to have increased with decreasing MgATP, as lowering MgATP prolongs the rigor state. This parameter, however, was not as sensitive as p in describing the prolongation of tatt with decreasing MgATP.Fluctuations and asymmetries in the elastic and viscous moduli, which are characteristic of the gamma distribution forPDFtatt(t), may not be easily seen in the averaged values for the moduli presented in Figure 5. To better demonstrate these fluctuations, an example pair of elastic and viscous moduli recorded from one β-MyHC muscle preparation at 1 mM MgATP is presented in Figures 6(c) and 6(d). The A- and B-processes have also been subtracted from the measured elastic and viscous moduli resulting in the C-process, which reflects the viscoelastic mechanics modeled here. The dotted lines shown in Figures 6(c) and 6(d) represent the expected C-process when p=1, the single exponential distribution. The recorded C-process, however, demonstrates fluctuations in the elastic modulus and an asymmetry in the viscous modulus different from the p=1 case and consistent with the fluctuations and asymmetries for p>1 illustrated in Figure 4.The meantatt was prolonged in the β-MyHC compared to α-MyHC at each MgATP examined and was also prolonged as MgATP was reduced (Figure 7(a)). These findings are consistent with known consequences of myosin isoform and MgATP availability [5, 9, 16]. The calculation of the mean tatt based on the gamma distributed tatt was generally higher than that calculated based on the exponential distributed tatt, showing that the shape of the PDFtatt(t) affects the calculated value of mean tatt. Notably, the correlation between tatt based on the gamma distribution and that based on the exponential distribution model was very strong, r2=0.932 (Figure 7(b)), which indicates that both methods are capable of making comparisons of tatt among multiple groups or conditions. The Pearson correlation coefficient between the recorded data and the fitted models was very high, but was higher in the gamma distribution model (r2=0.986±0.002) compared to single exponential model (r2=0.980±0.002). We would expect the estimate of mean tatt based on the gamma distribution to be somewhat more accurate than those based on the single exponential distribution. Unfortunately, at this time and without some other independent method to verify the actual mean tatt in the muscle strip, we cannot say with certainty what bias is introduced into our estimate of tatt using either model.Mean time of crossbridge attachment. (a) The mean time of crossbridge attachment using both models ofPDFtatt(t) demonstrated a longer tatt for β-MyHC and a prolongation in tatt as MgATP was reduced. The two models agreed qualitatively in these respects, but did not agree quantitatively as estimates using the gamma distribution (GD) model were ~40% higher than those using the single exponential distribution model (SED). (b) Values for tatt using the two models correlated strongly. Despite the relative differences in the estimates of tatt, either model produces similar results for comparing the effects of isoform and MgATP. P*<0.05 by t-test between MyHC isoforms at same MgATP.
(a)(b)The parameter estimates resulting from the fits to (16) were used to predict the distributions of tatt, that is, PDFtatt(t), for each MyHC isoform and at each MgATP condition examined (Figure 8). The distributions for the α-MyHC (Figure 8(a)) appear at shorter time periods compared to those of β-MyHC (Figure 8(b)). The distributions in both isoforms shift to longer time periods as MgATP is reduced. According to results from optical trapping [2, 4], in both isoforms inorganic phosphate Pi is released either prior to force production or very quickly after force production. Thus we depict Pi release in Figures 8(c)–8(f) to occur quickly. At saturating MgATP concentrations (i.e., 5 mM), the rigor state would be very short and the longest-lived state and the most significant contributor to the total tatt would be the ADP state (Figures 8(c) and 8(e)). Because lowering MgATP concentration would lengthen only the rigor state (Figures 8(d) and 8(f)), the mean tatt in both isoforms would have been prolonged due to longer rigor states. Particularly in the cases of lower MgATP, the value for tatt would be made up of the sum of multiple times periods, as indicated by the increased shape parameter, and the distribution of tatt is better represented by the gamma distribution.Figure 8
Predicted gamma distributedPDFtatt(t) for the two myosin isoforms under varying MgATP conditions at 37°C. (a) The parameter values estimated for the α-MyHC demonstrate an increasing shape parameter and prolonged tatt as MgATP was reduced. The resulting PDFtatt(t) looks more Gaussian in shape rather than exponential as MgATP was reduced. This would be consistent with a longer rigor state prior to MgATP binding and detachment of the myosin crossbridge. (b) The distributions for β-MyHC represent longer tatt compared to α-MyHC. Again, a longer rigor state would arise with lower MgATP. (c, d) In the α-MyHC at 5 mM MgATP, the release of Pi, release of ADP, and binding of ATP occur in relatively quick succession. The time periods that represent the Pi and ADP states, however, would not be affected by MgATP availability. With the reduction in MgATP concentration, the rigor state would be prolonged as detected in the viscoelastic mechanics. (e, f). In the β-MyHC, the ADP state is longer than that in the α-MyHC and the rigor state is more dramatically prolonged with decreasing MgATP compared to α-MyHC.


## 5. Discussion

It is generally understood that force recorded in response to a length perturbation reflects (i) the number of myosin crossbridges bound to actin in a force-generating postpower stroke confirmation and (ii) a viscoelastic response that is the consequence of the length perturbation having been applied to temporarily formed crossbridges. The viscoelastic response, however, is often overlooked as a significant component of the dynamic force response. By frequency domain analysis, however, it is clear that the viscoelastic response is substantial in magnitude, and the frequency characteristics of the viscoelastic response emerge as a consequence of crossbridge kinetics [14, 16, 19]. The macroscopic measurements of viscoelastic mechanics in striated muscle combined with the modeling results presented in this paper can be used to estimate the microscopic temporal parameters reflecting myosin crossbridge kinetics, specifically the distribution of myosin attachment times (tatt) in a muscle fiber.Our current modeling culminated in a specific mathematic representation of the mechanical transfer function that would arise from intermittently attached myosin crossbridges, (9). This representation was also provided previously [19], but without significant discussion regarding the normalized survival function, that is, the PDFτ(t) term, in the integrand. We have provided here an explanation of the survival function, which underlies PDFτ(t) and its relationship with the distribution of tatt, namely, PDFtatt(t) ((10a) and (10b)). Using the survival function and its relationship to the distribution of tatt was a very important step in demonstrating how the macroscopic viscoelastic mechanics emerge from molecular phenomena.The gamma distribution model provided a more generalized representation of the many possible distributions oftatt and in so doing provided an opportunity to discern the more subtle consequences of the tatt distribution. The shape parameter in particular permits this more generalized description of the tatt distribution. The effect of the shape parameter, which was always estimated to be near or greater than 2, indicates that the distributions of tatt that underlie that recorded elastic and viscous moduli in both isoforms and under various MgATP conditions were distinct from the single exponential distribution and with lower MgATP begin to resemble Gaussian distributions. The increasing shape parameter with decreasing MgATP is consistent with the prolongation of the rigor state of the myosin crossbridge. We interpret these results for p and β to indicate that the crossbridge tatt is best represented by two or more intermediate time periods, but we caution that a p approaching a value of 10 may simply reflect a longer tatt and not additional biochemical states. The production of additional models stipulating only 2- or 3-specific states may be warranted. Such a model may be useful in detecting which biochemical state of the crossbridge cycle is affected by an intervention such as that due to drug or posttranslational modification. We are aware, however, that increasing the number of fitted parameters may also lead to greater covariance among the parameters and less meaningful results.We believe that our results demonstrate a benefit and utility of the gamma distribution model, which allows for near-zero occurrences of very short values fortatt, as would be expected for some cases such as low MgATP concentrations and lower temperatures. The single exponential model, on the other hand, presumes that the most frequent values for tatt occur at the shortest times. Our results, particularly at lower MgATP, suggest that the single exponential model may not be sufficient to fully describe the distributions of tatt and the viscoelastic mechanical consequences as occuring in skinned muscle strips. It is notable, however, that while the gamma distribution model of tatt provides a generalized model of the mechanical consequences of the intermittent attachment of myosin to actin, we found that the single exponential model is sufficient for purposes of comparing estimated values of mean tatt based on the elastic and viscous moduli.The gamma distribution model of myosintatt leads to characteristics of elastic and viscous moduli that might not be obvious or intuitively predictable. Specifically, we did not expect the fluctuations and asymmetry in the elastic and viscous moduli that emerged from a gamma distributed tatt with p>1. As the shape parameter increases and the tatt distribution becomes more Gaussian in shape, the fluctuations and asymmetries observed in the elastic and viscous moduli in the frequency domain became more pronounced. Only through completion of the analytical derivation and computer modeling were we able to demonstrate these particular mechanical consequences.


### 5.1. Limitations

Myosin kinetics are known to demonstrate dependencies on load or strain [10, 12], which were not considered in the present work. The modeling in the present work focuses on explaining only the viscoelastic mechanics at higher frequencies and specifically ignores any strain dependencies thus providing a limited model of the macroscopic mechanical consequences of myosin kinetics. A more complete model of viscoelastic mechanics would need to predict and explain the dips in the elastic and viscous moduli like those observed at the lower frequencies. The negative viscous modulus, which corresponds to mechanical energy production by the muscle preparation, may well be one particular outcome of a strain dependency on myosin kinetics. For example, the myosin time attached may be shortened during muscle lengthening at these lower frequencies, perhaps in a Pi-dependent manner [31–33], and the number of force-generating crossbridges attached during shortening versus lengthening would therefore contribute to an observation of mechanical work production.


### 5.2. Conclusion

The importance of experimental observations of the single myosin molecule using optical techniques cannot be overstated. Without these methods, we would not likely be able to discern the nature of muscle force production at the molecular level. Discerning the performance of myosin within the context of myofilament lattice structure like that foundin vivo is also important. Mathematical and computer models of muscle performance, such as that presented in this work, are valuable in their providing the opportunity to detect a measure of myosin performance, that is, the distribution of tatt, within the context of in vivo sarcomeric structure that is not currently technologically possible by other means.


## 5.1. Limitations

Myosin kinetics are known to demonstrate dependencies on load or strain [10, 12], which were not considered in the present work. The modeling in the present work focuses on explaining only the viscoelastic mechanics at higher frequencies and specifically ignores any strain dependencies thus providing a limited model of the macroscopic mechanical consequences of myosin kinetics. A more complete model of viscoelastic mechanics would need to predict and explain the dips in the elastic and viscous moduli like those observed at the lower frequencies. The negative viscous modulus, which corresponds to mechanical energy production by the muscle preparation, may well be one particular outcome of a strain dependency on myosin kinetics. For example, the myosin time attached may be shortened during muscle lengthening at these lower frequencies, perhaps in a Pi-dependent manner [31–33], and the number of force-generating crossbridges attached during shortening versus lengthening would therefore contribute to an observation of mechanical work production.


## 5.2. Conclusion

The importance of experimental observations of the single myosin molecule using optical techniques cannot be overstated. Without these methods, we would not likely be able to discern the nature of muscle force production at the molecular level. Discerning the performance of myosin within the context of myofilament lattice structure like that foundin vivo is also important. Mathematical and computer models of muscle performance, such as that presented in this work, are valuable in their providing the opportunity to detect a measure of myosin performance, that is, the distribution of tatt, within the context of in vivo sarcomeric structure that is not currently technologically possible by other means.



---

*Source: 592343-2011-11-17.xml*