Mechano-chemical Interactions in Cardiac Sarcomere Contraction: A Computational Modeling Study

We developed a model of cardiac sarcomere contraction to study the calcium-tension relationship in cardiac muscle. Calcium mediates cardiac contraction through its interactions with troponin (Tn) and subsequently tropomyosin molecules. Experimental studies have shown that a slight increase in intracellular calcium concentration leads to a rapid increase in sarcomeric tension. Though it is widely accepted that the rapid increase is not possible without the concept of cooperativity, the mechanism is debated. We use the hypothesis that there exists a base level of cooperativity intrinsic to the thin filament that is boosted by mechanical tension, i.e. a high level of mechanical tension in the thin filament impedes the unbinding of calcium from Tn. To test these hypotheses, we developed a computational model in which a set of three parameters and inputs of calcium concentration and sarcomere length result in output tension. Tension as simulated appeared in good agreement with experimentally measured tension. Our results support the hypothesis that high tension in the thin filament impedes Tn deactivation by increasing the energy required to detach calcium from the Tn. Given this hypothesis, the model predicted that the areas with highest tension, i.e. closest to the Z-disk end of the single overlap region, show the largest concentration of active Tn’s.


Introduction
It has been widely known since the 1883 work of Ringer [1] that calcium ions (Ca 2+ ) regulate muscle contraction. Although the mechanism was unknown at the time, Ringer found that rat hearts did not contract forcefully unless Ca 2+ was added to the solution perfusing the muscle. Since Ringer's discovery the relationship between Ca 2+ concentration ([Ca 2+ ]) and tension in muscle has been widely studied. Experiments have shown a steep rise in steady state tension with relatively little increase in intracellular [Ca 2+ ]. The desire to understand the steepness of the [Ca 2+ ]-tension relationship has led to the hypothesis that cardiac muscle contraction is highly cooperative [2,3]. Cooperativity implies that a single event encourages subsequent similar events. In the case of muscle contraction, a single Ca 2+ binding invoking the generation of tension might enhance additional Ca 2+ bindings to generate more tension [4]. Various hypotheses about the mechanism of cooperativity in cardiac muscle have been proposed.
Before discussing cooperativity, we explain the structures and molecular mechanisms involved in cardiac contraction. The cardiac muscle is striated, meaning that it is composed of a repeating organization of contractile units known as sarcomeres (Fig 1A). Within the sarcomeres thick and thin filaments can interact with one another to generate tension [5,6].
The sarcomere is bounded by two opposing Z-disks with the thin filaments protruding from each ( Fig 1A). A thin filament contains a double helical arrangement of spherical actin monomers (Fig 1B), each with a myosin binding site on its surface facing the outside of the helix [7]. A thin filament is surrounded by thick filaments, each having myosin heads protruding every 14.3 nm that can bind to actin to form a cross-bridge (XB) that can generate force [8]. The thick and thin filaments range in length from~1.5 to~1.7 μm and~1.1 to 1.3 μm, respectively [9]. The distance between Z-disks is equal to the sarcomere length (SL). Opposing thin filaments can mutually overlap in the so-called double overlap region ( Fig 1A). The single overlap region is the section of the thick filament outside of the double overlap region that overlaps with the thin filament. We assume that XB's can form solely in the single overlap region [10]. As SL increases up to about 2.3 μm, the length of the double overlap region decreases, thereby lengthening the single overlap region, thus providing more opportunities for XB formation and subsequently higher force generation (Fig 1A).
Ca 2+ activates the thin filament, allowing XB's to form and force to develop (Fig 1B). Seven actin monomers, one tropomyosin (Tm) molecule, and one troponin complex (Tn) form a repeating structure called a regulatory unit (RU) along the thin filament (Fig 2A) [11]. Tm molecules are double helical molecules that follow the same helical path as the actin monomers, overlapping one another head to tail under the Tn. At rest, the Tm molecule wraps around the actin monomers in a position that blocks the myosin binding sites, thereby hindering the formation of XB's and thus contraction. Each Tn is composed of three subunits, the inhibitory subunit (TnI), the Tm binding subunit (TnT) that attaches the Tn to Tm, and the Ca 2+ binding subunit (TnC) [12] (Fig 2B). TnI binds to actin monomers at rest and essentially anchors the Tm molecule in place, blocking the binding sites on actin [13]. When the cell depolarizes, Ca 2+ enters and may bind to TnC opening a hydrophobic patch on the N-terminus of TnC. The COOH terminus of TnI exhibits a strong affinity for the NH 2 terminus opened on TnC and binds strongly to it, and the inhibitory peptide of TnI releases actin [14]. The Tm strand moves laterally, thus freeing nearby XB binding sites on the thin filament ( Fig 2C). We refer to the latter movement of Tm as the activation of the Tn (Fig 2B). It has been shown that, along with Ca 2+ binding to TnC, the formation of a strong XB is necessary to fully unblock the binding  sites on actin [15]. During relaxation the [Ca 2+ ] in the cell decreases, TnC releases the bound Ca 2+ ions, and TnI rebinds to actin moving the Tm back over the XB-binding sites preventing further force generation [16].
Commonly the relation between intracellular [Ca 2+ ] and muscle tension is described by the Hill-relation, having zero active tension at low [Ca 2+ ] and approaching a level of saturation at high [Ca 2+ ]. If a single Ca 2+ ion would bind to a single Tn complex without cooperativity, for a physiological muscle contraction, the required increase in [Ca 2+ ] would be unphysiologically large. By comparing the relatively small, physiological change in [Ca 2+ ] with the resulting large change in muscle tension, Ca 2+ binding to Tn requires an increase of cooperativity by the socalled Hill coefficient n. In experiments by Dobesh et al. [18] it is found n%7.
On the basis of conventional chemical principles, such a high level of cooperativity requires full cooperation between a number of binding sites at least equal to the Hill coefficient. In cardiac muscle, each Tn contains three binding sites for Ca 2+ . Each of the two intertwined actin strands in the thin filament contains a Tn after every sequence of seven actin monomers. The Tn's of both actin strands located at equal lengths form repetitive pairs of Tn's, or in-register Tn's. Galińska-Rakoczy et al. showed that the C-terminal section of TnI extends over the actin monomer on the adjacent strand reaching to the Tm molecule [17]. Their results suggest that A RU is a functional unit of the sarcomere that is considered to work together. It is composed of 7 actin monomers, one tropomyosin (Tm) molecule, and one troponin complex (Tn). At rest, the C-terminal section of TnI extends over the actin monomer on the adjacent strand reaching to the Tm molecule, hindering its movement [17]. That hindrance is removed upon the conformational changes resulting from Ca 2+ binding to that Tn. Note: Tn's on adjacent actin strands are located in-register with one another, but this is not displayed in the figure. upon Ca 2+ binding the Tn may influence the position of the Tm molecule on the surface of the adjacent actin strand as well. Because of this interaction, we consider some degree of cooperativity of Ca 2+ binding to in-register Tn's. Each Tn complex has 3 binding sites for Ca 2+ , with one of them having a substantially lower affinity for Ca 2+ than the others [19]. So, potentially, each pair of in-register Tn complexes contains 6 binding sites, limiting the intrinsic Hill coefficient to a maximum value of 6 if each of the binding sites exhibits an equal degree of cooperativity. The term 'intrinsic' is used to indicate that the binding is considered purely chemical, excluding effects of mechanical tension. In experiments by Sun et al. [20], in the absence of muscle tension, Ca 2+ was bound to Tn with a Hill coefficient of 3, indicating that the intrinsic cooperativity was about half of the abovementioned maximum. Even in the case of full cooperativity between in-register Tn's, i.e. with an intrinsic cooperativity coefficient equaling 6, one more binding site on a Tn in the next RU would be necessary to reach a cooperativity represented by the Hill coefficient of 7 in experiments of Dobesh et al. [18]. However, Tn's are separated along the thin filament Tn's by about 35 nm, a distance too large for conventional mechanisms of chemical cooperativity.
In many current sarcomere models, the trend is to represent cooperativity along the thin filament as the interaction between neighboring XB's and RU's. Razumova and colleagues [21] developed a stiffness-distortion model of the XB and tested three different models of nearest neighbor cooperative interactions (RU-RU, XB-XB, XB-RU). They found that while none of these interactions were able to completely account for the cooperative force development observed in experiments, each mechanism had a large impact on the [Ca 2+ ]-tension relationship. XB-XB interactions had the highest impact on peak force, RU-RU interactions had the greatest impact on the steepness of the [Ca 2+ ]-tension relationship, and XB-RU interactions had the greatest impact on Ca 2+ sensitivity. The model as developed by Rice and colleagues [22] was based on the model of Razumova and colleagues [21] including assumed nearest neighbor cooperativity. Campbell and colleagues [23] have developed a model of cooperativity in which adjacent Tm molecules overlap one another at the Tn, thereby encouraging RU activation when the nearest neighbor RU is activated. It is not clear yet on what physical principles such RU-RU, RU-XB and XB-XB interactions would work.
It has further been proposed that at longer SL's the constant volume of the myofiber lattice causes the filaments to squeeze closer together, possibly encouraging binding events by moving a myosin head within reach of a binding site on the thin filament [24,25]. Similarly, Daniel et al. propose that the cooperative realignment of binding sites occurs when the thin filament is strained by up to 2 nm [26], a lateral movement that could possibly move a binding site within reach of a myosin head. Still others believe that the end-to-end interactions of the Tm's along the thin filament are responsible for the cooperative activation [27].
In the present study we introduce a mechanism based on the mechanochemical interaction of Ca 2+ binding to Tn. Binding of Ca 2+ to Tn results in conformational changes in the Tn complex. We postulate that mechanical stretch of the thin filament facilitates the change in conformation by reducing the energetic increment required for binding of Ca 2+ to Tn. Thus, mechanical tension in the thin filament will shift the equilibrium of Ca 2+ binding to Tn towards the bound state. The latter shift in equilibrium is quantified by the use of the general physical principle that the ratio of prevalence of two states depends on the difference in energy between the two states. We think that such a mechanism shows properties of cooperativity even though it differs from a pure chemical mechanism.
On the basis of the abovementioned mechanochemical principle, we designed a computational model of cardiac sarcomere mechanics that we call the MechChem model. The Mech-Chem model is intended to predict the tension in a thin filament as a function of [Ca 2+ ] and sarcomere length. We investigated whether the mechanochemical mechanism can boost cooperativity from the intrinsic level with a Hill coefficient of 3 described by Sun et al. [20] up to the level of 7 shown by Dobesh et al. [18] in skinned muscle preparations. The energy to bind Ca 2+ to Tn is assumed to decrease linearly with thin filament tension. Currently we focus on the static behavior, implying that thin filament tension was calculated while maintaining constant sarcomere length and [Ca 2+ ]. Under static conditions we considered the chemical reactions involved to be in equilibrium.
In investigating the static behavior of cardiac muscle contraction, muscle tension as simulated with the MechChem model is compared with simulations based on conventional chemical binding of Ca 2+ to Tn with a high degree of cooperativity, resulting in the typical S-shaped Hill curve. Experimental results to compare with were obtained from skinned muscle preparations, subject to various static levels of [Ca 2+ ] and sarcomere length as reported by Dobesh et al [18]. When fitting both the MechChem model and Hill-type model to experimental findings, a [Ca 2+ ]-tension relation is found for each sarcomere length, requiring the estimation of 3 parameters. Furthermore, it is investigated which parameters values of the MechChem model could remain while changing sarcomere length.

Model design
Formation of XB's is initiated by the binding of Ca 2+ to TnC which we will refer to as the binding of Ca 2+ to an in-register pair of Tn, referred to in the equations as Tn 2 . The related chemical equilibrium reaction is represented by Parameter n indicates the coefficient related to the intrinsic cooperativity of Ca 2+ binding to Tn 2 . When Ca 2+ binds to Tn, Tm moves away from the position in which it blocks the XBbinding sites on the actin monomers. We assume that the unblocking of XB-binding sites, or activation of a Tn, automatically implies XB formation in the environment of that Tn. After formation, the XB exerts a longitudinal force F XB on the thin filament which is guided as tension along the thin filament to the Z-disk (Fig 3). Thus, toward the Z-disk, the longitudinal tension S(x) in the thin filament is increasing by a step F XB at each location where a XB has formed. Variable x represents the distance from the beginning of the single overlap region moving toward the Z-disk (Fig 3). The maximum value (x max ) of x represents the location on the thin filament where the single overlap region ends. The single overlap length depends on SL and the assumed filament lengths as shown in S1 Appendix.
Looking from the middle of the sarcomere towards the Z-disk, for tension S(x j ) just distal to XB j at position x j it holds: The symbol j max indicates the number of XB's attached to the thin filament in the single overlap region. Eq (2) is graphically elucidated in Fig 3. Tension S(x) increases with each attached XB until distance x exceeds the end of the single overlap region where no XB's can be formed. The total tension of the thin filament equals the tension acting on the Z-disk from that thin filament.

Implementation of the MechChem model
The following assumptions are key to the MechChem model implementation presented in this section: • Ca 2+ binding to Tn moves Tm azimuthally on the thin filament, unblocking the XB-binding sites on nearby actin monomers.
• Unblocking of binding sites implies XB binding and force development.
• XB's can form solely in the single overlap region.
• The total tension in the thin filament at the Z-disk is the sum of the individual XB force acting on the thin filament.
• Tension hinders the deactivation of the thin filament. The energy required for Tn to release Ca 2+ increases linearly with the mechanical tension in the thin filament.
• Binding of Ca 2+ to an in-register Tn pair (Tn 2 ) is intrinsically cooperative, characterized by the base cooperativity coefficient n (according to [20]: n = 3).
• When simulating isometric tension development, all XB's generate equal force.
Referring to Eq 2, the principle of the model is explained by the existence of discrete locations where XB's can form. In living muscle at any given time during contraction, there are countless XB's. Thus, we decided to represent muscular force generation by the average of many discrete states, allowing us to describe tensions as continuous functions of distance x along the representative thin filament. Focusing on the static conditions only, the involved chemical reactions are considered in equilibrium. According to Eq 1, we consider binding of Ca 2+ to Tn as a chemical binding to an in-register Tn pair, named Tn 2 .
Equilibrium concentrations of Ca 2+ , Tn 2 , and Tn 2 Ca n are related by: The symbol K TnCa represents the equilibrium constant, and n is the base cooperativity coefficient representing the intrinsic cooperativity of Ca 2+ binding to Tn 2 . We assume that when Ca 2+ binds to the Tn 2 , the nearby XB-binding sites are unblocked. For the proportion of Solving Eq (3) for [Tn 2 Ca n ] and substituting the result into Eq (4) renders an expression for P(x) as a function of [Ca 2+ ] and K TnCa (Eq 5).
Because in the MechChem model binding of Ca 2+ to Tn 2 is assumed to depend on mechanical tension in the thin filament, the equilibrium constant K TnCa is considered to depend on the position x along the thin filament. The energy required to detach Ca 2+ from the Tn is assumed to increase linearly with the tension S(x) in the thin filament. Thus, the equilibrium constant is multiplied by the exponential of the product of tension S(x) and a constant C s , representing the effect of the added affinity for Ca 2+ by an increase of tension. The latter relation is based on the general physical principle that in equilibrium the ratio of state concentrations is proportional to the exponential of a constant multiplied by the energy difference between both states. Analogously, the ratio of ion concentrations on both sides of a membrane depends on the voltage difference across the membrane. So, we express the dependence of the equilibrium concentration constant, K TnCa on tension S(x) by: The symbol K TnCa0 represents the equilibrium constant in the absence of tension in the thin filament. The physical dimension of constant C s is the inverse of tension. Replacing K TnCa (x) in Eq (5) with Eq (6) yields the tension-dependent expression for P(x) in Eq (7).
We assume that the density of attached XB's is proportional to P(x). Assuming that all attached XB's exert the same force during steady state isometric contraction, we find that the XB-induced force density f XB (x) acting at location x is proportional to P(x) (Eq 8).
Constant C f represents force density with full Tn activation, having the physical dimension of force per unit length along the thin filament. After elimination of P(x) by substitution of Eq 7 into Eq 8 and using the property that the derivative dS(x)/dx of tension with respect to x equals force density f XB (x), we find the following differential equation for S(x) with boundary condition S(0) = 0: Assuming n = 3 according to the findings by Sun et al. [20] for the intrinsic thin filament cooperativity, Eq (9) contains 3 independent parameters, i.e., C f , C s and K TnCa0 .
We have fitted the Hill-type model with the experimental data of Dobesh et al. [18] for comparison, expressing tension S H as a function of [Ca 2+ ]. The three parameters EC 50 , n H , and S max , represent the [Ca 2+ ] at the 50% level of maximum tension, the Hill coefficient, and maximum tension, respectively. Thus we used Simulation protocol Under normal conditions, cardiac muscle cells are enclosed by membranes that regulate the influx of Ca 2+ ions from the extracellular space through channels and pumps [28]. Submersing the muscle in detergent causes perforations in the cellular membranes, a procedure known as muscle skinning. Due to these perforations, the membrane channels and pumps no longer regulate intracellular ion concentrations. Thus, it is assumed that when a skinned muscle cell preparation is submerged in a solution containing ions, the intracellular ion concentrations are equal to that of the immersing solution. Thus, the intracellular [Ca 2+ ] can be manipulated, and the muscle will contract and generate tension in response. Data points published by Dobesh and colleagues (Fig 2A of [22], the equations of which are presented in the S1 Appendix (Eq A1-A3). The parameters C s , C f and K CaTn0 were varied so that the sum of the squared differences between the experimental data by Dobesh and the solution of the differential equation (Eq 9) was minimal. The minimization was performed on each individual curve. The Hill-type model with parameters EC 50 , S max , and n H was fitted to the same experimental data. For the MechChem and Hill-type models, the residual errors were assessed to find systematic differences between model and measurement. Additionally, the root mean squared error (RMSE) was calculated for all individual curves (Eq 11). S model and S Dobesh represent the tension generated in the model and the tension measured experimentally by Dobesh et al., respectively. Because the tension values reported by Dobesh et al. referred solely to the active tension generated, the results of the MechChem model are also reported as active tension. Additionally, we calculate the tension solely in the thin filament. Hence, the passive tension component contributed by the extracellular matrix or titin is viewed as a separate component that would be additional to the calculated active tension. The number of points is represented by j.
The RMSE compares the results of our model to the data from the experiments of Dobesh et al. [18].
In the MechChem model, tension was obtained by solving Eq 9 numerically in Matlab (MathWorks, Natick, MA) with the ode23 solver.

Results
The experimental data of Dobesh and colleagues in skinned cardiac muscle ( Table 2 and  In Table 1, the triplets of parameter values for the simulated curves are shown per sarcomere length for both the MechChem model and the fits according to the Hill-type model. The bestfit values we obtained for the Hill-type model were similar to those obtained by Dobesh et al. [18]. The differences we report are likely a result of either a different fitting routine or slight errors when the data were extracted from a published figure. On average, the root mean square  Table 2  error RSME for the MechChem model is smaller than that for the Hill-type model. The parameter C f varied by a maximum of 6.5% at different SL's, whereas K TnCa0 and C s decreased by 43.2% and 70.6%, respectively. Interestingly, for the MechChem model, parameter values K TnCa0 and C s exhibit strong interdependency. Our analyses showed that for a change in K TnCa0 by a factor of a (a>1), C s changed according to Eq 12 while exhibiting little influence on the resulting [Ca 2+ ]-tension relationship.
According to the abovementioned interdependency, it is possible that either K TnCa0 or C s changes with SL. We assume that K TnCa0 is the parameter that changes with SL because the necessary increase in K TnCa0 is much less than that in C s . Additionally, a large body of research has shown that the sensitivity of Tn to bind Ca 2+ changes with SL [18,29,30]. Therefore, we performed a best fit by varying K TnCa0 while keeping C f and C s fixed at their average values. The resulting fit is shown in Fig 5. As shown in Table 1, the fit is nearly as good as the fit by variation of all three parameters. For the Hill-type model, all parameters vary with sarcomere length, albeit that cooperativity only slightly diminishes with sarcomere length.
In the results, obtained with the MechChem model, cooperativity is apparent from the steep upslope of the [Ca 2+ ]-tension relationship (Fig 4A) [18] showed a decrease in EC 50 from 4.28 μM to 3.16 μM with SL increasing from 1.85 μm to 2.25 μm (as determined by the Hill fit in our analyses) (Fig 4B, Table 1). In the MechChem simulations covering the same SL range, EC 50 decreased from 4.18 μM to 3.17 μM when fitting the model to the [Ca 2+ ]-tension curves per SL (Fig 4A). Fig 6A and 6B show the MechChem simulation of tension and degree of activation along the thin filament in the single overlap zone, respectively, at different [Ca 2+ ]'s (3, 4, 5 μM with   ]-tension relationship with K TnCa0 changing with SL, while setting Cf = 1.93Á10 6 m -1 S and Cs = 6.5 S -1 . The unit S refers to the unit of tension as shown by Dobesh et al. [18]. The dots on the curve are plotted from the data presented in Table 2 and Fig 2A of

Discussion
The methods presented here provide a novel view of a mechanochemical mechanism of cooperativity in cardiac sarcomeres. The intrinsic cooperativity of Ca 2+ binding to an in-register pair of Tn's, i.e., cooperativity in the absence of mechanical tension, is characterized by a base cooperativity coefficient n = 3 according to experimental data reported by Sun et al. [20]. In skinned muscle experiments of Dobesh et al. [18], cooperativity is considerably higher, as characterized by a Hill coefficient of 7. We propose that the latter boost occurs by a mechanism where mechanical tension in the thin filament strengthens the chemical binding of Ca 2+ to Tn.
With each XB bound, the tension in the thin filament increases toward the Z-disk. This tension increments the energy required to deactivate a Tn, thus increasing the affinity of Ca 2+ for Tn towards the Z-disk. The model-simulated [Ca 2+ ]-tension relation appeared in good agreement (Table 1) with the experimental data from Dobesh et al. [18], thus supporting the hypothesis that an important contribution to cooperativity in cardiac muscle is delivered by the mechanical tension in the thin filament acting over a long range, i.e. the entirety of the thin filament.

Comparing measured [Ca 2+ ]-tension relation with MechChem model simulations
Considering the MechChem simulation with parameters C f and C s fixed and K TnCa0 depending on SL (Fig 5, Table 1 (Fig 4B), as shown by the residuals in Fig 4C and 4D

Proposed mechanism of cooperativity in the MechChem model
Cooperativity is assumed to have two components. The first is the intrinsic cooperativity related to binding of Ca 2+ to Tn in absence of mechanical tension. The related intrinsic cooperativity coefficient is set to 3 according to the findings by Sun et al. [20]. Secondly, the mechanochemical component boosts the cooperativity from the intrinsic coefficient of 3 up to the experimentally measured value of about 7. Sun et al. [20] inhibited tension development in trabeculae with blebbistatin and found that in the absence of tension the Ca-dependent activation of Tn's along the thin filament remained cooperative with a Hill coefficient n H = 3. The Hill coefficient indicates the minimum number of binding sites that cooperate. Each Tn contains 3 binding sites for Ca 2+ [20], but they do not all have the same affinity and are not likely to fully cooperate. Galińska-Rakoczy et al. [17] showed that the C-terminal section of TnI touches the Tm on the adjacent actin strand, thus coupling both Tn's, forming an in-register pair Tn 2 . Consequently, for each Tn 2 , there are 6 binding sites available for inherent cooperation. Apparently, cooperativity is partial, resulting in n %3.
By introduction of the MechChem model of tension-driven cooperativity, we propose that mechanical tension in the thin filament results in the strengthening of the binding between Ca 2+ and Tn, thereby hindering deactivation, implying that tension lowers the energy associated with Ca 2+ being bound to Tn. Combining intrinsic cooperativity with the mechanochemical mechanism results in the MechChem model, which results display a striking resemblance to the experimentally measured tension as a function of [Ca 2+ ] in skinned muscle preparations (Figs 4 and 5). Many correlations between change of Ca 2+ binding properties and mechanical events have been reported.
Rieck et al. [31] showed that even in the absence of tension in the thin filament, the formation of strong XB's stabilized the open conformation of the Tn. Since XB's can only form if a nearby Tn is in the unblocked state, that Tn may help to unblock the in-register Tn, thus enhancing the formation of additional XB's attached to the paired actin strand. Furthermore, Isambert et al. [32] showed that the rigidity of the thin filament decreased when Tm was in the unblocked state as compared to being in the blocked state. This finding shows that there is a mechanical coupling between the conformational changes of Tn and elastic properties of the thin filament. More recently, Desai et al. [33] directly showed that myosin binding was necessary for complete activation of the thin filament upon partial activation due to calcium binding.
The MechChem simulated best-fit [Ca 2+ ]-tension curve was determined for intrinsic cooperativity coefficients n = 1,2,3 and shown in the supplementary material (S1 Fig). The fit is excellent when n = 3, while for n = 1 or 2, the curves are not as steep as in the physiological situation. In the experiments by Sun et al., development of mechanical tension was shown not to have an effect on cooperativity, i.e. the Hill coefficient remained as low as 3. As they mentioned in their article [20], they could not exclude the possibility that the applied fluorescent probes attached to the various structures of the Tn complex may change the properties of Tn to some degree. We think that these probes may inhibit the mechanochemical enhancement of Ca 2+ affinity by tension.

Cooperativity in the MechChem model vs. the Hill model
The [Ca 2+ ]-tension relation in skinned muscle has generally been characterized by a modified Hill curve that takes on a symmetric s-shape [34]. The [Ca 2+ ]-tension relation simulated by our model is asymmetric (Fig 4A and 4C); the initial steep rise in tension generated with additional Ca 2+ decreases closer to saturation. Whereas Dobesh et al. fitted symmetric curves, they admitted that the Hill curves consistently overestimated the tension developed in the sarcomeres as the curves rounded toward saturation (Fig 4B and 4D). In an attempt to account for the asymmetry, Dobesh et al. [18] proposed that the experimental data be fit to 2 Hill coefficients that meet at EC 50 for each curve. In the Mech Chem model, the [Ca 2+ ]-tension relation is already asymmetric, so the introduction of an additional parameter is not necessary to reshape the curve.
Both the Hill-type model and the MechChem model require 3 parameters per curve. By fitting the models to the experimental data for each sarcomere length separately, in the Hill-type model all 3 parameters appeared to depend on SL, albeit that dependency of the Hill coefficient n H appeared weak (Table 1). In the MechChem model, parameters C f and C s did not depend clearly on SL, thus hinting us to keep these parameter values fixed, while only reestimating the equilibrium constant K CaTn0 , expressing Ca 2+ affinity in absence of mechanical tension. The resulting fits were nearly as good (Fig 5), suggesting that Ca 2+ affinity apparently depends somehow on SL, while the other parameters were general to all SL's.
Our model comprises two mechanisms of cooperativity that generate the [Ca 2+ ]-tension relationship, i.e. tension in Tm strengthening the bond between Ca 2+ and Tn and intrinsic cooperative activation of in-register pairs of Tn. The Hill function is based on common chemical equilibrium and provides the basis for many current models of Ca 2+ -Tn binding. For example, Rice et al. [22] utilize a modified Hill function to model the transition between the active and inactive states of Tn. The peak intracellular [Ca 2+ ] reaches 1.45 μM in the Rice model, a concentration that the Hill function mimics well. The Hill model does not provide a physical explanation for the mechanism of cooperativity in cardiac muscle, but instead utilizes a coefficient that characterizes the steepness of the [Ca 2+ ]-tension relation. Additionally, the Hill model was initially developed to understand the cooperative binding of oxygen to hemoglobin [34], where it makes physical sense that a first chemical binding of oxygen on a hemoglobin molecule will facilitate subsequent bindings because oxygen-binding sites on hemoglobin are separated by only 2.5 to 3.5 nm [35]. To reach the high physiologic level of cooperativity, several Tn's must interact. Because Tn's are separated by about 35 nm along a single actin strand within the thin filament, it is physically difficult to explain that a chemical binding can influence the binding of a different molecule that far away by conventional chemical principles.
It has been shown that the sensitivity of Tn to bind Ca 2+ changes with SL [18,29], a characteristic implemented in the MechChem model. There is currently no clear consensus regarding the mechanism behind length dependent Ca 2+ sensitivity in cardiac muscle (for review, [36]). One proposed mechanism is that the lattice spacing decreases with longer SL, moving myosin heads closer to the thin filament and rendering the thin filament more sensitive to Ca 2+ by enhancing binding [37]. However, the lattice spacing hypothesis has been questioned after it has been shown that muscle length does not necessarily correlate with myofilament spacing [30]. It has also been proposed that phosphorylation of sites on TnI by protein kinase A and protein kinase C alters the Ca 2+ sensitivity and could prove significant in the regulation of length dependent activation [14]. Conversely, Lee et al. have shown that the increase of passive tension in titin leads to increased Ca 2+ sensitivity [38]. It is possible that the change in K TnCa0 in the MechChem model is due to a combination of the abovementioned mechanisms.
The MechChem model vs. nearest neighbor cooperativity models As already indicated in the introduction, various models have been developed on the basis of nearest neighbor cooperativity, based on RU-RU, XB-XB or XB-RU interactions. These models generally result in Hill-type [Ca 2+ ]-tension relations that fit accurately to experimental data. However, the curves are slightly, but systematically different from the MechChem curves ( Fig  4). It is not clear yet if these differences are sufficiently strong to make a choice between the two model types. We predict a clear difference to be expected, yet we do not currently have the means to test the model. In a tension bearing sarcomere, the MechChem model predicts that the hindrance to Tn deactivation imposed by high tension in the thin filament causes a higher concentration of activated Tn's, bound Ca 2+ , and XB's toward the Z-disk end of the single overlap region where tension is highest (Fig 5B). With the nearest neighbor hypothesis, no preference is to be expected on the location of activated Tn's. However, although not explicitly noted, the MechChem model does take into account some of the nearest neighbor cooperativity mechanisms. The development of tension in the thin filament begets more tension development. Thus, XB's recruit more XB's. Additionally, the tension in the thin filament at position x strongly determines whether the Tn is active at point x (XB-RU cooperativity). RU-RU cooperativity is accounted for in the MechChem model through the intrinsic cooperativity coefficient n.

Comparison with long range cooperativity hypotheses
Most of the models previously discussed have included only local cooperativity mechanisms, yet there are also models and hypotheses that incorporate cooperative mechanisms acting along the entirety of the thin filament. Brandt and colleagues [39] proposed that the mutual overlap of Tm molecules under Tn causes a simultaneous unblocking or blocking of all Tm molecules along the filament. Conversely, the "cooperative realignment of binding sites", a model developed by Daniel and colleagues [26], ignores the effects of Tm but looks instead at the possible impact of strain on the thin filament. Hence, this model can predict cooperativity in tension development but not in activation. Still others attribute the cooperative effect to the constant volume property of the myofiber matrix. The stretch of sarcomeres causes the thick and thin filaments to squeeze closer together increasing the probability of XB binding [24]. Like the cooperative realignment of binding sites, the lattice spacing hypothesis can account for cooperative tension generation but not activation.
The model designed by Land & Niederer [40] represented the entire thin filament, a compilation of 26 RU's. They hypothesized that the state (position) of the Tm molecule in each RU has a corresponding free energy determined by the state of neighboring XB's and RU's. The energy term is used to compute the probability of a RU being blocked or unblocked. While the simulated results are consistent with available experimental data and insight can be gained through the Land & Niederer model, it is composed of a system of 750 ordinary differential equations. The tension-driven cooperativity model we present consists of a single ordinary differential equation and three key parameters, so we present a highly simplified model that still captures behaviors shown in available experimental data. Our model, like the Land & Niederer [40] model, includes an energy term in the computations of state P. It is a mechanical energy term related to the tension within the thin filament that increases along the thin filament from the start of the single overlap region towards the Z-disk. Our model is different from others because the mechanics of the thin filament (S(x)) directly impact the related chemistry (Ca 2+ binding to Tn) whereas most of the other models view the mechanics and chemistry separately.
Izakov et al. have hypothesized that the number of strongly bound XB's along the thin filament affects the binding affinity of Ca 2+ to Tn and have implemented this idea in a computational model [41]. Landesberg and Sideman proposed a similar model, but it was a loosely coupled model meaning that Ca 2+ was not required to remain bound for the Tn to remain active, but bound XB's were adequate to do so [42]. Our model differs because we propose that the tension in the thin filament contributes to an increment energy required to unbind Ca 2+ from Tn. This energy increases along the thin filament toward the Z-disk as tension increases. In the models discussed above, the global affinity for Ca 2+ to bind to Tn will change based on the number of strong XB's. Within our model, however, the affinities increase with tension in the thin filament, being explained by a linear increase of binding energy with tension in the thin filament.

Mechanism of relaxation as calcium concentration decreases
Due to the interaction between mechanics and chemistry within the model, neither mechanics nor chemistry fully account for the activation or subsequent deactivation of the thin filament. While we propose that high tension in the thin filament hinders deactivation, a decrease in intracellular [Ca 2+ ] will trigger relaxation. Due to the buildup of tension within the thin filament as shown in Fig 3, there is always a loose end closer to the mid-line in which little tension is developed. The hindrance to deactivation imposed by high tension does not exist in these areas, thus promoting the deactivation of the Tn's closest to the mid-line first. We propose that the areas of highest tension (closest to the z-disk) are the latest to deactivate.

Heterogeneity of cross-bridge density
Our model has provided a potentially experimentally testable hypothesis. Model results suggest directionality in the dispersion of XB's along the thin filament with a higher concentration of XB's closer to the Z-disk. Desai and colleagues [33] were recently able to fluorescently label single myosin heads and observe single bindings to the thin filament. They observed that although myosin binding activated the RU, there was no directionality in the binding of XB's. However, tension was not developed in this model because the myosin heads were not tethered to the thick filament. Additional experimental evidence is needed to test the hypothesis of XB dispersion in the loaded thin filament. The difficulty lies in developing an experimental approach that enables viewing of each individual XB binding in a skinned muscle under tension.

Model limitations
The model presented here assumes that when an area on the thin filament becomes unblocked, XB's are automatically formed. Thus, the fraction of unblocking is proportional to the XBforce developed along the strand. The assumption that XB's form automatically when RU's are activated is not physiologically accurate. We expect this to be corrected upon the explicit incorporation of the different steps of the XB cycle. However, the experimentally measured tension in skinned muscle preparations is obtained after the steady state has been reached. Hence, the experimental data is time-independent.
The current model results are limited to the [Ca 2+ ]-tension relationship. The introduction of dynamics such as the XB cycle requires additional assumptions to be made and additional unknown parameters added. However, the next logical step for the model is the implementation of the XB cycle that will be studied in the isometric twitch. The MechChem model is a mean field approximation, so it is highly simplified. Specific spatial details such as individual binding sites are not accounted for, but the simplification reduces computational cost considerably.

Conclusions
A novel mechanochemical model of tension generation by the sarcomere has been developed based on long range cooperativity imposed by mechanical tension in the thin filament and intrinsic cooperativity resulting from the interaction between the calcium-binding sites on the in-register troponin complexes. Simulated [Ca 2+ ]-tension curves resembled those obtained in steady state isometric muscle experiments. Thus, our results support the hypothesis that high tension in the thin filament impedes deactivation by increasing the energy required to detach calcium from the troponin complex. Furthermore, we found that the tension in the thin filament was relatively low toward the beginning of the single overlap region close to the mid-line of the sarcomere but increased steeply in the overlap region closer to the Z-disk. Model simulations suggest that the concentration of calcium bindings to the troponin complexes and active XB's are low at the free end of the thin filament and saturated closer to the Z-disk. Future experimental studies are needed to test the latter property, indicating the validity of our hypothesis on the cooperative effect of tension in the thin filament on force generation by the cardiac sarcomere.