Mathematical Modeling in Computational Biology:
Focus on Excitable Media and Cardiac Electrophysiology
Modeling normal human ventricular cell and tissue electrophysiology. I have been involved in developing a minimal model of human ventricular cell electrophysiology based closely on published data for three anatomical regions of the ventricular wall (epicardial, endocardial, and midmyocardial) that experience gradients in certain electrophysiological properties. Previously, three other models human ventricular cell physiology had been published, none of which adequately matched newly obtained experimental data. In addition, the other existing models do not provide accurate descriptions of the regional gradients across the ventricular wall. Using the minimal model, my colleagues and I have been able to explore the contribution of these regional gradients to the morphological characteristics of the electrocardiogram (ECG) using a realistic three-dimensional ventricular anatomical model and specifically to match experimental data obtained using a ventricular slab similar to experimental setups (see Fig. 1).
Figure 1. (A) Action potential duration (APD) as a function of diastolic interval (DI) obtained from experiments (dots) and using the minimal model (solid) for epicardial (blue), endocardial (green), and midmyocardial (red) human ventricular cells, with corresponding simulated action potentials. (B) Anatomical ventricular slab matching experimental slab conditions including regional variations in electrophysiological properties in the epicardium (blue), midmyocardium (red), and endocardium (green), with corresponding action potentials. (C) Simulated ECG in the anatomical structure shown in (B) obtained using the minimal model.
Modeling Brugada Syndrome. I have been developing a model of diseased physiology reflecting differences believed to be present with a type of genetic disease called the Brugada syndrome, which can cause sudden cardiac death without any underlying structural defects. In many cases, Brugada Syndrome has been shown to occur in the presence sodium ion channel mutation that results in loss of function that decreases sodium current. Most experimental data on Brugada syndrome has been obtained using dogs, so I initially developed a model to reproduce important properties of these experiments. More recently, Brugada syndrome data in humans has become available, and accordingly I have been working to implement these characteristics within the context of the human ventricular cell model. Among other things, I have been able to demonstrate the importance of regional variations in cell properties in inducing reentry reentrant arrhythmias (see Fig. 2) and have identified the conditions necessary for a reentrant arrhythmia (spiral wave) to occur, including slow pacing rates and appropriate sodium channel recovery from inactivation. Following this preliminary work, it is still necessary to investigate how Brugada Syndrome arrhythmias develop in the full anatomical context of the human ventricles, complete with anatomical complexity and regional gradients in ion channel properties.
Figure 2. Arrhythmia development in an anatomically realistic three-dimensional model of the canine right ventricular free wall due to the presence of Brugada syndrome, which causes various regions of the epicardium to delay the action potential dome (A, lower left) or to lose the dome entirely (A, upper left) compared to the normal epicardial action potential (A, right). Because repolarization time varies by region (B-C), the late dome region is able to generate a new wave that propagates into the fully repolarized loss-of-dome region above it (D-E). This initiates a reentrant arrhythmia (F-G). Electrical wave fronts are shown in red, with quiescent tissue in dark blue.
Movie of Brugada arrhythmia development in the RV wedge.
Movies of Brugada arrhythmia development in 2D preparations with the same distribution of normal, late dome, and lost dome regions as in Fig. 2A. Stimulation from the left and from below result in arrhythmia initiation, while stimulation from the right does not induce an arrhythmia.
Modeling rabbit ventricular cell and tissue electrophysiology. To understand arrhythmia dynamics in rabbit ventricles, mathematical it was necessary to construct models of rabbit ventricular action potentials to match electrophysiological properties in the presence of two different motion uncoupling agents used for optical mapping experiments: diacetyl monoxime and cytochalasin-D. These pharmacological agents are used to prevent ventricular contractions that would interfere with image recording. However, they also affect electrophysiological properties. Starting from action potential recordings and restitution information recorded in tissue using both optical and microelectrode techniques, two models were developed and then simulated in two- and three-dimensional structures. In two-dimensional sheets, spiral waves initiated using either model remained stable with linear core trajectories. However, using an anatomically realistic representation of rabbit ventricles, the dynamics differed, matching what was recorded experimentally: with diacetyl monoxime, only a single spiral wave corresponding to ventricular tachycardia was observed (Fig. 3A), while with cytochalasin-D, an initiated spiral wave broke apart into multiple waves (Fig. 3B), corresponding to ventricular fibrillation. Analysis of these results showed that the cytochalasin-D model produced more complex dynamics because, unlike the diacetyl-monoxime model, spirals waves experienced tight turns, making them vulnerable to several breakup mechanisms. Specifically, fibrillation arose from an instability of the core induced by the anisotropic fiber rotation and by dynamically created regions of dispersion of refractoriness produced by wave front-back interaction due to the large size of the wavelength relative to the ventricular size.
Figure 3. Simulated ventricular tachycardia (A) and fibrillation (B) using an anatomically realistic model of the rabbit ventricles including fiber orientation information. Electrical potential is shown, with voltage ranging from –85 mV (blue) to 25 mV (red). Anterior, posterior and apical views are shown, along with a semi-transparent view showing the locations of transmural reentries in red against the endocardial surfaces in green.
Movies of ventricular tachycardia in the presence of
diacetyl monoxime and ventricular fibrillation
in the presence of cytochalasin-D.
Modeling normal canine ventricular cell and tissue electrophysiology. Over the last decade, transmural gradients in electrophysiological properties have been described. We have used published experimental data in the canine heart to develop models of the three dominant types of cells: epicardial, endocardial, and specialized midmyocardial “M” cells. In addition, we subdivided the canine ventricular anatomical model into spatial regions inhabited by each of these cell types, as shown in Figure 4A for a representative cross-section. Then we introduced a Purkinje network to reproduce a published activation sequence observed experimentally (see Figure 4B) and simulated a series of normal activations and the accompanying ECG. We compared the results with and without the heterogeneity, and also with and without an apex-base gradient in repolarization. We found that although the apex-base gradient increased the magnitude of an ECG deflection that indicates repolarization (the T-wave), the transmural heterogeneity was necessary for T-wave polarity to be correct (upright), as shown in Figure 4C-D, in agreement with a current hypothesis of the origin of the T-wave.
Figure 4. (A) Transmural heterogeneities in the canine ventricular model. The red region shows M cells are assigned in the midmyocardium, with endocardial cells in the remaining inside tissue and epicardial cells in the remaining outside tissue. (B) Activation sequence in simulated canine ventricles. The results agree well with the experimentally obtained activation sequence of Durrer. (C) and (D) Representative pseudo-ECGs obtained in simulations. Tissue is homogeneous in (C) and includes transmural heterogeneity in (D), whose upright deflection following each spike closely resembles healthy ECGs.
Movie of the activation sequence in the canine ventricles.
Incorporating realistic intracellular calcium dynamics and contraction. Over the last several years, it has been demonstrated that voltage and calcium dynamics can operate independently, and recent optical mapping experiments have shown that the two systems sometimes desynchronize during fibrillation. Only one recently published mathematical model has been shown to allow independent behavior of the voltage and calcium systems at fast pacing rates. I have been working to couple this calcium model with a model of the cardiac action potential. Preliminary results are shown in Fig. 5, where alternans in voltage (A) are shown in concert with alternans in the intracellular calcium concentration (B).
Figure 5. Simulation of voltage oscillations known as electrical alternans (A) driven by alternans in intracellular calcium (B). Calcium alternans occur only at relatively fast rates but generally disappear at extremely fast rates.