Modeling the effects of thin filament near-neighbor cooperative interactions in mammalian myocardium

In the mammalian heart, each beat is characterized by a series of pressure and volume changes that underlie the ventricular filling (i.e., relaxation) and the subsequent ejection (i.e., contraction) of an appropriate stroke volume. Contraction is considered a positive cooperative process, in which the rapid development of pressure arises from the synergistic activation of the thin filament by calcium (Ca²⁺) and strong-binding myosin crossbridges (Razumova et al., 2000; Campbell et al., 2001; Lehrer and Geeves, 2014; Desai et al., 2015; Moore et al., 2016; Solis and Solaro, 2021). During relaxation, the dissociation of Ca²⁺ from troponin C (TnC) and the subsequent detachment of myosin crossbridges from actin results in the inactivation of the thin filament. This leads to an abrupt fall in pressure, which increases in ventricular compliance and augments ventricular filling. Thus, the sequential activation and inactivation of the myocardial thin filament due to cooperative mechanisms have profound implications for beat-to-beat ventricular systolic and diastolic function (Hanft et al., 2008; Solís and Solaro, 2021), particularly when considered in the context of a submaximal release of Ca²⁺ during a cardiac twitch (Bers, 2002). Evidence for the modulation of force development by cooperative processes is inferred from the observation of the nearly 10-fold variation in crossbridge kinetics as the [Ca²⁺] is increased from threshold to maximal levels in rodent myocardium (Stelzer et al., 2006a; Giles et al., 2019). The steep activation dependence of the kinetics of force development (ktr) can be explained, at least in part, by a conceptual model describing the cooperative effects associated with crossbridge binding and thin filament near-neighbor interactions (e.g., between adjacent regulatory units) that collectively act to slow crossbridge kinetics at low levels of Ca²⁺ activation (Campbell, 1997; Razumova et al., 2000; Campbell et al., 2001). While near-neighbor interactions can reasonably explain, at least in part, the cooperative development of force in mammalian myocardium, it should be emphasized that a number of mathematical models have been designed to investigate other potential determinants of thin filament activation. Each of these sophisticated models has examined the effects of various cooperative determinants, including end-to-end tropomyosin (Tm) interactions (Farman et al., 2010; Smith et al., 2003; Smith and Geeves, 2003), a crossbridge-mediated increase in the Ca²⁺-binding kinetics of TnC (Yadid et al., 2011; Kreutziger et al., 2011), and more recently, a role for thick filament activation (Caremani et al., 2022; Brunello and Fusi, 2024). For instance, in their elegant model, Smith and colleagues demonstrated that properties inherent in a flexible chain of actin-Tm-troponin provides a coherent basis for the cooperative binding of myosin to actin (Smith and Geeves, 2003). Of particular importance was their observation that the solution-based binding data could be fit without the overt need for incorporating near-neighbor cooperative factors. Although there is an assortment mathematical models, each utilizing a differing mathematical approach (e.g., inclusion or exclusion of parametric factors) and/or investigating a divergent array of cooperative determinants, an essential point has emerged. That is, these models collectively demonstrate that cooperativity is a fundamental property of mammalian myocardium, and its multifarious nature necessitates the continual investigation of its molecular and physical components.

It is well established that across mammalian species myocardial contractile kinetics are determined in part by the ventricular myosin heavy chain isoform (MyHC) expression. All mammalian species express both α- and β-MyHC isoforms in ventricular myocardium. However, the ratios of the two MyHC isoforms differ dramatically across species to account for the nearly 10-fold difference in beat frequency and contractile kinetics between small and large mammals. For example, with a resting heart rate of nearly 600 beats per minute, the murine ventricle primarily expresses the fast α-MyHC isoform, whereas large mammals (e.g., swine and human) predominately (i.e., ≥90%) express the comparatively slower β-MyHC isoform (Reiser et al., 2001; Sadayappan et al., 2009; Locher et al., 2011; Deacon et al., 2012; Walklate et al., 2016, 2021). Furthermore, myocardial twitch kinetics are also dictated by the responsiveness of the thin filament to the activating effects of Ca²⁺ and crossbridges, as mediated via local, near-neighbor (Razumova et al., 2000; Campbell et al., 2010; Kalda and Vendelin, 2020) and longer range interactions (Tanner et al., 2007, 2012; Moore et al., 2016; Dupuis et al., 2016). To examine whether molecular cooperativity underlies the species-dependent differences in contractile kinetics, we measured the activation dependence of the rate of force development in murine and porcine myocardium (Patel et al., 2023). We observed that the activation-dependent profile of ktr was significantly different between the two species. Although the rate of force redevelopment varied by a factor of three as the level of Ca²⁺ was raised from intermediate to maximal levels in porcine myocardium, near maximal ktr values were observed at very low levels of Ca²⁺ activation. Furthermore, in murine myocardium, the pCa₅₀ for ktr was significantly right-shifted compared with the pCa50 for steady-state force. That is, increases in steady-state force occurred well before increases in the rate of force development were observed. In contrast, porcine ventricular myocardium was characterized by a tighter coupling between steady-state force and ktr, such that Ca²⁺-dependent increases in steady-state force were accompanied by simultaneous increases in the kinetics of force development (Patel et al., 2023). Therefore, as an initial step to better understand the cooperative mechanisms underlying myocardial contractile dynamics in the hearts of large and small mammals (Fig. S2), we undertook a simplified mathematical approach to quantify the contribution of thin filament near-neighbor cooperative interactions. This information is critical because the relative strength of near-neighbor interactions will collectively contribute to the extent of cooperative spread along the thin filament and subsequent force development. Our mathematical analysis focused on the three types of near-neighbor interactions: (1) an XB–XB interaction (defined by parameter v), in which the binding of a crossbridge cooperatively recruits the binding of neighboring crossbridges; (2) an XB–RU interaction (defined by parameter w), in which a strongly bound crossbridge cooperatively activates a neighboring RU; and (3) an RU–RU interaction (defined by parameter u), whereby an activated RU cooperatively activates a neighboring RU. Our model contains several unique adaptations from previously published models (Razumova et al., 2000; Campbell et al., 2010; Kalda and Vendelin, 2020). These include the further extension of the cooperative parameter u to define the effect of RU–RU cooperative interactions more clearly. First, we introduce the use of the cooperative coefficients u₁ and u₂ to describe the strength of neighboring RU–RU interactions and their impacts on the transitions between the blocked and closed states. Second, we propose the use of the cooperative coefficients z₁ and z₂ to delineate how neighboring RU–RU interactions impact the transitions between the closed and open states. Third, we derived nearest neighbor interaction factors, designated as $α,α¯,β,β¯$⁠, to measure the extent to which each of the near-neighbor RU–RU, XB–XB, and XB–RU cooperative interactions influence the transition between the blocked-to-closed and closed-to-open states of the thin filament. Collectively, our model demonstrates that the near-neighbor RU–RU interaction plays a principal role in the cooperative activation of force in both murine and porcine myocardium. However, species-specific differences in the relative strengths of the RU–RU, XB–XB, and XB–RU cooperative interactions underlie the differences in myocardial contractile dynamics in the hearts of small and large mammals.

Our model combines the processes of thin filament activation (three states) and crossbridge cycling (two states) into a coupled four-state system (Fig. 1) that is adapted from previous models (Razumova et al., 2000; Campbell et al., 2010). These collections of states represent an ensemble of myosin heads, associated actins, and regulatory proteins of the thick and thin filament. The blocked state B represents the blocked state of a RU along thin filament that prevents the formation of strongly bound XBs. The closed state C stands for the closed state of a RU with the nearest XB assumed to be in either a detached or weakly bound state. State M₁ represents an open state of a RU in which nearest XB is strongly bound but is not generating force (i.e., contributes to the stiffness of the XB). The M₂ state is considered an open state of a RU in which the nearest XB is strongly bound and generating force. A RU in the M₂ state can return to the C state, where the XB head is detached or weakly bound. The M₂-to-C state transition is unidirectional. The crossbridge cycle is assumed to consume one ATP. In contrast, the other state transitions are bidirectional and do not involve ATP hydrolysis.

The model formulation is based on the following assumptions:

• (1)

  Ca²⁺ binding to TnC releases the inhibition of TnI on actin. The dissociation of Ca²⁺ from TnC reverses this process.

• (2)

  Tm transitions from the blocked to the closed state only when TnI is dissociated from actin (i.e., Ca²⁺ binding to TnC induces the transition of Tm from blocked-to-closed state).

• (3)

  Tm must be in the blocked state before TnI can rebind to actin and Ca²⁺ can dissociate from TnC.

• (4)

  Strong binding of myosin places and holds Tm in the open state.

• (5)

  Each RU regulates a single myosin-binding site.

• (6)

  Ca²⁺ binding within a RU cannot directly expose binding sites in adjacent RUs.

• (7)

  RU and XB interactions are considered along a single infinitely long thin filament (no end effects). We have not modeled interactions between multiple thick and thin filaments.

• (8)

  Sarcomere length-dependent activation is not considered.

• (9)

  States are assumed to be randomly distributed along the length of the myofilament (which is Bragg–Williams mean field approximation in statistical physics). This simplification helps avoid the requirement for probabilistic Monte–Carlo methods in computing solutions, which facilitates the use of an ordinary differential equation system, thereby reducing computation cost.

• (10)

  The total number of actin-myosin binding sites is fixed.

in which $Ca502+$ is the Ca²⁺ concentration of thin filament-binding sites at which the ratio in two above formulas equals 0.5. We only consider conditions of constant Ca²⁺ activation. By assumption 5, we also regard B as the number of XBs that are in the unbounded state and C as the number of XBs that are in the weakly bound state. M₁ and M₂ are also considered as the number of XBs that are in the strongly bound non-force–generating state and in the strongly bound force-generating state, respectively. Notice that, in some subsections below, we sometimes utilize M to indicate that an RU is in open state or a XB is in either M₁ or M₂ state.

The total crossbridge population is divided into two subpopulations of non-cycling (state B) and cycling crossbridges (i.e., states C, M₁, and M₂). State M₂ is the unique state that can generate force during isometric conditions, and hence, isometric muscle force is proportional to the number of cycling crossbridges in the M₂ state.

Each state and the combinations of states may be expressed as a fraction of R_T:

• $λB=BRT$—fraction of RUs that are blocked or fraction of XBs that are in the unbounded state.

• $λC=CRT$—fraction of RUs that are closed or fraction of XBs that are in the weakly bound state.

• $λM=M1+M2RT$—fraction of RUs that are open or fraction of XBs that are either in M₁ or M₂ state.

• $λM1=M1RT$—fraction of RUs that are in M₁ state or fraction of XBs that are in the strongly bound non-force–generating state.

• $λM2=M2RT$—fraction of RUs that are in M₂ state or fraction of XBs that are in the strongly bound force-generating state.

Due to assumption 9, the above fractions represent the probability that a given RU or XB is in any one state. For simplicity, from now on we assume that R_T = 1. Hence, $λB=B,λC=C,λM1=M1$ and $λM2=M2$⁠.

Below, we briefly present how nearest neighbor RU–RU, XB–XB, and XB–RU interactions affects the transition rates k_BC, k_CB, $kCM1$ and $kM1C$ and the subsequent combination of all three nearest neighbor interactions into our thin filament model. For detailed presentation, we refer readers to Data S1.

The coefficient u₂ is proportional to the activation energy needed to be overcome for the transition from the B state to the C state, and while the coefficient u₁ is proportional to the activation energy needed to be overcome for the transition from the C state to the B state. Because the effects of Ca²⁺ on k_BC and k_CB are independent of the effects of neighbor interactions, $kBC(1,1)$ and $kCB(1,1)$ incorporate the Ca²⁺ effect as given by Eqs. 5 and 6.

The coefficient z₂ is proportional to the activation energy needed to be overcome for the transition from the C state to the M₁ state and z₁ from the M₁ state to the C state.

Note that the parameter v is closely related to the activation energy required for the transition between the closed and open states.

Notice that the parameter w is closely related to the activation energy required for the transition between the blocked and closed states.

For convenience, we designate the parameters $(α,α¯,β,β¯)$ as nearest neighbor interaction factors. Model state variables and equations (Table S2) and transitions rate coefficients, near-neighbor cooperative coefficients, and nearest neighbor interaction factors (Table S3) are summarized in the Data S1.

A detailed description of the mathematical model and materials and methods used in the present study can be found in the Fig. S3; and Tables S4 and S5.

Fig. S1 shows the steady-state tension development in mammalian-permeabilized myocardium. Fig. S2 shows the Ca²⁺- and activation-dependence of the rate of force redevelopment in mammalian-skinned myocardium. Fig. S3 shows the model resolution matrices generated by fitted results for modeling the Ca²⁺ dependence of the rate of force redevelopment in mammalian myocardium. Fig. S4 shows modeling the effect of near-neighbor cooperative interactions on the force–pCa relationship in mammalian myocardium. Table S1 shows the summary of steady-state mechanical data in murine- and porcine-permeabilized ventricular myocardium. Table S2 shows the summary of state variables and model equations. Table S3 shows the summary of rate coefficients, cooperative coefficients, and nearest neighbor interaction factors. Table S4 shows the fitted parameters for modeling force–pCa relationships in murine- and porcine-permeabilized ventricular myocardium. Table S5 shows the parameter identification for model resolution matrices. Data S1 shows the steady-state contractile measurements.

The kinetics of myofilament activation and crossbridge cycling were described by a four-state model, utilizing a system of three ordinary differential equations (Eqs. 2, 3, and 4) under a conservative constraint that the total number of actin-myosin sites along the whole thin filament is fixed. The effect of Ca²⁺ to activate RUs is represented by the dependence of the rate constants k_BC and k_CB on Ca²⁺ given by Eqs. 5 and 6. Without any nearest neighbor interaction, $kBC=kBC(1,1)=fBC0$⁠, $kCB=kCB(1,1)=fCB0$⁠, $kCM1=kCM1(1,1)=fCM10$⁠, and $kM1C=kM1C(1,1)=fM1C0$⁠.

the $Ca2+/Ca502+$ ratio, 6 cooperative parameters (u₁, u₂, z₁, z₂, v, and w) measuring the strength of RU–RU (i.e., u₁, u₂, z₁, and z₂), XB–XB (i.e., v), and XB–RU (i.e., w) near-neighbor interactions, and four nearest neighbor interaction factors $(α,α¯,β,and β¯)$ measuring the extent to which each near-neighbor interaction effects the forward and reverse blocked-to-closed and closed-to-open transitions. We integrated the model (Eqs. 2, 3, and 4) using the fourth-order Runge–Kutta method to predict steady-state force during varying levels of constant Ca²⁺ activation (represented by changing the pCa = −log[Ca²⁺]) and the time course of force redevelopment (ktr) starting from zero-force initial conditions.

We designed eight discrete parameter sets using assorted values of nearest neighbor interaction factors and cooperative coefficients $(α,α¯,β,β¯,u1,u2,z1,z2,v,and w)$ to quantify the unitary, tandem, and ensemble effects of the RU–RU, XB–XB, and XB–RU cooperative interactions (Table 1).

We consistently observed that the use of those parameter sets that excluded the effects of the RU–RU cooperative interaction failed to adequately fit our in vitro contractility data. This observation is illustrated by the failure of parameter set 1 (i.e., the elimination of RU–RU, XB–XB, and XB–RU cooperative interactions) to adequately fit the force–pCa (Fig. S4 A) and ktr–pCa (Fig. 5) relationships in murine and porcine myocardium. In a similar manner, parameter set 3 (i.e., unitary effects of XB–XB cooperative interactions), parameter set 4 (i.e., unitary effects of XB–RU cooperative interactions), and parameter set 6 (i.e., the tandem effects of XB–XB and XB–RU cooperative interactions) failed to fit the force–pCa and ktr–pCa relationships (data not shown). In contrast, our in vitro contractility data obtained in murine- and porcine-permeabilized myocardium was well fit using parameter sets that included the RU–RU cooperative interaction, i.e., the unitary effects of RU–RU interaction (set 2), the tandem effects of RU–RU and XB–XB (set 5) and RU–RU and XB–RU interactions (set 7), and the ensemble effects of RU–RU, XB–XB, and RU–XB interactions (set 8). For simplicity, here we present the modeling results using parameter sets 2 and 8. To examine the unitary effect of the RU–RU cooperative interaction (i.e., parameter set 2), we fit our in vitro contractility data by allowing the RU–RU cooperative coefficients describing the blocked-closed transition (u₁ and u₂) and the closed-open transition (z₁ and z₂) to become more cooperative (i.e., attain values >1) and eliminated the cooperativity associated with near-neighbor XB–XB (i.e., v = 1) and XB–RU (i.e., w = 1) interactions. Further, we fixed the nearest neighbor interaction factors (i.e., $α=α¯=β=β¯=1$⁠) to ensure that the RU–RU cooperative interaction was the sole component effecting the blocked-closed and closed-open transitions. The fitted parameters derived from the unitary RU–RU cooperative interaction are summarized in Table 2, and the in silico–derived ktr–pCa and rel ktr–pCa relationships are illustrated in Fig. 6. To provide a more physiologic framework for thin filament near-neighbor communication, we examined the ensemble effects of the RU–RU, XB–XB, and XB–RU cooperative interactions. To do so, we fit the in vitro contractility data using parameter set 8 in which the accompanying RU–RU, XB–XB, and XB–RU cooperative coefficients (i.e., u₁, u₂, z₁, z₂, v, and w) were permitted to attain values >1 (i.e., introducing cooperativity into these near-neighbor interactions). To examine the extent to which the RU–RU and XB–RU interactions effected the blocked-closed transition and the RU–RU and XB–XB interactions effected the closed-open transition, we let the nearest neighbor interaction factors vary between 0 and 1 (i.e., 0 < α < 1, $0<α¯<1$⁠, 0 < β < 1, and $0<β¯<1$⁠). The fitted parameters derived from the ensemble effects of the RU–RU, XB–XB, and XB–RU interactions are summarized in Table 2, and the in silico–derived ktr–pCa and rel ktr–pCa relationships are illustrated in Fig. 7. The use of parameter set 2 and parameter set 8 produced remarkably similar fits of the in vitro contractility data. Both parameter sets produced large values for the cooperative coefficient u₂ in murine and porcine myocardium (i.e., 2.2659 versus 17.0820; Table 2). A summary of the fraction of crossbridges in the non-cycling and cycling populations for parameter sets 1, 2, and 8 are presented in Table 3. Collectively, these results demonstrate that the interaction between near-neighbor RUs is the predominant cooperative mechanism for the spread of thin filament activation in mammalian myocardium as previously predicted (Campbell et al., 2010; Kalda and Vendelin, 2020). Furthermore, the larger fitted value of u₂ in porcine myocardium indicates a relatively stronger RU–RU cooperative coupling exists within the porcine thin filament.

The regulation of contraction in mammalian cardiac muscle is a multivariable process triggered by the initial binding of Ca²⁺ to TnC, which in turn induces a switch-like effect on thin filament regulatory proteins to permit crossbridge binding to actin and the subsequent development of force (Solis and Solaro, 2021). It is important to stress that Ca²⁺-activated force is not a simple linear function of the myoplasmic [Ca²⁺]. Rather, the effect of Ca²⁺ binding to TnC is amplified by positive cooperativity in myosin crossbridge binding to actin, in which the initial binding of Ca²⁺ to TnC increases the likelihood of subsequent crossbridge binding to actin (Moss et al., 2004, 2015). These effects are mediated via local interactions along a segment of the thin filament within and between troponin-Tm regulatory units (Moore et al., 2016; Solís and Solaro, 2021). As a result, once a myosin crossbridge binds, the activation profile expands away from the central Ca²⁺-bound TnC, permitting the cooperative recruitment of additional crossbridge binding within the activated RU and among neighboring RUs (Moss et al., 2015). However, the relative strength of near-neighbor interactions (i.e., the degree of cooperativity between near-neighbors) will affect, at least in part, the extent of cooperative spread along the thin filament due to initial Ca²⁺ and crossbridge binding (Razumova et al., 2000; Campbell et al., 2010; Kalda and Vendelin, 2020). In murine ventricular myocardium, the rate constant of force redevelopment (ktr) exhibits a steep activation dependence, increasing nearly 10-fold across varying levels of Ca²⁺ activation (Stelzer et al., 2006a; Giles et al., 2019). In contrast, the activation dependence of ktr in porcine ventricular myocardium is significantly reduced, as manifested by the near-maximal rates of force redevelopment at low levels of Ca²⁺ activation. A critical question is what molecular mechanism(s) can account, at least in part, for the differing activation-dependent profiles of ktr in the hearts of small and large mammals?

At all levels of submaximal [Ca²⁺], porcine myocardium exhibited significantly faster relative rates of force development than murine myocardium, with the greatest difference in crossbridge kinetics observed at very low [Ca^2+]. One possible way to explain the markedly divergent ktr values between porcine and murine myocardium at very low levels of Ca²⁺ activation is the binding of Ca²⁺ to a few discrete regions of the thin filament and the isolation of active RUs with Ca²⁺ bound from adjacent near-neighbor RUs having no Ca²⁺ bound (i.e., RUs in “blocked” state). The first condition may arise if the binding of Ca²⁺ to the thin filament is not uniform along the thin filament. This would result in the nonuniform (i.e., stochastic) activation of the thin filament (Risi et al., 2021), which would limit the number of strongly bound crossbridges and account for the small forces developed in both species. The second condition may result from species-dependent differences in the activation energy required to cooperatively recruit near-neighbor RUs. In porcine myocardium, the relatively higher activation energy required to cooperatively recruit near-neighbor RUs would effectively isolate and perhaps even eliminate near-neighbor RU communication at very low [Ca^2+]. As a result, the small number of activated RUs in the porcine thin filament would have very few or no adjacent near-neighbor RUs with either Ca²⁺ or a crossbridge bound, thereby making the cooperative recruitment of crossbridges into the adjacent near-neighbor RUs less likely. If the length of the functional (i.e., cooperative) unit is essentially the length of a structural unit spanning seven actin monomers (Gillis et al., 2007), these effects would appreciably increase ktr.

At intermediate [Ca^2+], relative ktr values measured in porcine myocardium were less than that observed at low and maximal [Ca^2+] but remained significantly faster than murine ktr values measured at similar levels of Ca²⁺. Once again, if we assume that Ca²⁺ binding occurs within randomly distributed RUs along the thin filament (Risi et al., 2021), then some regions of the thin filament will have a greater collection of activated RUs and other regions will have fewer numbers of activated RUs. Thus, the probability of activated RUs being in the near vicinity of other RUs in the closed or open state will be increased. As the [Ca^2+] rises to intermediate levels, the effective size of the functional group appears to increase as more activated RUs are incorporated. The central region of an activated functional group will have the greatest amounts of Ca²⁺ and crossbridges bound. Because of this, the rate of force redevelopment will be progressively slowed due to cooperative recruitment of crossbridges into the end regions of the functional group or adjacent functional groups.

As the [Ca^2+] continues to rise toward maximal levels, the number of activated RUs with Ca²⁺ bound increases in both murine and porcine myocardium. Thus, fewer inactivated RUs remain, and the concomitant impact of near-neighbor cooperative interactions will be reduced. As a result, the time course of force redevelopment in porcine and murine myocardium becomes progressively faster, and ktr converges to the respective values obtained during maximal activation.

Here, we propose that the relative contribution of the ensemble near-neighbor RU–RU, XB–XB, and XB–RU cooperative interactions vary in a species-specific manner to modulate the activation dependence of the rate of force redevelopment. As a result of the model, two RU–RU cooperative coefficients (i.e., u₂ and z₂) were derived and shown to play a critical role in fitting the in vitro contractility data obtained in porcine and murine ventricular myocardium. Collectively, these cooperative coefficients demonstrate the dynamic near-neighbor cooperative behavior of mammalian myocardium. The cooperative coefficients u₂ and z₂ measure the strength of the near-neighbor RU–RU cooperative interaction, i.e., the degree of cooperative coupling between neighboring RUs. The cooperative coefficient u₂ is closely related to the activation energy needed to be overcome for a central RU to make the transition from the blocked to the closed state. When u₂ = 1, there is no cooperative RU–RU interaction. However, as the value of u₂ increases (i.e., u₂ > 1), a greater activation energy is encountered for the transition of a central RU from the blocked to the closed state (Kalda and Vendelin, 2020). In both the murine and porcine thin filament, a central RU located between two RUs in the blocked state (BB, Fig. 8) will have a lower success frequency of transitioning compared with a central RU located between two RUs in the closed state (CC, Fig. 8). But the greater value of u₂ modeled for porcine myocardium compared with murine myocardium (e.g., u₂ = 17.082 versus 2.265; Table 2) suggests that a higher activation energy exists between near-neighbor RUs, making it harder for these neighboring RUs in the porcine thin filament to successfully transition to the closed state. It has been proposed that the introduction of stronger cooperative coupling between near-neighboring RUs should effectively slow the time course for force redevelopment (Campbell, 1997; Razumova et al., 1999). In both the murine and porcine thin filament, the initial Ca²⁺ binding to TnC would activate a RU, causing it to transition from the blocked to the closed state. At very low Ca²⁺, this activated RU would be located among neighboring RUs positioned in the blocked state. The relative strength of near-neighbor RU–RU interactions will initially define the extent of cooperative spread along the thin filament due to Ca²⁺ and subsequent crossbridge binding within the activated RU (Moss et al., 2004, 2015). As such, stronger near-neighbor RU–RU interactions, as modeled for the porcine thin filament, will tend to hold adjacent RUs in the blocked state, whereas weaker RU–RU interactions, as in the murine thin filament, would permit adjacent RUs to transition to the closed state. Since the strength of the RU–RU cooperative interactions is related to activation energy (Razumova et al., 2000; Kalda and Vendelin, 2020), the higher the value of u₂ the greater the activation energy barrier to cooperatively recruit neighboring RUs to the closed state. Thus, the neighboring RUs in the porcine thin filament will exhibit a lower success frequency of transitioning to the closed state. If the spread cooperative activation is limited to a shorter segment of the thin filament (e.g., within a single RU versus within two to three RUs) at very low [Ca^2+] (Gillis et al., 2007), then the effect of a strong cooperative coupling to slow the time course of force redevelopment will be minimized, resulting in rates of force redevelopment that approach values seen at maximal [Ca^2+] (Fig. 6 B and Fig. 7 B). In a similar manner to that proposed for u₂, the cooperative parameter z₂ is closely related to the activation energy needed to be overcome for a central RU to make the transition from the closed to the open state. In both the murine and porcine thin filament, a central RU located between two RUs in the closed state (CC, Fig. 8) will have a lower success frequency of transitioning to the open state compared with a central RU located between two RUs in the open state (MM, Fig. 8). As the value of z₂ increases, the degree of cooperative coupling between neighboring RUs becomes greater. In our model, the greater value of z₂ derived for murine myocardium compared with porcine myocardium (e.g., z₂ = 1.0118 versus 1.000; Table 2) suggests that a slightly stronger cooperative interaction governs the transition of a central RU from the closed to the open state in the murine thin filament. This suggests that the murine myocardial thin filament requires a slightly greater number of strongly bound crossbridges for activation, making it the more cooperative of the two muscles in this regard. Further elucidation of the cooperative mechanisms underlying the activation dependence of ktr can be ascertained by examining the model-derived nearest neighbor interaction factors $(α,α¯,β,and β¯)$⁠. The extent to which the RU–RU and XB–RU interactions effect the B-to-C transition (i.e., k_BC) can be described by α and 1 − α, while the C-to-B transition (i.e., k_CB) can be described by $α¯$ and $1−α¯$⁠. In a similar fashion, the effects of the RU–RU interaction and XB–XB interaction on the C-to-M₁ transition (i.e., $kCM1$⁠) can be described by β and 1 − β, while the M₁-to-C transition (i.e., $kM1C$⁠) can be described by $β¯$ and $1−β¯$⁠. The nearest neighbor interaction factors for porcine myocardium $(α,α¯,β,and β¯)=(1.0,0.8296,1.0,and 1.0)$ indicate that the process of force redevelopment is predominantly modulated by RU–RU cooperative interactions with only a slight involvement of XB–RU interactions. However, the nearest neighbor interaction factors for murine myocardium $(α,α¯,β,and β¯)=(0.9996,0.8897,0.9466,and 0.1376)$ signify that while force redevelopment is largely modulated by the RU–RU interactions, there is also additional significant involvement of the near-neighbor XB–RU and XB–XB cooperative interactions.

Our modeling data suggest that the responsiveness of the myocardial thin filament to the activating effects of Ca²⁺ and myosin crossbridges varies in a species-dependent manner. To further investigate the role of the RU–RU cooperative interactions on the contractile dynamics in murine and porcine myocardium, we modeled the effects of independently varying either u₂ or z₂ on the activation dependence of ktr (Fig. 9 and Table 4). By using Eqs. 15 and 16, we can define an “activation factor” as K, where K = k_BC/k_CB. Therefore, K becomes a function of the Ca²⁺ concentration, the RU–RU and XB–RU cooperative coefficients (u₁, u₂, and w), the RU–RU and XB–RU nearest neighbor interaction factors of the transitions between the B and C state (α, 1 − α, $α¯$⁠, $1−α¯$⁠), and the state variables B, M₁, and M₂, i.e., K = K(Ca²⁺, u₁, u₂, w, α, 1 − α, $α¯$⁠, $1−α¯$⁠, B, M₁, and M₂). Changes in K reflect the sensitivity of myocardial thin filament to the activating effects of Ca²⁺ binding when any of the cooperative coefficients (u₁, u₂, and w) are altered. By using Eqs. 17 and 18, we can define a “cross-bridge recruitment factor” as N, where N = k_CM1/k_M1C. Thus, N is a function of the RU–RU and XB–XB cooperative coefficients (z₁, z₂, and v), the RU–RU and XB–XB nearest neighbor interaction factors of the transitions between the C and M₁ state (β, 1 − β, $β¯$⁠, $1−β¯$⁠), and the state variables B, M₁, and M₂, i.e., N = N(z₁, z₂, v, β, 1 − β, $β¯$⁠, $1−β¯$⁠, B, M₁, and M₂). Changes in N indicate the sensitivity of the myocardial thin filament to the activating effects of crossbridge binding when any of the cooperative coefficients (z₁, z₂, and v) are altered.

We modeled the effect of varying u₂ on the blocked-to-closed RU transition using parameter set 8 with the cooperative coefficients u₁, z₁, z₂, v, and w held constant. In murine myocardium, we examined the effects on the activation dependence of ktr by progressively strengthening the degree of cooperative coupling between near-neighbor RUs by increasing u₂ (Fig. 9 A and Table 4). Increasing u₂ from 2.2659 to 4.0 significantly steepened the activation dependence of ktr. The nearly 24-fold increase in ktr from low to maximal levels of Ca²⁺ activation was due to the positive cooperative feedback between M₂ and the activation factor K. By Eqs. 15 and 16, since α = 0.9996 and $α¯=0.8897$⁠, the relatively small increase in u₂ increased the value of K, implying that the forward transition rate k_BC increases faster than the reverse transition rate k_CB. This results in a larger number of crossbridges available for cooperative recruitment into cycling population, which increases M₂ at all levels of Ca²⁺ activation. Because z₁ and z₂ (i.e., 1.0047 and 1.0118) are close to 1.0, the value of N remains constant as u₂ is increased. In porcine myocardium, we modeled the effects of decreasing u₂ from 17.082 to 10.0 (Fig. 9 B and Table 4), which results in a progressive weakening of the cooperative coupling between neighboring RUs. Because z₁ = 1.0001 and z₂ = 1.0, the value of crossbridge recruitment factor N remains constant (N = 0.1315) as u₂ is decreased. Since u₁ = 1.1055 (which is close to 1.0), the gradual decrease in u₂ reduces the activation factor K, which decreases the value of M₂, leading to a further decrease in K. This suggests there is negative feedback between the activation factor K and the force-generating crossbridge population M₂, when u₂ is reduced. As a result, the rate of force redevelopment at low and intermediate levels of Ca²⁺ activation increases, thereby reducing the activation dependence of ktr. The modeling result for porcine myocardium is intriguing when considered in the context of hypertrophic cardiomyopathy (HCM)-related mutations in Tm (Loong et al., 2012; Moore et al., 2016; Matyushenko et al., 2017, 2023; Tsaturyan et al., 2022; Berry et al., 2024) and troponin T (TnT) (Sommese et al., 2013; Wang et al., 2018; Landim-Vieira et al., 2023). In the human myocardial thin filament, structural changes induced by HCM-related mutations in either TnT or Tm could have a substantial impact on the responsiveness of the thin filament to the activating effects of Ca²⁺ and crossbridge binding. If the TnT or Tm HCM mutations reduce the cooperative coupling between neighboring RUs (i.e., decreasing the activation energy), the level of thin filament activation would be increased, resulting in a concomitant reduction in the activation dependence of ktr (Table 4).

To investigate the species specificity of the closed-to-open RU transition, we examined the effects of increasing z₂ on the activation dependence of ktr. To do so, we used the ensemble parameter set 8 with all cooperative coefficients held constant and allowed z₂ to increase. In both murine (Fig. 9 C) and porcine myocardium (Fig. 9 D), increasing z₂ from ∼1.0 to 2.0 resulted in a significant decrease in ktr at low and intermediate levels of Ca²⁺ activation. This effect arises due to the positive feedback between M₂ state and the crossbridge recruitment factor N as seen in Eqs. 17 and 18. Likewise, with the RU–RU and XB–RU cooperative coefficients (i.e., u₁, u₂, and w) held constant, Eqs. 15 and 16 imply that a slight increase in M₂ will increase the activation factor K, resulting in a net movement of crossbridges from the non-cycling (B) to the cycling (C, M₁, and M₂) population. This will effectively reduce ktr as the [Ca²⁺] increases from threshold to maximal levels.

In summary, our mathematical model accurately fit the in vitro contractility data measured in murine- and porcine-permeabilized myocardium. Our model quantified the effects of near-neighbor interaction on the cooperative development of force, thereby allowing for elucidation of the contribution of the RU–RU, XB–XB, and XB–RU near-neighbor interactions in mammalian myocardium. Our modeling result demonstrated that in both murine and porcine thin filament, the RU–RU cooperative interaction was the predominant near-neighbor effect, consistent with that proposed previously (Campbell et al., 2010). From an energetic perspective, the near-neighbor interaction requiring the largest free energy change would be predicted to be the rate-limiting cooperative factor. Our results agree with the model of Kalda and Vendelin (2020), who proposed that the changes in free energy associated with the RU–RU interaction are approximately fivefold larger than that for the XB–XB interaction. Furthermore, our model revealed species-specific differences in the ensemble effects of the combined RU–RU, XB–XB, and XB–RU cooperative interactions. The species-specific differential aggregate contribution of near-neighbor cooperative mechanisms underlies the differences in myocardial contractile dynamics and confers a physiologic adaptation to synchronize ventricular contractility with intrinsic beat frequency and twitch kinetics.

It is well established that cooperative mechanisms modulate both force development and force relaxation in mammalian cardiac muscle. As such, molecular cooperativity is not described by a single factor but rather is defined by a wide array of mechanisms operating from the molecular level to more globally across the sarcomere. Given the multivariant nature of cooperativity, we chose to simplify our modeling efforts to initially examine near-neighbor interactions (i.e., RU–RU, RU–XB, and XB–XB) to establish a baseline understanding of how these cooperative interactions modulate force development. Our sensitivity analysis, represented by resolution matrices, indicated that while our model is sensitive to both forward and reverse cooperative parameters, there is a slightly higher sensitivity to the forward parameters under our experimental conditions.

Furthermore, we readily acknowledge that additional forms of cooperativity most likely play significant roles in modulating cardiac contraction and relaxation. For example, our model did not examine cooperative mechanisms that exert sarcomere-wide effects (Dupuis et al., 2016; Mijailovich et al., 2021), including myofilament compliance and lattice geometry (Tanner et al., 2007, 2012), as well as communication between adjoining thin filaments (Risi et al., 2021). Furthermore, our model did not explicitly examine the contributions of the thin filament regulatory proteins Tm (Smith and Geeves, 2003; Campbell et al., 2010; Geeves et al., 2011; Moore et al., 2016; McConnell et al., 2017) and troponin (Farman et al., 2010; Kreutziger et al., 2011; Land and Niederer, 2015; Landim-Vieira et al., 2023) or the effects of phosphorylation of the myosin regulatory light chain (Toepfer et al., 2013; Karabina et al., 2015; Turner et al., 2023) and myosin-binding protein-C (Moss et al., 2015; Previs et al., 2016; Risi et al., 2022; Wong et al., 2024). The incorporation of these mechanisms, particularly the phosphorylation of cardiac-binding protein-C (cMyBP-C), is likely to significantly impact model outcomes and provide a more comprehensive understanding of the regulation of ventricular function. Cardiac-binding protein-C phosphorylation increases myocardial contractility (Stelzer et al., 2006b; Sadayappan et al., 2011; Gresham et al., 2017) at least in part, by promoting actin–MyBP-C interactions and augmenting the level of thin filament activation (Mun et al., 2014; Kampourakis et al., 2014; Risi et al., 2022; Wong et al., 2024), thereby potentially interacting with the near-neighbor cooperative interactions that we have modeled. While computational constraints necessitated our initial focus on near-neighbor cooperative interactions, we are actively working to refine our model. We recognize that the model-derived cooperative parameters provide a qualitative interpretation of the factors underlying myocardial contractile kinetics. Subsequent experimentation will be required to validate the quantitative contribution of these factors. To accomplish this, our immediate plans include incorporating relaxation transients into our experimental model, which will allow us to better explore mechanisms underlying cooperative inactivation of the thin filament. Subsequently, we aim to systematically integrate additional cooperative mechanisms, prioritizing those that our preliminary analyses suggest may have the most significant impact on force development and relaxation kinetics. Working with the framework of near-neighbor cooperative interactions, we will continually refine our model to incorporate additional cooperative mechanisms, particularly given how these mechanisms may be altered in disease states, such as cardiomyopathy (Sewanan et al., 2016; Farman et al., 2018; Mertens et al., 2024). We aim to enhance its predictive power and provide more comprehensive insights into cardiac function in health and disease.

All in vitro contractility data underlying this study are available in the published article and its online supplemental material. The MATLAB code for our model has been deposited in Github: (https://github.com/tphan86/TF_model).

Henk L. Granzier served as editor.

Research reported in this publication was supported by the National Institute of General Medical Sciences of the National Institutes of Health (NIH) under Award Number P20GM104420. The content is solely the responsibility of the authors and does not necessarily represent the official views of the NIH.

Author contributions: T.A. Phan: formal analysis, methodology, software, validation, visualization, and writing—original draft, review, and editing. D.P. Fitzsimons: conceptualization, investigation, project administration, resources, supervision, validation, visualization, and writing—original draft, review, and editing.