Bifurcations and Proarrhythmic Behaviors in Cardiac Electrical Excitations

Abstract

The heart is a hierarchical dynamic system consisting of molecules, cells, and tissues, and acts as a pump for blood circulation. The pumping function depends critically on the preceding electrical activity, and disturbances in the pattern of excitation propagation lead to cardiac arrhythmia and pump failure. Excitation phenomena in cardiomyocytes have been modeled as a nonlinear dynamical system. Because of the nonlinearity of excitation phenomena, the system dynamics could be complex, and various analyses have been performed to understand the complex dynamics. Understanding the mechanisms underlying proarrhythmic responses in the heart is crucial for developing new ways to prevent and control cardiac arrhythmias and resulting contractile dysfunction. When the heart changes to a pathological state over time, the action potential (AP) in cardiomyocytes may also change to a different state in shape and duration, often undergoing a qualitative change in behavior. Such a dynamic change is called bifurcation. In this review, we first summarize the contribution of ion channels and transporters to AP formation and our knowledge of ion-transport molecules, then briefly describe bifurcation theory for nonlinear dynamical systems, and finally detail its recent progress, focusing on the research that attempts to understand the developing mechanisms of abnormal excitations in cardiomyocytes from the perspective of bifurcation phenomena.

Keywords: action potential; afterdepolarizations; bifurcation theory; cardiac arrhythmias; nonlinear dynamical system.

Figure 1

Schematic representation of a typical action potential (AP) in the ventricular myocyte and the respective ion channel currents that contribute to AP formation. It also indicates the major molecules that constitute ion channels and transporters (for the meanings of individual symbols, see text).

Figure 2

Relationships between membrane potential changes in response to depolarization and repolarization in a normal AP without EAD (A) and an abnormal AP with EADs (B) and voltage-dependent activation/inactivation of the L-type Ca²⁺ channel. Excessive APD prolongation due to slow repolarization leads to a long-lasting stay of the membrane potential in the L-type Ca²⁺ channel current (I_CaL) window current region (blue areas), resulting in large reactivation of I_CaL. SS-act: steady-state activation curve (black lines); SS-inact: steady-state inactivation curve (red lines); Depo: depolarization; Repo: repolarization.

Figure 3

Schematic diagram showing how the dynamics of individual ion channels, transporters, [Ca²⁺]_I, and [Na⁺]_i contribute to EAD generation and regulation. The diagram depicts the major functional components (I_Ks, I_Kr, I_CaL, I_NCX, and I_NKA), factors involved in Ca²⁺- and Na⁺-handling, heart rate (HR), and their interactions related to EAD generation and regulation. The upward and downward arrows represent an increase (or enhancements) and decrease (or attenuations) in each factor, respectively. A decrease in HR (bradycardia) decreases Na⁺ influx via the Na⁺ channel activation, resulting in a decrease in [Na⁺]_i. This reduction in [Na⁺]_i facilitates EAD generation through a decrease in outward I_NKA during AP phase 2. On the other hand, a decrease in HR prolongs AP phase 4, i.e., diastolic interval (DI), and the prolongation of DI increases the amount of Na⁺ influx through NCX and thus causes [Na⁺]_i elevation. This [Na⁺]_i elevation increases Na⁺ efflux through NAK and outward I_NAK during AP phase 2. This may counteract the decrease in [Na⁺]_i due to reduced I_Na and the resulting reduction in I_NKA (modified from Figure 11 in [85]). For details of other depicted interactions that affect EAD formation, see text.

Figure 4

Schematic diagram of the structure and electrical properties of a ventricular myocyte. (A): Ion channel and transporter molecules embedded in the cell membrane, sarcoplasmic reticulum (Ca²⁺ cycling), and Ca²⁺-binding molecules related to excitation and contraction. (B): An equivalent circuit model of the cell membrane, composed of the ion channel conductance (g), membrane capacitance (C_m), and electromotive force (E) representing ion concentration gradients. The current (I) flowing through each voltage-gated ion channel is represented by the product of the time-dependent variable conductance and driving force as the difference between the membrane potential (V_m) and individual reversal potential. g_x: maximum conductance of ion channels x, for x = Na⁺, Ca²⁺, K⁺, etc.; I_Na: fast sodium channel current; I_to: transient outward K⁺ channel current; I_Kr and I_Ks: fast and slow components, respectively, of delayed rectifier K⁺ channel currents; I_CaL: L-type Ca²⁺ channel current; I_CaT: T-type Ca²⁺ channel current; I_K1: inward rectifier K⁺ channel current; I_NKA: Na⁺/K⁺ ATPase current; I_NCX: Na⁺/Ca²⁺ exchanger current; I_p(Ca): Ca²⁺ pump current in the sarcolemmal membrane; J_rel: SR Ca²⁺ release flux by Ca²⁺-induced Ca²⁺-release; J_up: Ca²⁺ uptake flux via Ca²⁺ pump (SERCA) in the SR; J_leak: Ca²⁺ leakage flux from the SR; JSR: junctional SR; NSR: network SR; CMDN: calmodulin; TRPN: troponin; SS: a subspace of the myoplasm. For details of other symbols, refer to references [99,114,115,116,117].

Figure 5

Examples of a dynamical system with continuous time and with discrete time. (A): An action potential (AP) train evoked in a cardiomyocyte with stimuli (I_stim) applied at a cycle length of CL, which is a typical example of the response of a dynamical system with continuous time. AP duration of the AP evoked by the ith stimulus is denoted as APD_i for i = 0, 1, 2, …, n, …, k − 1, k, …. (B): A point sequence consisting of APD values determined for the AP train and plotted on the (APD_n, APD_n+1)-plane. If a mapping P (Poincaré map) is obtained to represent the relationship between a point (APD_n) and the next point (APD_n+1), i.e., the dynamics of the point sequence, then the differential system with continuous time is transformed into a difference system with discrete time. The dynamics of such the difference system with discrete time can then be studied by projecting the point sequence dynamics into the state space, e.g., (APD_n, APD_n+1)-plane. Thereby, it is possible to examine the dynamics of the original differential system with continuous time more efficiently. In general, the Poincaré map is difficult to obtain analytically and is mostly obtained numerically. For examples of experimental and numerical methods for obtaining Poincaré maps, see [124,125].

Figure 6

Hopf bifurcation. Examples of the dynamic responses to changes in the parameter λ in the 2-dimensional (x, y)-state space (A) and 3-dimensional (x, y, λ)-space (B). Schematic examples of membrane potential changes at the parameter λ₀, λ₁, λ₂ in the 3-dimensional (x, y, λ)-space of panel (A,C). Changing the system parameter λ, the stability of an equilibrium point (EP) changes via a Hopf bifurcation. Before and after the occurrence of a Hopf bifurcation (λ*), a stable EP becomes unstable with the emergence of a limit-cycle (LC) oscillation.

Figure 7

Saddle-node (SN) bifurcation. Schematic diagrams of the SN bifurcations of equilibrium (A) and fixed (B) points. As simple examples, in the dynamical systems dx/dt = f(x, λ) (A) and x_n+1 = P(x_n, λ) (B), the functions f(x) and P(x_n) are varied up and/or down by changing system parameter λ. The SN bifurcation of equilibrium points occurs when the curve of f(x) touches the x-axis (dx/dt = 0 axis). On the other hand, the SN bifurcation of fixed points occurs when the curve of P(x_n) touches the diagonal line (x_n+1 = x_n) at a value λ^* of parameters. After that, a node and a saddle point emerge or disappear. A schematic diagram of Arnold’s tongue structure is also shown for a non-autonomous system with stimuli of various amplitudes and frequencies (C). Colored regions indicate parameter regions in which a periodic oscillation observed in an autonomous system can be entrained by external periodic stimuli, e.g., the orange region represents harmonic synchronized oscillation observed in the non-autonomous system. In general, the harmonic synchronized region and even the sub-harmonic synchronous region (e.g., blue and green regions) are divided in the parameter space by SN bifurcation sets. Oscillations that are asynchronous to the periodic stimuli, such as quasi-periodic oscillations, appear when set to parameters outside the synchronization region. f₀: an intrinsic frequency of the limit-cycle oscillation in the autonomous system.

Figure 8

Bifurcation phenomena associated with stable and unstable manifolds of an equilibrium point. (A) A schema of stable and unstable manifolds of the equilibrium point x₀ in the state space. W^s(x₀), stable manifold; W^u(x₀), unstable manifold. (B) Homoclinic bifurcation. When the parameter λ becomes the value of Hom, the stable and unstable manifolds of the equilibrium point x₀ coincide, giving rise to a closed-loop orbit, which is referred to as “homoclinic orbit”. This closed-loop orbit is vulnerable to changes in the parameters and breaks down immediately. (C) Saddle-node invariant cycle bifurcation. At the bifurcation parameter value λ*, saddle-node and homoclinic bifurcations occur simultaneously.

Figure 9

Schematic diagrams of period-doubling and Neimark–Sacker bifurcations. (A): A typical example of the period-doubling (PD) bifurcation in non-autonomous systems is AP alternans (i). As the cycle length (CL) becomes shorter, the stable fixed point on the (APD_n, APD_n+1)-plane becomes unstable at the period-doubling bifurcation point (λ*), and a pair of periodic points is generated (ii). (B): The change in a limit-cycle (LC) oscillation through the Neimark–Sacker (NS) bifurcation in autonomous systems. The NS bifurcation occurs at λ* as the parameter λ changes. Then, the LC oscillation becomes unstable and a torus-like oscillatory response, so-called “quasi-periodic oscillation”, occurs around the unstable LC oscillation. When the quasi-periodic oscillation is discretized as a point sequence via the Poincaré mapping, a closed curve, which is referred to as an “invariant closed curve (ICC)”, appears around the unstable fixed point that reflects the unstable LC oscillation.

Figure 10

Bifurcation phenomena and action potential behaviors observed in the Kurata model [112] of a human ventricular myocyte. (A) A two-parameter bifurcation diagram on the (G_Kr, G_Ks)-parameter plane (top) obtained using MATCONT [154] and a phase diagram (bottom) obtained by AP simulations. The maximum conductances of the rapid (G_Kr) and slow (G_Ks) components of delayed rectifier K⁺ channel currents are expressed as normalized values, i.e., ratios to the control values. In the two-parameter bifurcation diagram (top), the symbols H and PD represent the loci of parameter sets that cause the Hopf bifurcation of an equilibrium point and period-doubling bifurcation of a limit cycle (LC), respectively. The gray region indicates the area of parameters in which a stable LC can be observed. On the other hand, in the phase diagram (bottom), the orange and white regions represent parameter regions in which an AP with and without EAD, respectively, can be observed in the Kurata model. The gray point with the symbol N indicates the control condition with the normal G_Kr and G_Ks. (B) A merged phase diagram. (C) A one-parameter bifurcation diagram of APD at 90% repolarization (APD₉₀) in each AP response as a function of G_Kr; see the blue arrow in panel (B), which indicates the change in G_Kr with the fixed normal G_Ks (1). The solid and dashed lines in (C) represent stable and unstable AP responses, respectively. SN: saddle-node bifurcation; EAD1–3: AP with 1–3 EAD(s); PD: period-doubling bifurcation. NS: Neimark–Sacker bifurcation. (D) An example of tetra-stable AP dynamics in the Kurata model at 0.39 × G_Kr with 0.70 × G_Ks. Colored and grey lines indicate the steady-state and transient responses, respectively. Dots indicate the application of current pulses. Pacing cycle length = 2 s. Each panel was modified from [85,130].

Figure 11

Examples of the slow–fast decomposition analysis. AP simulations and slow–fast decomposition analyses were performed using the Kurata model [112] with G_Kr = 100% (A) and G_Kr = 0.5 × control (B). (Left) The membrane potential (V_m) and the time course of I_Ks open probability (n²) in the full system. (Center) One-parameter bifurcation diagrams of the fast subsystem, depicting quasi-equilibrium potentials (qEQ_1–3) and the potential extrema of quasi limit cycles (qLCs) as a function of the slow variable n². The trajectory of the full system projected onto the (V_m, n²)-plane is superimposed on the one-parameter bifurcation diagram of the fast subsystem. Solid and dashed thin lines of qEQ_1–3 indicate stable and unstable equilibrium potentials, respectively. The black and gray thick lines are stable and unstable qLCs, respectively. H: Hopf bifurcation; SN: saddle-node bifurcation of limit cycles; PD: period-doubling bifurcation of limit cycles; hom: homoclinic bifurcation. (Right) Schematic diagrams representing the relationships between the dynamic behavior projected onto the 3-dimensional state space of the full system and the bifurcation structure of the fast subsystem. Modified from [197].