Introduction

Sarcomeres are the basic contractile cytoskeletal units of striated muscle, including cardiac muscle. In most species, sarcomeres are ∼2 µm in length and consist of highly ordered bipolar bundles of myosin motor proteins interdigitated with actin filaments, tied together in the Z-bands. (Fig. 1A). Force generation by myosin interacting with actin is regulated by intracellular Ca2+ and fueled by chemical energy (ATP)1. Sarcomeres are serially connected into myofibrils. A cardiomyocyte (CM) contains tens of parallel myofibrils. Overall CM contractility is a mesoscopic process that emerges from a very large number of coupled, non-linear, stochastic, elementary force generators. A fundamental question that has not been answered during more than a century of muscle research is how molecular-scale stochastic dynamics in this non-equilibrium, excitable system eventually results in the emergent, organized dynamics at the cell and organ scale.

Sarcomere tracking in genetically engineered Z-line labeled iPSC-derived cardiomyocytes on micropatterned soft gels.

(A) Sketch of a human cardiomyocyte (CM) on soft gel (top) and sarcomere structure in relaxed and contracted state (bottom). (B) ACTN2-Citrine cardiomyocytes (culture day 20) on a polyacrylamide gel substrate (Young’s modulus: 15 kPa), patterned with rectangular stripes of Synthemax (70 × 10 µm). More than 50% of the stripes were typically occupied by single cardiomyocytes. Inset: zoomed-in view of CMs on pattern. (C) Workflow for measurement of sarcomere motion: (1) individual CM adherent to a 15 kPa substrate, first frame of a 1,500-frame confocal time-lapse recording. (2) Deep-learning (Siam-U-Net)-based segmentation of sarcomere Z-bands in the CMs depicted above. (3) Kymograph from the region of interest (ROI) labeled with yellow lines in B and C, demonstrating sarcomere z-band motion. Inset shows Z-band intensity profiles of one time frame with and without Siam-U-Net. (D-E) Confocal images (top row) of representative ACTN2-Citrine-labeled CMs with corresponding deep-learning Z-band segmentation (bottom row) and kymographs (dashed lines mark the start of contraction cycles) recorded from the automatically selected ROI (red lines). (Bottom panels) Overlay plot of single sarcomere length change ΔSL(t) for all tracked sarcomeres in the marked ROIs (first 4 seconds). Thin colored lines are individual sarcomere length changes, the thick black lines display the average of the sarcomere length changes. Contraction intervals are highlighted in gray. Sarcomere popping events are marked with asterisks. Conditions (substrate stiffness): 5 kPa (D); 15 kPa (physiological); E); 85 kPa (F).

Macroscopic muscle contraction, on the one hand, as well as single-molecule kinetics of acto-myosin, on the other hand, have been well studied25. In connecting the microscopic to the macroscopic level, a critical assumption commonly used is that the force that individual sarcomeres exert are primarily a function of sarcomere length (force-length relation)4, i.e., determined by the overlap between actin and myosin filaments, which are protected from slipping apart by passive mechanical elements, such as titin6. A technical limitation of whole-fiber length measurements or diffraction experiments is that they do not report the length of individual sarcomeres, but instead only an average length7,8. Individual sarcomere dynamics and heterogeneity among sarcomeres is thereby averaged out. Furthermore, it is uncertain if the classical model of sarcomere contraction, governed by a steady-state force-length relationship, applies to the rapid periodic dynamics of cardiac muscle, with cyclic activation of acto-myosin by ∼300 ms long Ca2+ transients9.

Individual cardiac sarcomeres have indeed been shown to exhibit spontaneous oscillatory contractions (SPOCs) when sub-maximally activated, although the physiological relevance of this phenomenon remains unclear10. In simple dynamic models of motor proteins, it has been shown by numerical simulations that coupled motors can exhibit complex non-linear phenomena, such as dynamic instabilities and limit-cycle oscillations11,12, some of which have also been observed in in-vitro experiments13. In general, it is known that complex dynamic phenomena can emerge when non-linear dynamic systems are coupled14,15.

In cardiac muscle, deterministic dynamics are vital for efficient and robust organ function. It is thus important to understand at what level in the hierarchical construction of muscle microscopic stochasticity merges into deterministic dynamics. While some intrinsic stochasticity on the single-sarcomere level and heterogeneity among coupled sarcomeres might be unavoidable or even favorable to reduce local stresses, it might under pathological conditions, such as in a fibrotic (stiffened) ventricle, contribute to heart muscle dysfunction. Theoretical studies have explored functional implications of (potentially pathological) static sarcomere non-uniformity, indicating that already small degrees of sarcomere non-uniformity could dramatically affect muscle function 1618. If individual sarcomeres would vary in strength due to damage or disease, weak sarcomeres would continuously be overpowered by their neighbors, which could lead to muscle damage and failure19. Experiments on isolated myofibrils, mostly from skeletal muscle at low temperatures (<20°C) and with constant steady-state calcium activation, found evidence of sarcomere non-uniformity, particularly under isometric and overstretched conditions2023. It remained unclear what these results mean for sarcomere dynamics during in-vivo cardiac contraction at 37°C under transient calcium activation. Recent in-vivo measurements of sarcomere dynamics in mouse hearts found evidence for some asynchronicity24. The physiology of mice hearts with much higher beating frequencies (470 bpm) and expression of distinct sarcomere protein isoforms25 differs substantially from human hearts. It is not known if heterogeneous stochastic sarcomere dynamics can also occur in the human heart at physiological contraction rates (∼60 bpm).

Human pluripotent stem-cell (PSC)-derived CMs are, despite their intrinsic immaturity, well-established as an accessible and scalable platform to investigate human CM structure and function 26. To enhance structural uniformity, to establish a “physiological” aspect ratio with anisotropic alignment of myofibrils, and to support auxotonic contractions, PSC-derived CMs can be cultured on micro-printed cell-adhesive patterns on soft hydrogels27,28. To facilitate structural and functional analyses in living CMs, fluorescent labeling by genome editing using CRISPR/Cas9 in combination with high-resolution microscopy techniques can be employed 2932. We created a human iPSC reporter line with endogenous fluorescent labeling (Citrine) of the myocyte specific z-band protein alpha-actinin 2 (ACTN2) to study sarcomere dynamics in real-time with high spatial and temporal resolution (Härtter et al., submitted33). Culturing of ACTN2-Citrine CMs on micropatterned soft gels with various Young’s moduli (5 - 85 kPa) enabled us to record movies of sarcomere dynamics in single spontaneously beating CMs under different mechanical boundary conditions. For an unbiased analysis of sarcomere dynamics, we employed the deep learning-based Sarcomere Analysis Multi-tool (SarcAsM)33 that can segment and track the motion of individual sarcomere z-bands with 27 nm spatial and 16 ms temporal resolution. Surprisingly, we found that the motility of individual sarcomeres during contraction cycles was almost independent of the external load imposed by the elastic substrates, while heterogeneity between sarcomeres increased strongly with substrate stiffness, leading to a decrease of overall CM contractions. Our findings show that myofibril-level contractility under different load conditions is not regulated within single sarcomeres but emerges in a tug-of-war-like competition among highly dynamic sarcomeres that switch from synchronous contraction to asynchronous dynamics under more rigid constraints. Furthermore, we could show that heterogeneity among sarcomeres is largely stochastic, i.e., not predetermined by static non-uniformity.

Results

Individual ACTN2-Citrine hiPSC-derived CMs on micro-patterned elastic substrates

To resolve sarcomere dynamics with high spatial and temporal resolution in PSC-derived cardiomyocytes, we used a CRISPR-engineered induced pluripotent (iPSC) stem cell line, expressing a yellow fluorescing protein as ACTN2-Citrine N-terminal fusion protein after cardiomyocyte differentiation 33. Seeding of differentiated cardiomyocytes on soft elastic polyacrylamide gels functionalized with a micro-printed pattern promoted anisotropic myofibril assembly in uniformly elongated cardiomyocytes (70 × 10 µm; Fig. 1B). Culturing cells on gels with defined elasticities (Young’s moduli: 5 kPa, 9 kPa, 15 kPa, 29 kPa, 49 kPa, 85 kPa; Table S1) allowed us to impose auxotonic loads on the cells on a scale considered to be relevant in vivo under physiological (10-20 kPa) and pathological (e.g., fibrosis; ≥30 kPa) conditions34. After a maturation period of 20-35 days, we recorded 20-30 s long movies of, in total, 1,362 spontaneously beating CMs at a frame rate of 66 Hz (Movie S1, Fig. S1).

Automated tracking of sarcomere motion at high spatial and temporal resolution

We used SarcAsM, a machine-learning (ML) software tool for automated segmentation of z-bands and analysis of sarcomere structure and dynamics33 for data evaluation. SarcAsM automatically identified 3,985 regions of interest (ROIs), up to 4 per cell, with well-organized registered sarcomeres (≥10 sarcomeres; Fig. 1C, S1). Along each ROI line (∼800 nm in width), an intensity kymograph of Z-band motion was extracted from deep-learning processed (SarcAsM) movies (Fig. 1C). This approach provided more accurate and robust localization and tracking of Z-band trajectories Zi(t) of individual sarcomeres than using raw microscopy data (∼28 nm and 15 ms resolution). From Zi(t), we obtained sarcomere length SLi(t), sarcomere length change ΔSLi(t) and sarcomere velocity Vi(t) of each sarcomere i as well as multi-sarcomere averages and ,(t) for each ROI (details see Ref 33) (Fig. 1D-F).

Cardiomyocyte and individual sarcomere contractility as a function of substrate stiffness

We found that hiPSC-CMs beat spontaneously at a frequency of 0.91 ± 0.38 Hz, which was largely independent of the substrate stiffness (Fig. S2A). Beat-to-beat intervals, though, showed increasing irregularities with increasing substrate stiffness (Fig. S2B). Contraction durations Tcshortened with increasing substrate stiffness (Fig. S2C).

In line with previous reports27,35,36, total cell contraction amplitudes, quantified by the inward motion of the outermost z-bands in myofibrils, substantially decreased with increasing substrate stiffness (Fig. 1D-F). A surprising finding was that the length changes of individual sarcomeres ΔSLifell out of synchronization and became distinctly heterogenous, deviating strongly from the average length change with increasing substrate stiffness. Sarcomeres also showed more and more singular large-amplitude extensions (“sarcomere popping”) far beyond their resting lengths, often very pronounced at the end of contractions (marked with asterisks in bottom row of Fig. 1D-F and large excursions in phase-space plots in Fig. 2A-C). Despite the rapid and heterogenous motion of individual sarcomeres, the emergent contraction at the myofibril scale remained regular and smooth on all substrates.

Analysis of sarcomere dynamics in length change vs. velocity phase-space.

(A-C) Phase-space plots of sarcomere length change ΔSL vs. velocity V for three representative ROIs from CMs on substrates of increasing stiffness (same ROIs as in Fig. 1). Thin black lines show individual sarcomere dynamics, red lines average dynamics. Annotations highlight selected maximal and minimal values of single and average ΔSL and V. (D) Box plot of maximal average contractions as function of substrate stiffness. Maximal average extensions are always close to 0, and not shown. (E) Maximal shortening and lengthening amplitudes ΔSL+/−of individual sarcomeres quantified in each contraction cycle. (F) Maximal average sarcomere lengthening and shortening velocities V+/−. (G) Maximal individual sarcomere lengthening and shortening velocities V+/−. Boxes show quartiles, red lines the median, green triangles the mean and whiskers the 5th and 95th percentile of the distribution per condition. Each data point corresponds to the extremal value within one contraction cycle. To weigh each ROI equally, only the first 10 contraction cycles in each recording were considered. D-G show data of 1,652 ROIs (5 kPa: 122, 9 kPa: 361, 15 kPa: 306, 29 kPa: 306, 49 kPa: 343, 85 kPa: 214). Statistical analysis was performed using the Kruskal-Wallis and Dunn’s posthoc tests, with significance set at p < 0.01. All differences were significant, unless marked (n.s.). (H,I) Time-series and Morlet Wavelet Scalogram of average and single sarcomere length changes ΔSL(t). The top plot displays average (H) and representative single (I) sarcomere length changes over time, with purple areas indicating contraction periods. The bottom plot presents the wavelet scalograms, depicting the evolution of frequency content in the signal over time, with the blue dashed line signifying the cell’s beating rate. (J) Comparison of time-averaged oscillation frequencies of average (black) and single (red) sarcomere length changes of one representative ROI, showing high-frequency intrinsic oscillatory motion of individual sarcomeres with frequencies of 3-4 Hz, which cancel out on the myofibril scale. For time-averaging, only the contraction periods were included. The black curve shows mean ± S.D. of frequencies of 16 sarcomeres in one representative ROI.

To analyze the deviations of individual sarcomere dynamics from the myofibril average during contraction cycles, we collected minima and maxima of ΔSL, i.e., maximal contraction and elongation amplitudes, and of V for each contraction cycle for both individual sarcomeres and averages over all sarcomeres in a given ROI (Fig. 2A-G). In total, we analyzed a data set of 3,985 ROIs recorded from 1,362 cells, each over at least 10 contraction cycles. We excluded irregularly beating (beat-to-beat variability >0.1 s) and very rapidly (>2 Hz) or slowly (<0.5 Hz) beating cells. We further limited the analysis to ROIs spanning whole myofibrils, where the average length change was strictly negative, i.e., contracting. Out of the full data set of 3,985 ROIs, 2,062 ROIs met these criteria. This selection ensured that we focused on rhythmically and uniformly contracting cells. As expected, maximal average contraction amplitudes , measured at the peaks of shortening in each contraction cycle, were the largest (0.18 ± 0.05 µm) on the softest (5 kPa) substrates and were close to zero (0.04 ± 0.03 µm) on the hardest substrates (85 kPa; Fig. 2D). Maximal contraction amplitudes ΔSL of individual sarcomeres in each contraction cycle were also largest on 5 kPa substrates at 0.23 ± 0.09 µm and declined to 0.14 ± 0.07 µm on 85 kPa (Fig. 2E). Note that the decline of maximal individual sarcomere shortening (ΔSL) by ∼30% from 5 to 85 kPa was far less than the ∼75% decline of the average sarcomere contraction (). In addition, while myofibrils never elongated beyond resting length, as expected for auxotonic contractions, individual sarcomeres frequently elongated well beyond their resting length during contraction cycles (median ΔSL+ = 0.05 µm). The distributions of maximal sarcomere extension amplitudes (ΔSL+) were remarkably similar between substrate conditions (Fig. 2E).

Next, we evaluated the effect of substrate elasticity on average and individual sarcomere contraction and extension velocities. The average extremal velocities +and followed the same trend as the extremal average length changes and and were the largest on 5 kPa and declined steadily and strongly with increasing substrate stiffness (from 5 to 85 kPa by 58% for contraction and by 54% for elongation; Fig. 2F). The single extremal contraction and extension velocities V and V+, in contrast, were far less affected by substrate elasticity and differed only by 3% - 6% for the various substrate conditions (Fig. 2G). Note that the maximal extension velocities were 29 - 38% larger than the maximal absolute contraction velocities, showing a strong asymmetry between single sarcomere contraction and extension dynamics. Unlike the average sarcomere motion with clearly distinguishable contraction and relaxation phases, single sarcomeres exhibited rapid switching between slow contractile and fast extensile motion patterns with sometimes multiple shortening and lengthening phases within one contraction cycle (Fig, 2H,I, Fig. 1D-F). To better quantify this phenomenon, we analyzed the oscillation frequencies of the average and individual length changes ΔSL(t) over time, utilizing Morlet wavelet analysis, a method that makes it possible to extract the instantaneous oscillation frequencies of a signal. The analysis revealed a narrow range of frequencies surrounding the beating rate in the scalogram of average sarcomere length changes (Fig. 2H,I). Conversely, for individual sarcomeres, we found a broader range of frequencies, with notable peaks occurring at the dominant and macroscopically visible whole-cell beating rate as well as a surprising, large secondary peak at 3 ± 0.5 Hz (2,062 ROIs, Fig. 2J). This distinctive high-frequency oscillation was consistently observed across numerous ROIs and appeared to be independent of the beating rate. The intrinsic oscillatory motion of individual sarcomeres suggests a mechanism of active acto-myosin powered contraction coupled with rapid lengthening in a relaxation oscillator-like behavior 15.

Correlation analysis of sarcomere motion identifies stochastic and static heterogeneity

The observed heterogeneity of sarcomere dynamics, in particular on rigid substrates, might have been predetermined by static non-uniformities among sarcomeres (e.g., by structural differences in myosin numbers). To distinguish this possibility from stochastically occurring heterogeneity, we extracted the motion ΔSLi,k of each sarcomere i during each contraction cycle k (Fig. 3A). We then calculated the Pearson correlation coefficients r(i,j),(k,l) of each pair of motion patters ΔSLi,k and ΔSLj,l:

and examined the 4-dimensional array of correlation coefficients r(i,j),(k,l). The entries (i = j, kl) are the serial correlations of each sarcomere across different contraction cycles, while the entries (ij, k = l) are the mutual correlations of all sarcomere pairs (i, j) in a given contraction cycle. For each ROI, we calculated the average serial and mutual correlation coefficients rm = ⟨ri≠j,k=l and rs = ⟨ri=j,k≠l for ΔSL and V and examined the distributions of rmand rsfor different substrate conditions (Fig. 3B). The serial and mutual correlations for ΔSL and V were maximal on 5 kPa substrates and declined steadily with increasing substrate stiffness. The decrease of mutual correlation, on the one hand, reflects the decrease of synchrony between sarcomeres due to the tug-of-war competition between sarcomeres imposed by the rigid constraints on the myofibril level (Fig. 3B). The decrease of serial correlation, on the other hand, indicates an increased beat-to-beat variability of the motions of individual sarcomeres (Fig. 3B). Interestingly, the mutual correlation declined more strongly than the serial correlation, by 50% to more than 90% from 5 - 85 kPa. For sarcomere velocity V, we observed a similar trend, however, with a lower base level correlation than ΔSL.

Static versus stochastic heterogeneity of sarcomere motion.

(A) Concatenated contraction periods from a representative ROI in one beating cell on a 15 kPa substrate showing sarcomere length changes during activation, labeled as ΔSLi,k with sarcomere number i and contraction cycle number k. Each box is 0.6 s in width and −0.5 to 0.5 µm in height. Insets show correlation scatter plots from the ΔSLi,k-pairs marked in red (serial correlation) and blue (mutual correlation) respectively in A. The respective Pearson correlation coefficient r(i,j),(k,l) (see Eq. 1) is noted in the graphs. (B) (left column) Average serial correlation rsbetween ΔSL and V (i = j) of individual sarcomere from different contraction cycles (kl), quantifying the variability of motion patterns of individual sarcomere across contraction cycles. (Middle column) Average mutual correlation rmof ΔSL and V between all sarcomeres in a myofibril (ij) in each contraction cycle (k = l), quantifying the variability between motions of neighboring sarcomeres in each contraction cycle. (Right column) Ratio R between average mutual and serial correlation of sarcomere ΔSL and V respectively. R is a measure for the degree of stochasticity in the motions, purely stochastic when R = 1, or largely static when R ≪ 1. (C) Average Pearson correlation coefficients of ΔSL/Vi,k and ΔSL/Vj,l. The x-axis shows the mutual correlation of the motion of different sarcomeres (ij), the y-axis the serial correlation of the different cycles of one sarcomere (i = j), while (kl). The set of dashed lines mark regions of stochastic heterogeneity (right) and static heterogeneity (left). The yellow star marks data of the ROI in A. B,C show data of 2,062 ROIs (5 kPa: N = 134, 9 kPa: 442, 15 kPa: 339, 29 kPa: 393, 49 kPa: 449, 85 kPa: 305). Boxes show quartiles, red lines the median, green triangles the mean and whiskers the 5th and 95th percentile of the distribution per condition. Statistical analysis was performed using the Kruskal-Wallis and Dunn’s posthoc tests, with significance set at p < 0.01. All differences were significant. (D,E) Illustrative sketches of sarcomere length changes ΔSL in different contraction cycles, assuming purely static heterogeneity (D) and purely stochastic heterogeneity (E). Color maps denotes sarcomere contractile strength.

We introduce the ratio R of serial to mutual correlation coefficients to distinguish static and stochastic heterogeneity (Fig. 3B,C) as:

If R = 1, the heterogeneity is perfectly stochastic. A random shuffling of sarcomeres (i, j) and contraction cycles (k, l) would not affect R. If R ≪ 1, the heterogeneity is largely static. In this case, mutual correlations are much smaller than serial correlations, indicating that the motion varies much more between sarcomeres than for a given sarcomere from one contraction cycle to another (Fig. 3D,E). The broad distribution of R values for all examined ROIs shows that heterogeneity was neither fully stochastic nor fully static, but rather distributed across the full spectrum from stochastic to static (Fig. 3C). We examined R for different substrate conditions and found a strong dependence on substrate elasticity for both length change and velocity correlations (Fig. 3B). Heterogeneity is mostly stochastic on 5 kPa and 9 kPa substrates and increasingly static at stiffer, non-physiological elasticities.

Sarcomere popping events are stochastically independent and not only occur in structurally “weak” sarcomeres

The large spread in lengthening amplitudes and velocities (Fig. 2E,G) we observed on all substrates was caused by isolated large elongations of some sarcomeres in a given myofibril beyond their resting lengths during phases of contraction, mostly at the end of systoles (Figs. 1D-F). This all-or-nothing phenomenon of strong elongation, to up to 0.5 µm beyond resting length, is hereafter denoted as sarcomere “popping”. A similar phenomenon has been described in skeletal muscle (38), but has, to our knowledge, not been shown experimentally in CMs. Since sarcomeres can only actively contract, such elongation is likely a passive relaxation after a threshold tension causes an avalanche-like release of all myosin heads from the actin filaments in an individual sarcomere while connected sarcomeres keep contracting. To analyze when and where sarcomeres pop, we identified popping events in all sarcomere ROIs by detecting extensions beyond 0.25 µm (Fig. 4A,B). For each ROI we then created a binary 2D array (popping = 1, no popping = 0; Fig. 4C). Based on that data, we analyzed the temporal and spatial distribution of popping events. We calculated the overall popping frequency n(P) = #total events/(#cycles x #sarcomere), as well as the frequency for each contraction cycle nc(P) = #events in a given cycle/#sarcomeres and for each individual sarcomere ns(P) = #events in a given sarcomere/#cycles. Substrate elasticity did affect the popping frequency (Fig. 4D), with popping getting the more frequent the stiffer the substrate was. On 9 kPa and 49 kPa substrates, popping events were 40% and 95% more frequent than on a 5 kPa substrate, respectively. Because popping frequency of single sarcomeres might be affected by sarcomere-specific characteristics, such as the equilibrium length SL0, we assessed the effect of SL0 on popping frequency ns(P) of single sarcomeres (Fig. 4E). “Longer” sarcomeres with SL0 ≈ 2 µm popped significantly more often than “shorter” sarcomeres with SL0 = 1.5 µm, on 5 and 9 kPa by a factor of ∼5, and on 15 - 85 kPa by a factor of ∼2. Popping frequencies were never smaller than 5% under any condition (Fig. 4D,E).

Analysis of sarcomere popping.

(A) Representative ΔSL time-series of one section of a myofibril containing 16 sarcomeres in one representative CM on a 29 kPa substrate. Popping events, the elongation of a sarcomere within one contraction cycle beyond 0.25 µm, are marked in red. (B) Zoomed plot of 4 consecutive contraction cycles showing two popping events (shaded red). (C) Spatio-temporal map of popping events (black), from data shown in A. Upper and right bar graphs show marginal popping frequencies per contraction cycle nc(P) or per sarcomere ns(P). (D) Overall popping frequencies for different substrate conditions. Boxes show quartiles, red lines the median, green triangles the mean and whiskers the 5th and 95th percentile of the distribution per condition. Statistical analysis was performed using the Kruskal-Wallis and Dunn’s posthoc tests, with significance set at p < 0.01. All differences were significant. (E) Popping frequency as a function of sarcomere equilibrium length SL0. Lines show average of all ROIs for the respective substrate stiffness. (F) Randomly generated events with probability p = 0.3, demonstrating that apparent clustering of popping events shows up even in purely random event sequences. (G,H) Distributions of distance d and time τ between popping events for single ROI in comparison with corresponding geometric distribution (red line). D and E show data of 2,062 ROIs.

Popping events appeared to occur randomly in time and in uncorrelated locations along the myofibrils. Note that apparent accumulations of popping events also appear with simulated Bernoulli-distributed random events (Fig. 4F). If popping events were indeed random in time and location (i.e., be stochastically independent), both the distances between popping events within a myofibril d and the time gaps between two popping events in a given sarcomere τ, would follow a geometric distribution, equivalent to the time interval between successes in a number of Bernoulli trials (1 = success, 0 = failure) with probability p 37. The geometric distribution G(k) for k ∈ {d,τ} reads:

We assessed the distributions of d and τ for all ROIs, and compared them with the geometric distributions G(k) with respective event probability p = n(P) (Fig. 4G,H). Using the Kolmogorov-Smirnov test, we could not reject the hypothesis, that popping events are stochastically independent, for 49% of ROIs for distances d, and 47% of ROIs for time gaps τ (p-value = 0.01). Given the considerable number of regions of interest (ROIs) exhibiting stochastically independent popping events, we hypothesize that while the prior probabilities for popping can indeed be influenced by individual sarcomere characteristics, such as length, strengths, or damaged morphology, popping also happens in myofibrils with uniformly structured sarcomeres across all substrate stiffness conditions.

Discussion

Culturing human stem-cell derived CMs, genome edited to produce fluorescent Z-band markers, made it possible to study sarcomere dynamics in regularly beating CMs under close-to-physiological conditions. Hidden behind the regular contractions of the whole cells, we discovered rich dynamic behavior of individual sarcomeres, including rapid oscillatory motions, stochastic heterogeneity and large-amplitude sarcomere popping.

While whole cells contracted less on stiffer substrates, substrate stiffness had only a weak inhibitory effect on single sarcomere motion amplitudes. Instead, increasing external constraints by more rigid substrates forced sarcomeres into a tug-of-war like competition leading to loss of synchrony. A possible cause of this behavior is an intrinsic dynamical instability of sarcomeres where, at a threshold force, myosin attachment to actin is lost and sarcomeres extend rapidly. This can be formally described by a non-monotonicity of the force-velocity relation, as suggested by biophysical theories for coupled molecular motors1113. Each CM sarcomere contains ∼105-106 myosin motors and actin filaments. Therefore, extension events are a collective, avalanche-like phenomenon in a given sarcomere, i.e. a dynamic phase transition between active shortening and fast passive lengthening/slipping. Once the tension on the acto-myosin system is partially released, myosin heads gain traction again and the sarcomere transitions again to the shortening mode, resulting in high-frequency relaxation oscillations. The stiffer the underlying substrate, the more and earlier the individual sarcomeres pass the force threshold and switch to lengthening mode.

Through detailed analysis of sarcomere trajectories and the occurrence of popping events, our data suggest that both static and stochastic heterogeneity among sarcomeres contribute to the observed dynamics. Static non-uniformities among sarcomeres, perhaps due to differences in the number of molecular components or structural defects, might predispose the motion patterns of individual sarcomeres. We did, however, under all conditions, observe a substantial number of ROIs that showed a high degree of stochastic heterogeneity (R ≈ 1) and sarcomere popping events occurring in a stochastically independent manner. Randomly occurring extensions would in that case ultimately be caused by stochastic acto-myosin kinetics, generating fluctuating forces. It is characteristic for dynamic instabilities in complex systems that small stochastic fluctuations get amplified and are thus fully sufficient to drive highly divergent dynamics even of structurally uniform sarcomeres.

The popping of individual sarcomeres towards the end of contractions indicates that sarcomere dynamic instabilities are especially pronounced during the contraction cycle phase in which acto-myosin coupling is reduced and thus the sarcomere is less stable. Interestingly, popping appeared to be fully reversible and was largely independent of beating frequencies and occurred on all substrates. We believe, that popping events are not aberrations but intrinsic aspects of sarcomere dynamics, manifesting stochastically in structurally uniform sarcomeres in healthy and well-ordered cardiomyocytes.

In summary, these findings suggest that stochastic heterogeneity and popping are intrinsic non-pathological phenomena, while static non-uniformities partially destroy this “stochastic balance” among sarcomeres. In healthy adult cardiomyocytes with more ordered sarcomere morphologies than in stem-cell derived cardiomyocytes38, static non-uniformities are expected to be much less pronounced, but stochastic heterogeneity could persist. This is indeed suggested by recent in-situ measurements of single sarcomere dynamics in mice24.

The (patho-)physiological roles and potential benefits of sarcomere stochastic heterogeneity and popping for cardiac muscle function are unexplored. Sarcomere popping has been previously reported in skeletal muscle in the context of residual force enhancement during active overstretch39. We hypothesize that popping of some sarcomeres in cardiac muscle towards the end of contractions might speed up the return to equilibrium length on myofibril level at the onset of quiescent period, while the individual sarcomeres are still non-uniform in length and relax more slowly. One role of the stochasticity might be that possibly unavoidable extreme stresses would be uniformly distributed among sarcomeres, which would enhance system robustness and reduce the likelihood of mechanical damage. Non-uniformities at a low intrinsic level would be compensated by the inherent stochasticity and not translate to static heterogeneity, e.g., making one marginally weaker sarcomere always pop. Diseases affecting sarcomeric structure and function could perturb this balance by introducing additional static non-uniformity and thus lead to damage19,38,40. In future research, an in-depth analysis of the interplay between local sarcomere-level and cell-level morphology and function might create further insight into the relevance of heterogeneity and stochasticity in healthy and diseased cardiac muscle.

Materials and Methods

Generation and culturing of hiPSC ACTN2-Citrine-derived cardiomyocytes

Cardiomyocyte differentiation of a ACTN2-Citrine reporter hiPSC line33 was performed according to Tiburcy et al.41. hiPSC-ACTN2-Citrine-derived cardiomyocytes were cultured in 6 well plates (Cat 3516, Corning) in serum-free “cardio” medium (0.4 mM Ca2+, RPMI 1640 with GlutaMAX (Cat 61870, Invitrogen), 1% Penicillin/Streptomycin (Cat 15140, Invitrogen), 2% B27 supplement (Cat 17504-044, Invitrogen) at 37 °C in a 5% CO2 incubator with culture medium changes at every other day. For re-seeding, cells were detached using Accutase® digestion medium (StemPro® Accutase® cell dissociation reagent (Cat A11105-01, Gibco), 0.025% Trypsin (Cat 15090-046, Gibco), 20 µg/mL DNaseI (Cat 260913, Calbiochem) for 15-20 min at 37°C. Digestion was stopped using “cardio” medium supplemented with 5 µmol/L Rock Inhibitor (Stemolecule Y27632, Cat 04-0012-10, Reprocell) three times the volume of the Accutase® mix. Cell clumps were separated using a 100 µm cell strainer. Cells were seeded using ∼150,000 cells per micropatterned patterned substrate of ∼1 cm2 size and cultured 24 hours in “cardio” medium with 5 µmol/L Rock Inhibitor. Cardiomyocytes on soft gels were maintained in serum-free “cardio” medium with daily medium changes at 37 °C in a 5% CO2 incubator for up to 30 days. Cardiomyocytes were imaged after a maturation period of 20-30 days post seeding on the soft gels.

Sample preparation and measurement

For live-cell video imaging, the soft substrates were mounted in a custom-built holder for round Ø25 mm glass cover slides with serum-free medium. Movies of beating cardiomyocytes were obtained with a confocal microscope (TCS SP5 II, Leica, Germany) at 37°C and 5% CO2. We used a 8,000 Hz resonant scanner with bidirectional scanning mode and recorded up to 20-30 s long movies of 1024 × 200 pixels and a temporal resolution of 15 ms / 67 frames per second.

Fabrication of cell-adhesive micropatterns on polyacrylamide soft gels

Silicon wafers (Microchemicals GmbH, Ulm, Germany) were coated with SU-8 photoresist (Series 3005, MicroChem, Newton, USA) using a two-step spin-coating process. The wafers were then exposed to UV light through custom photomasks (Compugraphics, Jena, Germany) and developed to create photoresist masters (details in Ref33). Polydimethylsiloxane (PDMS) stamps were produced by mixing PDMS and curing agent 10:1 (Sylgard 184 kit, Dow Corning) and pouring it onto the photoresist masters. After degassing and curing, the PDMS was cut and peeled off to create the stamps. Micropatterned polyacrylamide gels were prepared by treating PDMS stamps with plasma to make them hydrophilic, then incubating them with 0.1 mg/ml Synthemax™ (Cat 3535, Corning). The protein-coated stamps were placed on plasma-cleaned glass coverslips and weighted to transfer the protein. Gel solutions of acrylamide and bis-acrylamide were prepared to achieve different elastic moduli (Table S1), measured with a rheometer (Physica MCR 501 rheometer, Anton Paar, Austria). The gel solution was polymerized with the Synthemax-patterned glass on top. The gels were stored in PBS, then washed before use.

Microscopy data analysis and statistical analysis

For processing of the live-cell confocal movies, and tracking and analysis of sarcomere trajectories, our Python package SarcAsM (Sarcomere Analysis Multitool) was used (https://github.com/danihae/SarcAsM, repository private until publication)33. All additional analyses and illustrations were created using custom scripts written in Python. All data is displayed as mean ± standard deviation unless indicated otherwise. Whenever applicable, we employed the Kruskal-Wallis test, a non-parametric method, for hypothesis testing, and proceeded to a post-hoc multi-comparison using Dunn’s test, setting the threshold for statistical significance at a p-value of less than 0.01. All condition differences were deemed significant except where explicitly marked as not significant (n.s.).

Data and materials availability

All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Raw microscopy data and code for data analysis are available upon request.

Acknowledgements

We thank Florian Rehfeldt, Andrei Vilfan, Lev Truskinowsky, David Brückner, Chase Broedersz and Pierre Ronceray for helpful discussions. We gratefully acknowledge the use of the microscopy facility of the Max Planck Institute for Multidisciplinary Sciences for access to cell culture and high-speed confocal microscopy. CFS and DH would like to thank the Isaac Newton Institute for Mathematical Sciences for support and hospitality during the programme “New statistical physics in living matter: non equilibrium states under adaptive control”. DH acknowledges the support from the German Academic Foundation (Studienstiftung des Deutschen Volkes) for providing a doctoral fellowship and the Campus Institute for Data Science (CIDAS) at the University of Göttingen for awarding a postdoctoral fellowship. The research of CFS was supported by the European Research Council under the European Union’s Seventh Framework Programme (FP7/2007-2013) with the ERC grant agreement n°340528. WHZ acknowledges the support from the DZHK (German Center for Cardiovascular Research), the German Federal Ministry of Education and Research (IndiHEART; 161L0250A), the German Research Foundation (DFG SFB 1002 C04/S01, IRTG 1816, RTG 2824, EXC 2067-1), and the Fondation Leducq (20CVD04).

Author contributions

D.H., T.D., W.H.Z. and C.F.S. conceived and designed the study. D.H., L.H., and T.D. developed the methodology. WHZ and CFS supplied the essential resources. DH, LH, and TD performed micropatterning and cardiomyocyte experiments. DH and KN analyzed data. DH created the visualizations. The original draft was written by DH, and all authors contributed to the review and editing of the manuscript.

Competing interests

WHZ is founder and holds equity of myriamed GmbH and Repairon GmbH. All other authors declare that they have no competing interests.

Supplementary Information

Representative CM selection with automatically and individually detected ROIs.

CMs are randomly selected from the data sets for three substrate stiffnesses (9, 15 and 85 kPa); automatically determined ROIs depicted as red lines.

Effects of substrate elasticity on spontaneous beating of cardiomyocytes.

(A) Spontaneous beating frequencies. (B) Beating irregularity: relative standard deviation of beating periods. (C) Average length of contraction cycles Tc (identified as systole). A-C include data from 3,985 ROIs. Boxes show quartiles, red lines the median, green triangles the mean and whiskers the 5th and 95th percentile of the distributions.

Composition of polyacrylamide soft gels used in the study.

All gels used here were made from the same respective stock solution (10 ml). The Young’s modulus was measured for a gel sample from each prepared stock solution using a rheometer (Physica MCR 501 Rheometer, Anton Paar, Germany) using a 25 mm, 2° cone plate with a sample volume of 140 μl. A time sweep (1 h, spacing 30 s, 1% strain, 1 Hz), a frequency sweep (3 measurements per decade, 1% strain, 0.01 – 100 Hz), and an amplitude sweep (3 measurements per decade, 0.01 – 100% strain, 1 Hz) were performed consecutively.