A three-dimensional simulation model of cardiomyocyte integrating excitation-contraction coupling and metabolism

Multi-Scale Computational Modeling of Spatial Calcium Handling From Nanodomain to Whole-Heart: Overview and Perspectives

Regulation of intracellular calcium is a critical component of cardiac electrophysiology and excitation-contraction coupling. The calcium spark, the fundamental element of the intracellular calcium transient, is initiated in specialized nanodomains which co-locate the ryanodine receptors and L-type calcium channels. However, calcium homeostasis is ultimately regulated at the cellular scale, by the interaction of spatially separated but diffusively coupled nanodomains with other sub-cellular and surface-membrane calcium transport channels with strong non-linear interactions; and cardiac electrophysiology and arrhythmia mechanisms are ultimately tissue-scale phenomena, regulated by the interaction of a heterogeneous population of coupled myocytes. Recent advances in imaging modalities and image-analysis are enabling the super-resolution reconstruction of the structures responsible for regulating calcium homeostasis, including the internal structure of nanodomains themselves. Extrapolating functional and imaging data from the nanodomain to the whole-heart is non-trivial, yet essential for translational insight into disease mechanisms. Computational modeling has important roles to play in relating structural and functional data at the sub-cellular scale and translating data across the scales. This review covers recent methodological advances that enable image-based modeling of the single nanodomain and whole cardiomyocyte, as well as the development of multi-scale simulation approaches to integrate data from nanometer to whole-heart. Firstly, methods to overcome the computational challenges of simulating spatial calcium dynamics in the nanodomain are discussed, including image-based modeling at this scale. Then, recent whole-cell models, capable of capturing a range of different structures (such as the T-system and mitochondria) and cellular heterogeneity/variability are discussed at two different levels of discretization. Novel methods to integrate the models and data across the scales and simulate stochastic dynamics in tissue-scale models are then discussed, enabling elucidation of the mechanisms by which nanodomain remodeling underlies arrhythmia and contractile dysfunction. Perspectives on model differences and future directions are provided throughout.

Isolation and reconstruction of cardiac mitochondria from SBEM images using a deep learning-based method

Mitochondrial morphological defects are a common feature of diseased cardiac myocytes. However, quantitative assessment of mitochondrial morphology is limited by the time-consuming manual segmentation of electron micrograph (EM) images. To advance understanding of the relation between morphological defects and dysfunction, an efficient morphological reconstruction method is desired to enable isolation and reconstruction of mitochondria from EM images. We propose a new method for isolating and reconstructing single mitochondria from serial block-face scanning EM (SBEM) images. CDeep3M, a cloud-based deep learning network for EM images, was used to segment mitochondrial interior volumes and boundaries. Post-processing was performed using both the predicted interior volume and exterior boundary to isolate and reconstruct individual mitochondria. Series of SBEM images from two separate cardiac myocytes were processed. The highest F1-score was 95% using 50 training datasets, greater than that for previously reported automated methods and comparable to manual segmentations. Accuracy of separation of individual mitochondria was 80% on a pixel basis. A total of 2315 mitochondria in the two series of SBEM images were evaluated with a mean volume of 0.78 µm³. The volume distribution was very broad and skewed; the most frequent mitochondria were 0.04-0.06 µm³, but mitochondria larger than 2.0 µm³ accounted for more than 10% of the total number. The average short-axis length was 0.47 µm. Primarily longitudinal mitochondria (0-30 degrees) were dominant (54%). This new automated segmentation and separation method can help quantitate mitochondrial morphology and improve understanding of myocyte structure-function relationships.

Geometric Understanding of Local Fluctuation Distribution of Conduction Time in Lined-Up Cardiomyocyte Network in Agarose-Microfabrication Multi-Electrode Measurement Assay

We examined characteristics of the propagation of conduction in width-controlled cardiomyocyte cell networks for understanding the contribution of the geometrical arrangement of cardiomyocytes for their local fluctuation distribution. We tracked a series of extracellular field potentials of linearly lined-up human embryonic stem (ES) cell-derived cardiomyocytes and mouse primary cardiomyocytes with 100 kHz sampling intervals of multi-electrodes signal acquisitions and an agarose microfabrication technology to localize the cardiomyocyte geometries in the lined-up cell networks with 100–300 μm wide agarose microstructures. Conduction time between two neighbor microelectrodes (300 μm) showed Gaussian distribution. However, the distributions maintained their form regardless of its propagation distances up to 1.5 mm, meaning propagation diffusion did not occur. In contrast, when Quinidine was applied, the propagation time distributions were increased as the faster firing regulation simulation predicted. The results indicate the “faster firing regulation” is not sufficient to explain the conservation of the propagation time distribution in cardiomyocyte networks but should be expanded with a kind of community effect of cell networks, such as the lower fluctuation regulation.

Dominant rule of community effect in synchronized beating behavior of cardiomyocyte networks

Exploiting the combination of latest microfabrication technologies and single cell measurement technologies, we can measure the interactions of single cells, and cell networks from "algebraic" and "geometric" perspectives under the full control of their environments and interactions. However, the experimental constructive single cell-based approach still remains the limitations regarding the quality and condition control of those cells. To overcome these limitations, mathematical modeling is one of the most powerful complementary approaches. In this review, we first explain our on-chip experimental methods for constructive approach, and we introduce the results of the "community effect" of beating cardiomyocyte networks as an example of this approach. On-chip analysis revealed that (1) synchronized interbeat intervals (IBIs) of cell networks were followed to the more stable beating cells even their IBIs were slower than the other cells, which is against the conventional faster firing regulation or "overdrive suppression," and (2) fluctuation of IBIs of cardiomyocyte networks decreased according to the increase of the number of connected cells regardless of their geometry. The mathematical simulation of this synchronous behavior of cardiomyocyte networks also fitted well with the experimental results after incorporating the fluctuation-dissipation theorem into the oscillating stochastic phase model, in which the concept of spatially arranged cardiomyocyte networks was involved. The constructive experiments and mathematical modeling indicated the dominant rule of synchronization behavior of beating cardiomyocyte networks is a kind of stability-oriented synchronization phenomenon as the "community effect" or a fluctuation-dissipation phenomenon. Finally, as a practical application of this approach, the predictive cardiotoxicity is introduced.

Clinical and pharmacological application of multiscale multiphysics heart simulator, UT-Heart

A heart simulator, UT-Heart, is a finite element model of the human heart that can reproduce all the fundamental activities of the working heart, including propagation of excitation, contraction, and relaxation and generation of blood pressure and blood flow, based on the molecular aspects of the cardiac electrophysiology and excitation-contraction coupling. In this paper, we present a brief review of the practical use of UT-Heart. As an example, we focus on its application for predicting the effect of cardiac resynchronization therapy (CRT) and evaluating the proarrhythmic risk of drugs. Patient-specific, multiscale heart simulation successfully predicted the response to CRT by reproducing the complex pathophysiology of the heart. A proarrhythmic risk assessment system combining in vitro channel assays and in silico simulation of cardiac electrophysiology using UT-Heart successfully predicted druginduced arrhythmogenic risk. The assessment system was found to be reliable and efficient. We also developed a comprehensive hazard map on the various combinations of ion channel inhibitors. This in silico electrocardiogram database (now freely available at http://ut-heart.com/) can facilitate proarrhythmic risk assessment without the need to perform computationally expensive heart simulation. Based on these results, we conclude that the heart simulator, UT-Heart, could be a useful tool in clinical medicine and drug discovery.

Three dimensional coupled reaction–diffusion modeling of calcium and inositol 1,4,5-trisphosphate dynamics in cardiomyocytes

Nanoparticles have shown great promise in improving cancer treatment efficacy by changing the intracellular calcium level through activation of intracellular mechanisms. One of the mechanisms of the killing of the cancerous cell by a nanoparticle is through elevation of the intracellular calcium level. Evidence accumulated over the past decade indicates a pivotal role for the IP₃ receptor mediated Ca²⁺ release in the regulation of the cytosolic and the nuclear Ca²⁺ signals. There have been various studies done suggesting the role of IP₃ receptors (IP₃R) and IP₃ production and degradation in cardiomyocytes. In the present work, we have proposed a three-dimensional unsteady-state mathematical model to describe the mechanism of cardiomyocytes which focuses on evaluation of various parameters that affect these coupled dynamics and elevate the cytosolic calcium concentration which can be helpful to search for novel therapies to cure these malignancies by targeting the complex calcium signaling process in cardiomyocytes. Our study suggests that there are other factors involved in this signaling which can increase the calcium level, which can help in finding treatment for cancer. The cytosolic calcium level may be controlled by IP₃ signaling, leak, source influx of calcium (σ) and maximum production of IP₃ (V_P). We believe that the proposed model suggests new insight into finding treatment for cancer in cardiomyocytes through elevation of the cytosolic Ca²⁺ concentration by various parameters like leak, σ, V_P and especially by other complex cell signaling dynamics, namely IP₃ dynamics.

We propose a three-dimensional unsteady-state mathematical model to describe the mechanism of cardiomyocytes.

Computational Model of Calcium Signaling in Cardiac Atrial Cells at the Submicron Scale

In cardiac cells, calcium is the mediator of excitation-contraction coupling. Dysfunctions in calcium handling have been identified as the origin of some cardiac arrhythmias. In the particular case of atrial myocytes, recent available experimental data has found links between these dysfunctions and structural changes in the calcium handling machinery (ryanodine cluster size and distribution, t-tubular network, etc). To address this issue, we have developed a computational model of an atrial myocyte that takes into account the detailed intracellular structure. The homogenized macroscopic behavior is described with a two-concentration field model, using effective diffusion coefficients of calcium in the sarcoplasmic reticulum (SR) and in the cytoplasm. The model reproduces the right calcium transients and dependence with pacing frequency. Under basal conditions, the calcium rise is mostly restricted to the periphery of the cell, with a large concentration ratio between the periphery and the interior. We have then studied the dependence of the speed of the calcium wave on cytosolic and SR diffusion coefficients, finding an almost linear relation with the former, in agreement with a diffusive and fire mechanism of propagation, and little dependence on the latter. Finally, we have studied the effect of a change in RyR cluster microstructure. We find that, under resting conditions, the spark frequency decreases slightly with RyR cluster spatial dispersion, but markedly increases when the RyRs are distributed in clusters of larger size, stressing the importance of RyR cluster organization to understand atrial arrhythmias, as recent experimental results suggest (Macquaide et al., 2015).

Insights on the impact of mitochondrial organisation on bioenergetics in high-resolution computational models of cardiac cell architecture

Recent electron microscopy data have revealed that cardiac mitochondria are not arranged in crystalline columns but are organised with several mitochondria aggregated into columns of varying sizes spanning the cell cross-section. This raises the question—how does the mitochondrial arrangement affect the metabolite distributions within cardiomyocytes and what is its impact on force dynamics? Here, we address this question by employing finite element modeling of cardiac bioenergetics on computational meshes derived from electron microscope images. Our results indicate that heterogeneous mitochondrial distributions can lead to significant spatial variation across the cell in concentrations of inorganic phosphate, creatine (Cr) and creatine phosphate (PCr). However, our model predicts that sufficient activity of the creatine kinase (CK) system, coupled with rapid diffusion of Cr and PCr, maintains near uniform ATP and ADP ratios across the cell cross sections. This homogenous distribution of ATP and ADP should also evenly distribute force production and twitch duration with contraction. These results suggest that the PCr shuttle and associated enzymatic reactions act to maintain uniform force dynamics in the cell despite the heterogeneous mitochondrial organization. However, our model also predicts that under hypoxia activity of mitochondrial CK enzymes and diffusion of high-energy phosphate compounds may be insufficient to sustain uniform ATP/ADP distribution and hence force generation.

Mammalian cardiomyocytes contain a high volume of mitochondria, which maintains the continuous and bulk supply of ATP to sustain normal heart function. Previously, cardiac mitochondria were understood to be distributed in a regular, crystalline pattern, which facilitated a steady supply of ATP at different workloads. Using electron microscopy images of cell cross sections, we recently found that they are not regularly distributed inside cardiomyocytes. We created new spatially accurate computational models of cardiac cell bioenergetics and tested whether this heterogeneous distribution of mitochondria causes non-uniform energy supply and contractile force production in the cardiomyocyte. We found that ATP and ADP concentrations remain uniform throughout the cell because of the activity of creatine kinase (CK) enzymes that convert ATP produced in the mitochondria into creatine phosphate. Creatine phosphate rapidly diffuses to the myofibril region where it can be converted back to ATP for the contraction cycle in a timely manner. This mechanism is called the phosphocreatine shuttle (PCr shuttle). The PCr shuttle ensures that different areas of the cell produce the same amount of force regardless of the mitochondrial distribution. However, our model also shows that when the cellular oxygen supply is limited—as can be the case in conditions such as heart failure—the PCr shuttle cannot maintain uniform ATP and ADP concentrations across the cell. This causes a non-uniform acto-myosin force distribution and non-uniform twitch duration across the cell cross section. Our study suggests that mechanisms other than the PCr shuttle may be necessary to maintain uniform supply of ATP in a hypoxic environment.

A 3D diffusional-compartmental model of the calcium dynamics in cytosol, sarcoplasmic reticulum and mitochondria of murine skeletal muscle fibers

Variations of free calcium concentration ([Ca2+]) are powerful intracellular signals, controlling contraction as well as metabolism in muscle cells. To fully understand the role of calcium redistribution upon excitation and contraction in skeletal muscle cells, the local [Ca2+] in different compartments needs to be taken into consideration. Fluorescent probes allow the determination of [Ca2+] in the cytosol where myofibrils are embedded, the lumen of the sarcoplasmic reticulum (SR) and the mitochondrial matrix. Previously, models have been developed describing intracellular calcium handling in skeletal and cardiac muscle cells. However, a comprehensive model describing the kinetics of the changes in free calcium concentration in these three compartments is lacking. We designed a new 3D compartmental model of the half sarcomere with radial symmetry, which accounts for diffusion of Ca2+ into the three compartments and simulates its dynamics at rest and at various rates of stimulation in mice skeletal muscle fibers. This model satisfactorily reproduces both the amplitude and time course of the variations of [Ca2+] in the three compartments in mouse fast fibers. As an illustration of the applicability of the model, we investigated the effects of Calsequestrin (CSQ) ablation. CSQ is the main Ca2+ buffer in the SR, localized in close proximity of its calcium release sites and near to the mitochondria. CSQ knock-out mice muscles still preserve a near-normal contractile behavior, but it is unclear whether this is caused by additional SR calcium buffering or a significant contribution of calcium entry from extracellular space, via stored-operated calcium entry (SOCE). The model enabled quantitative assessment of these two scenarios by comparison to measurements of local calcium in the cytosol, the SR and the mitochondria. In conclusion, the model represents a useful tool to investigate the impact of protein ablation and of pharmacological interventions on intracellular calcium dynamics in mice skeletal muscle.

Quantitative systems models illuminate arrhythmia mechanisms in heart failure: Role of the Na+ -Ca2+ -Ca2+ /calmodulin-dependent protein kinase II-reactive oxygen species feedback

Quantitative systems modeling aims to integrate knowledge in different research areas with models describing biological mechanisms and dynamics to gain a better understanding of complex clinical syndromes. Heart failure (HF) is a chronic complex cardiac disease that results from structural or functional disorders impairing the ability of the ventricle to fill with or eject blood. Highly interactive and dynamic changes in mechanical, structural, neurohumoral, metabolic, and electrophysiological properties collectively predispose the failing heart to cardiac arrhythmias, which are responsible for about a half of HF deaths. Multiscale cardiac modeling and simulation integrate structural and functional data from HF experimental models and patients to improve our mechanistic understanding of this complex arrhythmia syndrome. In particular, they allow investigating how disease-induced remodeling alters the coupling of electrophysiology, Ca²⁺ and Na⁺ handling, contraction, and energetics that lead to rhythm derangements. The Ca²⁺ /calmodulin-dependent protein kinase II, which expression and activity are enhanced in HF, emerges as a critical hub that modulates the feedbacks between these various subsystems and promotes arrhythmogenesis. This article is categorized under: Physiology > Mammalian Physiology in Health and Disease Models of Systems Properties and Processes > Mechanistic Models Models of Systems Properties and Processes > Cellular Models Models of Systems Properties and Processes > Organ, Tissue, and Physiological Models.

An automated workflow for segmenting single adult cardiac cells from large-volume serial block-face scanning electron microscopy data

This paper presents a new algorithm to automatically segment the myofibrils, mitochondria and nuclei within single adult cardiac cells that are part of a large serial-block-face scanning electron microscopy (SBF-SEM) dataset. The algorithm only requires a set of manually drawn contours that roughly demarcate the cell boundary at routine slice intervals (every 50th, for example). The algorithm correctly classified pixels within the single cell with 97% accuracy when compared to manual segmentations. One entire cell and the partial volumes of two cells were segmented. Analysis of segmentations within these cells showed that myofibrils and mitochondria occupied 47.5% and 51.6% on average respectively, while the nuclei occupy 0.7% of the cell for which the entire volume was captured in the SBF-SEM dataset. Mitochondria clustering increased at the periphery of the nucleus region and branching points of the cardiac cell. The segmentations also showed high area fraction of mitochondria (up to 70% of the 2D image slice) in the sub-sarcolemmal region, whilst it was closer to 50% in the intermyofibrillar space. We finally demonstrate that our segmentations can be turned into 3D finite element meshes for cardiac cell computational physiology studies. We offer our large dataset and MATLAB implementation of the algorithm for research use at www.github.com/CellSMB/sbfsem-cardiac-cell-segmenter/. We anticipate that this timely tool will be of use to cardiac computational and experimental physiologists alike who study cardiac ultrastructure and its role in heart function.

Dynamic Model for Characterizing Contractile Behaviors and Mechanical Properties of a Cardiomyocyte

Studies on the contractile dynamics of heart cells have attracted broad attention for the development of both heart disease therapies and cardiomyocyte-actuated micro-robotics. In this study, a linear dynamic model of a single cardiomyocyte cell was proposed at the subcellular scale to characterize the contractile behaviors of heart cells, with system parameters representing the mechanical properties of the subcellular components of living cardiomyocytes. The system parameters of the dynamic model were identified with the cellular beating pattern measured by a scanning ion conductance microscope. The experiments were implemented with cardiomyocytes in one control group and two experimental groups with the drugs cytochalasin-D or nocodazole, to identify the system parameters of the model based on scanning ion conductance microscope measurements, measurement of the cellular Young's modulus with atomic force microscopy indentation, measurement of cellular contraction forces using the micro-pillar technique, and immunofluorescence staining and imaging of the cytoskeleton. The proposed mathematical model was both indirectly and qualitatively verified by the variation in cytoskeleton, beating amplitude, and contractility of cardiomyocytes among the control and the experimental groups, as well as directly and quantitatively validated by the simulation and the significant consistency of 90.5% in the comparison between the ratios of the Young's modulus and the equivalent comprehensive cellular elasticities of cells in the experimental groups to those in the control group. Apart from mechanical properties (mass, elasticity, and viscosity) of subcellular structures, other properties of cardiomyocytes have also been studied, such as the properties of the relative action potential pattern and cellular beating frequency. This work has potential implications for research on cytobiology, drug screening, mechanisms of the heart, and cardiomyocyte-based bio-syncretic robotics.

Community effect of cardiomyocytes in beating rhythms is determined by stable cells

The community effect of cardiomyocytes was investigated in silico by the change in number and features of cells, as well as configurations of networks. The theoretical model was based on experimental data and accurately reproduced recently published experimental results regarding coupled cultured cardiomyocytes. We showed that the synchronised beating of two coupled cells was tuned not to the cell with a faster beating rate, but to the cell with a more stable rhythm. In a network of cardiomyocytes, a cell with low fluctuation, but not a hight frequency, became a pacemaker and stabilised the beating rhythm. Fluctuation in beating rapidly decreased with an increase in the number of cells (N), almost irrespective of the configuration of the network, and a cell comes to have natural and stable beating rhythms, even for N of approximately 10. The universality of this community effect lies in the fluctuation-dissipation theorem in statistical mechanics.

Integrate and fire model with refractory period for synchronization of two cardiomyocytes

We investigate an integrate and fire model for two cardiomyocytes interacting with each other. A feature of the model is to incorporate the refractory periods of the cardiomyocytes as well as the influence of firing of adjacent cells. The present model predicts that, if refractory periods of the two cells are nearly equal, the beating rhythms of the two cells always synchronize and their beating rate is tuned to the faster rate between the two cells. On the other hand, if their refractory periods significantly differ, they exhibit various kinds of harmonious beating rhythms. These results successfully explain the well known characteristics of synchronized beating of cultured cardiomyocytes. We also discuss effects of a delay time of cell-to-cell interaction, that gives further complicated phase diagrams for the beating rhythms.

A two dimensional electromechanical model of a cardiomyocyte to assess intra-cellular regional mechanical heterogeneities

Experimental studies on isolated cardiomyocytes from different animal species and human hearts have demonstrated that there are regional differences in the Ca2+ release, Ca2+ decay and sarcomere deformation. Local deformation heterogeneities can occur due to a combination of factors: regional/local differences in Ca2+ release and/or re-uptake, intra-cellular material properties, sarcomere proteins and distribution of the intracellular organelles. To investigate the possible causes of these heterogeneities, we developed a two-dimensional finite-element electromechanical model of a cardiomyocyte that takes into account the experimentally measured local deformation and cytosolic [Ca2+] to locally define the different variables of the constitutive equations describing the electro/mechanical behaviour of the cell. Then, the model was individualised to three different rat cardiac cells. The local [Ca2+] transients were used to define the [Ca2+]-dependent activation functions. The cell-specific local Young's moduli were estimated by solving an inverse problem, minimizing the error between the measured and simulated local deformations along the longitudinal axis of the cell. We found that heterogeneities in the deformation during contraction were determined mainly by the local elasticity rather than the local amount of Ca2+, while in the relaxation phase deformation was mainly influenced by Ca2+ re-uptake. Our electromechanical model was able to successfully estimate the local elasticity along the longitudinal direction in three different cells. In conclusion, our proposed model seems to be a good approximation to assess the heterogeneous intracellular mechanical properties to help in the understanding of the underlying mechanisms of cardiomyocyte dysfunction.

A computational study of the role of mitochondrial organization on cardiac bioenergetics

All cells in the body have a specific shape and internal organization which is specific to that cell's function. Heart cells are rod-shaped, and contain arrays of contractile protines (myofibrils) and mitochondria (organelles that produce energy) that are aligned along the length of the rod. This arrangement is presumed to allow the cell to generate maximal contractile force for each heartbeat and for energy metabolites to be readily available to generate this force. Heart disease phenotypes, such as diabetic cardiomyopathy and heart failure, exhibit altered organization of mitochondria. However, physiological and computational studies have predominantly investigated the effect of the biochemical changes that accompany the disease alone, such as reduced rates of ATP production by mitochondria. We present a modeling study that examines the effect of mitochondrial organization on energy metabolite distribution during the heartbeat. A 2D micrograph of the cell cross-section was selected from a 3D image stack of structural data of a cardiac cell. The image was used to generate a 2D finite element model, on which mitochondrial oxidative phosphorylation and energy metabolite diffusion was modelled. Results illustrate that mitochondrial density can induce heterogeneity in the distribution of metabolites across the cell area. We discuss the implications of these findings and avenues for future, more indepth studies.

Stokes flow patterns induced by a single cardiac cell

Stokes flow motions induced by a beating single cardiac cell (cardiomyocyte) are obtained numerically using the method of fundamental solutions (MFS). A two-dimensional meshfree-Stokeslets computational framework is used to solve the Stokes governing equations around an isolated cardiomyocyte. An approximate beating kinematical model is derived and used to approximate the cell-length shortening over a complete cardiac cycle. The induced flow patterns have been found to be characterized by the presence of counter-rotating vortices at both cell's edges. These vortical flow structures are clearly shown by rendering the velocity streamlines. The static pressure contours are also calculated at different time snapshots during both contraction and relaxation phases of the beating motion. The pressure signal is calculated at a point in the neighborhood of cell surface to capture the induced normal stress (traction) by the cell morphological motions to the surrounding fluid medium. The presented results have shown that, cells with a slightly different shortening/beating profile can induce different flow field. This implies that, each cell is characterized by a unique flow pattern "signature", which potentially can be correlated to the sub-cellular excitation-contraction processes of cardiac cells.

Sarcoplasmic reticulum Ca2+, Mg2+, K+, and Cl- concentrations adjust quickly as heart rate changes

During systole, Ca²⁺ is released from the sarcoplasmic reticulum (SR) through ryanodine receptors (RyRs) while, simultaneously, other ions (specifically K⁺, Mg²⁺, and Cl^-) provide counter-ion flux. These ions move back into the SR during diastole through the SERCA pump and SR K⁺ and Cl^- channels. In homeostasis, all ion concentrations in different cellular regions (e.g., junctional and non-junctional SR, dyadic cleft, and cytosol) are the same at the beginning and end of the cardiac cycle. Here, we used an equivalent circuit compartment model of the SR and the surrounding cytoplasm to understand the heart rate dependence of SR ion homeostasis. We found that the Ca²⁺, Mg²⁺, K⁺, and Cl^- concentrations in the SR and the cytoplasm self-adjust within just a few heartbeats with only very small changes in Mg²⁺, K⁺, and Cl^- concentrations and membrane voltages (just a few percent). However, those small changes were enough to compensate for the large heart-rate-dependent changes in SR and cytoplasmic Ca²⁺ concentrations in the new steady state. The modeling suggests that ion adaptation to increases in heart rate is inherent to the system and that physiological changes that increase contractility and cardiac output are accommodated by the same self-adjusting mechanism of producing small changes in ion driving forces. Our findings also support the long-held hypothesis that SR membrane potentials are small (~1-2mV).

Distinct functional roles of cardiac mitochondrial subpopulations revealed by a 3D simulation model

Experimental characterization of two cardiac mitochondrial subpopulations, namely, subsarcolemmal mitochondria (SSM) and interfibrillar mitochondria (IFM), has been hampered by technical difficulties, and an alternative approach is eagerly awaited. We previously developed a three-dimensional computational cardiomyocyte model that integrates electrophysiology, metabolism, and mechanics with subcellular structure. In this study, we further developed our model to include intracellular oxygen diffusion, and determined whether mitochondrial localization or intrinsic properties cause functional variations. For this purpose, we created two models: one with equal SSM and IFM properties and one with IFM having higher activity levels. Using these two models to compare the SSM and IFM responses of [Ca(2+)], tricarboxylic acid cycle activity, [NADH], and mitochondrial inner membrane potential to abrupt changes in pacing frequency (0.25-2 Hz), we found that the reported functional differences between these subpopulations appear to be mostly related to local [Ca(2+)] heterogeneity, and variations in intrinsic properties only serve to augment these differences. We also examined the effect of hypoxia on mitochondrial function. Under normoxic conditions, intracellular oxygen is much higher throughout the cell than the half-saturation concentration for oxidative phosphorylation. However, under limited oxygen supply, oxygen is mostly exhausted in SSM, leaving the core region in an anoxic condition. Reflecting this heterogeneous oxygen environment, the inner membrane potential continues to decrease in IFM, whereas it is maintained to nearly normal levels in SSM, thereby ensuring ATP supply to this region. Our simulation results provide clues to understanding the origin of functional variations in two cardiac mitochondrial subpopulations and their differential roles in maintaining cardiomyocyte function as a whole.

An integrated finite element simulation of cardiomyocyte function based on triphasic theory

In numerical simulations of cardiac excitation-contraction coupling, the intracellular potential distribution and mobility of cytosol and ions have been mostly ignored. Although the intracellular potential gradient is small, during depolarization it can be a significant driving force for ion movement, and is comparable to diffusion in terms of net flux. Furthermore, fluid in the t-tubules is thought to advect ions to facilitate their exchange with the extracellular space. We extend our previous finite element model that was based on triphasic theory to examine the significance of these factors in cardiac physiology. Triphasic theory allows us to study the behavior of solids (proteins), fluids (cytosol) and ions governed by mechanics and electrochemistry in detailed subcellular structures, including myofibrils, mitochondria, the sarcoplasmic reticulum, membranes, and t-tubules. Our simulation results predicted an electrical potential gradient inside the t-tubules at the onset of depolarization, which corresponded to the Na(+) channel distribution therein. Ejection and suction of fluid between the t-tubules and the extracellular compartment during isometric contraction were observed. We also examined the influence of t-tubule morphology and mitochondrial location on the electrophysiology and mechanics of the cardiomyocyte. Our results confirm that the t-tubule structure is important for synchrony of Ca(2+) release, and suggest that mitochondria in the sub-sarcolemmal region might serve to cancel Ca(2+) inflow through surface sarcolemma, thereby maintaining the intracellular Ca(2+) environment in equilibrium.

A multiscale computational model of spatially resolved calcium cycling in cardiac myocytes: from detailed cleft dynamics to the whole cell concentration profiles

Mathematical modeling of excitation-contraction coupling (ECC) in ventricular cardiac myocytes is a multiscale problem, and it is therefore difficult to develop spatially detailed simulation tools. ECC involves gradients on the length scale of 100 nm in dyadic spaces and concentration profiles along the 100 μm of the whole cell, as well as the sub-millisecond time scale of local concentration changes and the change of lumenal Ca²⁺ content within tens of seconds. Our concept for a multiscale mathematical model of Ca²⁺ -induced Ca²⁺ release (CICR) and whole cardiomyocyte electrophysiology incorporates stochastic simulation of individual LC- and RyR-channels, spatially detailed concentration dynamics in dyadic clefts, rabbit membrane potential dynamics, and a system of partial differential equations for myoplasmic and lumenal free Ca²⁺ and Ca²⁺-binding molecules in the bulk of the cell. We developed a novel computational approach to resolve the concentration gradients from dyadic space to cell level by using a quasistatic approximation within the dyad and finite element methods for integrating the partial differential equations. We show whole cell Ca²⁺-concentration profiles using three previously published RyR-channel Markov schemes.

Examination of the Effects of Heterogeneous Organization of RyR Clusters, Myofibrils and Mitochondria on Ca2+ Release Patterns in Cardiomyocytes

Spatio-temporal dynamics of intracellular calcium, [Ca²⁺]_i, regulate the contractile function of cardiac muscle cells. Measuring [Ca²⁺]_i flux is central to the study of mechanisms that underlie both normal cardiac function and calcium-dependent etiologies in heart disease. However, current imaging techniques are limited in the spatial resolution to which changes in [Ca²⁺]_i can be detected. Using spatial point process statistics techniques we developed a novel method to simulate the spatial distribution of RyR clusters, which act as the major mediators of contractile Ca²⁺ release, upon a physiologically-realistic cellular landscape composed of tightly-packed mitochondria and myofibrils. We applied this method to computationally combine confocal-scale (~ 200 nm) data of RyR clusters with 3D electron microscopy data (~ 30 nm) of myofibrils and mitochondria, both collected from adult rat left ventricular myocytes. Using this hybrid-scale spatial model, we simulated reaction-diffusion of [Ca²⁺]_i during the rising phase of the transient (first 30 ms after initiation). At 30 ms, the average peak of the simulated [Ca²⁺]_i transient and of the simulated fluorescence intensity signal, F/F₀, reached values similar to that found in the literature ([Ca²⁺]_i ≈1 μM; F/F₀≈5.5). However, our model predicted the variation in [Ca²⁺]_i to be between 0.3 and 12.7 μM (~3 to 100 fold from resting value of 0.1 μM) and the corresponding F/F₀ signal ranging from 3 to 9.5. We demonstrate in this study that: (i) heterogeneities in the [Ca²⁺]_i transient are due not only to heterogeneous distribution and clustering of mitochondria; (ii) but also to heterogeneous local densities of RyR clusters. Further, we show that: (iii) these structure-induced heterogeneities in [Ca²⁺]_i can appear in line scan data. Finally, using our unique method for generating RyR cluster distributions, we demonstrate the robustness in the [Ca²⁺]_i transient to differences in RyR cluster distributions measured between rat and human cardiomyocytes.

Calcium (Ca²⁺) acts as a signal for many functions in the heart cell, from its primary role in triggering contractions during the heartbeat to acting as a signal for cell growth. Cellular function is tightly coupled to its ultra-structural organization. Spatially-realistic and biophysics-based computational models can provide quantitative insights into structure-function relationships in Ca²⁺ signaling. We developed a novel computational model of a rat ventricular myocyte that integrates structural information from confocal and electron microscopy datasets that were independently acquired and includes: myofibrils (protein complexes that contract during the heartbeat), mitochondria (organelles that provide energy for contraction), and ryanodine receptors (RyR, ion channels that release the Ca²⁺ required to trigger myofibril contraction from intracellular stores). Using this model, we examined [Ca²⁺]_i dynamics throughout the cell cross-section at a much higher resolution than previously possible. We estimated the size of structural maladaptation that would cause disease-related alterations in [Ca²⁺]_i dynamics. Using our methods for data integration, we also tested whether reducing the density of RyRs in human cardiomyocytes (~40% relative to rat) would have a significant effect on [Ca²⁺]_i. We found that Ca²⁺ release patterns between the two species are similar, suggesting Ca²⁺ dynamics are robust to variations in cell ultrastructure.

Rapid changes in NADH and flavin autofluorescence in rat cardiac trabeculae reveal large mitochondrial complex II reserve capacity

KEY POINTS

A photometry-based technique was developed to measure nicotinamide adenine dinucleotide (NADH) and flavin adenine dinucleotide (FAD) autofluorescence and contractile properties simultaneously in intact rat trabeculae at a high time resolution. This provides insight into the function of mitochondrial complex I and II. Maximal complex I and complex II activities were determined in saponin-permeabilized right ventricular tissue by respirometry. In trabeculae, complex II function was considerably smaller than the maximal complex II activity, suggesting large complex II reserve capacity. Up-down asymmetry in NADH and FAD kinetics suggests a complex interaction between mitochondrial and contractile function. These data show that simultaneous measurement of contractile properties and NADH and FAD kinetics in cardiac trabeculae provides a mean to study the differences in complex I and II function in intact preparations in health and disease.

ABSTRACT

The functional properties of cardiac mitochondria in intact preparations have been mainly studied by measurements of nicotinamide adenine dinucleotide (NADH) autofluorescence, which reflects mitochondrial complex I function. To assess complex II function, we extended this method by measuring flavin adenine dinucleotide (FAD)-related autofluorescence in electrically stimulated cardiac trabeculae isolated from the right ventricle from the rat at 27°C. NADH and FAD autofluorescence and tension responses were measured when stimulation frequency was increased from 0.5 Hz to 1, 2 or 3 Hz for 3 min, and thereafter decreased to 0.5 Hz. Maximal complex I and complex II activity in vitro were determined in saponin-permeabilized right ventricular tissue by respirometry. NADH responses upon an increase in stimulation frequency showed a rapid decline, followed by a slow recovery towards the initial level. FAD responses followed a similar time course, but in the opposite direction. The amplitudes of early rapid changes in the NADH and FAD concentration correlated well with the change in tension time integral per second (R(2)  = 0.833 and 0.660 for NADH and FAD, respectively), but with different slopes for the up and down transient. Maximal velocity of the increase in FAD concentration (16 ± 4 μm s(-1) ), measured upon an increase in stimulation frequency from 0.5 to 3 Hz was considerably smaller than that of the decrease in NADH (78 ± 13 μm s(-1) ). The respiration measurements indicated that the maximal velocity of NADH utilization (143 ± 14 μm s(-1) ) was 2 times smaller than that of FADH2 (291 ± 19 μm s(-1) ). This indicates that in cardiac mitochondria considerable complex II activity reserve is present.

Computational modeling of subcellular transport and signaling

Numerous signaling processes in the cell are controlled in microdomains that are defined by cellular structures ranging from nm to μm in size. Recent improvements in microscopy enable the resolution and reconstruction of these micro domains, while new computational methods provide the means to elucidate their functional roles. Collectively these tools allow for a biophysical understanding of the cellular environment and its pathological progression in disease. Here we review recent advancements in microscopy, and subcellular modeling on the basis of reconstructed geometries, with a special focus on signaling microdomains that are important for the excitation contraction coupling in cardiac myocytes.

Sensitivity of rabbit ventricular action potential and Ca²⁺ dynamics to small variations in membrane currents and ion diffusion coefficients

Little is known about how small variations in ionic currents and Ca²⁺ and Na⁺ diffusion coefficients impact action potential and Ca²⁺ dynamics in rabbit ventricular myocytes. We applied sensitivity analysis to quantify the sensitivity of Shannon et al. model (Biophys. J., 2004) to 5%-10% changes in currents conductance, channels distribution, and ion diffusion in rabbit ventricular cells. We found that action potential duration and Ca²⁺ peaks are highly sensitive to 10% increase in L-type Ca²⁺ current; moderately influenced by 10% increase in Na⁺-Ca²⁺ exchanger, Na⁺-K⁺ pump, rapid delayed and slow transient outward K⁺ currents, and Cl⁻ background current; insensitive to 10% increases in all other ionic currents and sarcoplasmic reticulum Ca²⁺ fluxes. Cell electrical activity is strongly affected by 5% shift of L-type Ca²⁺ channels and Na⁺-Ca²⁺ exchanger in between junctional and submembrane spaces while Ca²⁺-activated Cl⁻-channel redistribution has the modest effect. Small changes in submembrane and cytosolic diffusion coefficients for Ca²⁺, but not in Na⁺ transfer, may alter notably myocyte contraction. Our studies highlight the need for more precise measurements and further extending and testing of the Shannon et al. model. Our results demonstrate usefulness of sensitivity analysis to identify specific knowledge gaps and controversies related to ventricular cell electrophysiology and Ca²⁺ signaling.

Mitochondrial colocalization with Ca2+ release sites is crucial to cardiac metabolism

In cardiomyocyte subcellular structures, colocalization of mitochondria with Ca2+ release sites is implicated in regulation of cardiac energetics by facilitating Ca2+ influx into mitochondria to modulate the tricarboxylic acid (TCA) cycle. However, current experimental techniques limit detailed examination of this regulatory mechanism. Earlier, we developed a three-dimensional (3D) finite-element cardiomyocyte model featuring a subcellular structure that integrates excitation-contraction coupling and energy metabolism. Here, using this model, we examined the influence of distance between mitochondria and Ca2+ release sites by comparing a normal (50-nm) distance model and a large (200-nm) distance model (LD). The influence of distance was minimal under a low pacing rate (0.25 Hz), but under a higher pacing rate (2 Hz), lower levels of mitochondrial Ca2+ and NADH, elevated phosphate, and suppressed force generation became apparent in the LD model. Such differences became greater when functional impairments (reduced TCA cycle activity, uncoupling effect, and failing excitation-contraction coupling) were additionally imposed. We concluded that juxtaposition of the mitochondria and the Ca2+ release sites is crucial for rapid signal transmission to maintain cardiac-energy balance. The idealized 3D model of cardiac excitation-contraction and metabolism is a powerful tool to study cardiac energetics.

Modeling for cardiac excitation propagation based on the Nernst-Planck equation and homogenization

The bidomain model is a commonly used mathematical model of the electrical properties of the cardiac muscle that takes into account the anisotropy of both the intracellular and extracellular spaces. However, the equations contain self-contradiction such that the update of ion concentrations does not consider intracellular or extracellular ion movements due to the gradient of electric potential and the membrane charge as capacitive currents in spite of the fact that those currents are taken into account in forming Kirchhoff's first law. To overcome this problem, we start with the Nernst-Planck equation, the ionic conservation law, and the electroneutrality condition at the cellular level, and by introducing a homogenization method and assuming uniformity of variables at the microscopic scale, we derive rational bidomain equations at the macroscopic level.

Decoding myocardial Ca²⁺ signals across multiple spatial scales: a role for sensitivity analysis

Numerous studies have employed mathematical modeling to quantitatively understand release of Ca(2+) from the sarcoplasmic reticulum (SR) in the heart. Models have been used to investigate physiologically important phenomena such as triggering of SR Ca(2+) release by Ca(2+) entry across the cell membrane and spontaneous leak of Ca(2+) from the SR in quiescent heart cells. In this review we summarize studies that have modeled myocardial Ca(2+) at different spatial scales: the sub-cellular level, the cellular level, and the multicellular level. We discuss each category of models from the standpoint of parameter sensitivity analysis, a common simulation procedure that can generate quantitative, comprehensive predictions about how changes in conditions influence model output. We propose that this is a useful perspective for conceptualizing models, in part because a sensitivity analysis requires the investigator to define the relevant parameters and model outputs. This procedure therefore helps to illustrate the capabilities and limitations of each model. We further suggest that in future studies, sensitivity analyses will aid in simplifying complex models and in suggesting experiments to differentiate between competing models built with different assumptions. We conclude with a discussion of unresolved questions that are likely to be addressed over the next several years.

Modeling effects of L-type ca(2+) current and na(+)-ca(2+) exchanger on ca(2+) trigger flux in rabbit myocytes with realistic T-tubule geometries

The transverse tubular system of rabbit ventricular myocytes consists of cell membrane invaginations (t-tubules) that are essential for efficient cardiac excitation-contraction coupling. In this study, we investigate how t-tubule micro-anatomy, L-type Ca²⁺ channel (LCC) clustering, and allosteric activation of Na⁺/Ca²⁺ exchanger by L-type Ca²⁺ current affects intracellular Ca²⁺ dynamics. Our model includes a realistic 3D geometry of a single t-tubule and its surrounding half-sarcomeres for rabbit ventricular myocytes. The effects of spatially distributed membrane ion-transporters (LCC, Na⁺/Ca²⁺ exchanger, sarcolemmal Ca²⁺ pump, and sarcolemmal Ca²⁺ leak), and stationary and mobile Ca²⁺ buffers (troponin C, ATP, calmodulin, and Fluo-3) are also considered. We used a coupled reaction-diffusion system to describe the spatio-temporal concentration profiles of free and buffered intracellular Ca²⁺. We obtained parameters from voltage-clamp protocols of L-type Ca²⁺ current and line-scan recordings of Ca²⁺ concentration profiles in rabbit cells, in which the sarcoplasmic reticulum is disabled. Our model results agree with experimental measurements of global Ca²⁺ transient in myocytes loaded with 50 μM Fluo-3. We found that local Ca²⁺ concentrations within the cytosol and sub-sarcolemma, as well as the local trigger fluxes of Ca²⁺ crossing the cell membrane, are sensitive to details of t-tubule micro-structure and membrane Ca²⁺ flux distribution. The model additionally predicts that local Ca²⁺ trigger fluxes are at least threefold to eightfold higher than the whole-cell Ca²⁺ trigger flux. We found also that the activation of allosteric Ca²⁺-binding sites on the Na⁺/Ca²⁺ exchanger could provide a mechanism for regulating global and local Ca²⁺ trigger fluxes in vivo. Our studies indicate that improved structural and functional models could improve our understanding of the contributions of L-type and Na⁺/Ca²⁺ exchanger fluxes to intracellular Ca²⁺ dynamics.

Toward an integrative computational model of the Guinea pig cardiac myocyte

The local control theory of excitation-contraction (EC) coupling asserts that regulation of calcium (Ca²⁺) release occurs at the nanodomain level, where openings of single L-type Ca²⁺ channels (LCCs) trigger openings of small clusters of ryanodine receptors (RyRs) co-localized within the dyad. A consequence of local control is that the whole-cell Ca²⁺ transient is a smooth continuous function of influx of Ca²⁺ through LCCs. While this so-called graded release property has been known for some time, its functional importance to the integrated behavior of the cardiac ventricular myocyte has not been fully appreciated. We previously formulated a biophysically based model, in which LCCs and RyRs interact via a coarse-grained representation of the dyadic space. The model captures key features of local control using a low-dimensional system of ordinary differential equations. Voltage-dependent gain and graded Ca²⁺ release are emergent properties of this model by virtue of the fact that model formulation is closely based on the sub-cellular basis of local control. In this current work, we have incorporated this graded release model into a prior model of guinea pig ventricular myocyte electrophysiology, metabolism, and isometric force production. The resulting integrative model predicts the experimentally observed causal relationship between action potential (AP) shape and timing of Ca²⁺ and force transients, a relationship that is not explained by models lacking the graded release property. Model results suggest that even relatively subtle changes in AP morphology that may result, for example, from remodeling of membrane transporter expression in disease or spatial variation in cell properties, may have major impact on the temporal waveform of Ca²⁺ transients, thus influencing tissue level electromechanical function.