Computationally efficient simulation of calcium signaling in cardiomyocytes

The heart pumps blood by the orchestrated contraction of ∼ 3 billion muscle cells called cardiomyocytes. Calcium ions provide an essential link between electrical stimulation of the cardiomyocyte membrane and their contraction as well as sustaining normal heart rhythm. Release of calcium from intracellular stores can be transient and localised near release sites (observed as Calcium sparks by fluorescent microscopy) or can occur as a propagating wave supported by ∼ 20000 release sites throughout the cell. The latter in single cells are recognised as a hallmark of arrhythmia susceptible whole hearts. Development of Ca2+ wave models that incorporated the known properties of calcium release sites has proved challenging and no model has yet achieved this. Moreover, these models are extremely computationally intensive because accurate simulations are usually required at spatial and temporal scales covering several orders of magnitude. The model developed here is composed of more than two million stiff and nonlinear systems with discrete stochastic state transitions. Each system is composed of 62 state equations. The methods used to simulate this model include diagonally iterated fully implicit Runge-Kutta integration with both L-stability and step-size control, constrained Newton’s root finding method with automatic preconditioner generation, and preconditioned Conjugate Residual Squared. These methods are implemented in a parallelised form with load balancing to distribute computation across multiple GPU and CPU computing systems. The simulated Ca2+ activity is similar to that observed in an isolated cardiomyocyte. Ca2+ spark rising time is 10 − 20 ms and falls to baseline within 50 − 100 ms. Waves occur with an approximate period of ∼ 3 seconds, amplitudes rising to 2−10 μm, and propagate with speeds around ∼ 100 μms−1. This model presents a new hypothesis for wave initiation, where calcium embers remain open for orders of magnitude longer than sparks, to supply the increased Ca2+ observed in waves. Sensitivity analysis was performed on the model, using a one-factor-at-a-time method to determine a partial derivative of each model output with respect to each model input. The key wave properties used as output variables were the peak cytoplasmic [Ca2+] during a wave, the minimum cytoplasmic [Ca2+] between waves and the time between wave peaks. The model indicated that wave period was strongly influenced by 1) RyR2 activity, 2) the degree of Ca2+ buffering by CSQ and SERCa rate, 3) rate of entry of Ca2+ into permeabilised cells and 4) rate of diffusion between the LSR and TSR. The last two factors have not been considered in experimental studies and may account for cell-to-cell variations in wave period.