The theoretical basis accompanied by the examples of implementation of the discussed algorithms with the help of MATLAB

Повний текст

(1)MINISTRY OF EDUCATION AND SCIENCE OF UKRAINE NATIONAL TECHNICAL UNIVERSITY OF UKRAINE “IGOR SIKORSKY KYIV POLYTECHNIC INSTITUTE”. K.Ivanko, N.Ivanushkina. Advances in Digital Processing of Low-Amplitude Components of Electrocardiosignals. Kyiv, 2017 - 1-.

(2) UDC 621.391:004.93:616.12-073.97(075.8)621.38:61(0.75.8). Reviewers: O.Kovalenko, Doctor of Medical Sciences, Professor of International Research and Training Center for Information Technologies and Systems of NASU O.Sharpan, Doctor of Technical Sciences, Professor of Igor Sikorsky Kyiv Polytechnic Institute Responsible Editor: V. Timofeev, Doctor of Technical Sciences, Professor of Igor Sikorsky Kyiv Polytechnic Institute Advances in Digital Processing of Low-Amplitude Components of Electrocardiosignals: Teaching book/K.Ivanko, N.Ivanushkina.─ Kyiv, 2017.─198 p. ISBN ______________________ The teaching book is devoted to development and research of methods and tools for non-invasive detection of subtle manifistations of heart electrical activity. Particular attention is paid to the improvement of information and algorithmic support of high resolution electrocardiography for early diagnosis of myocardial electrical instability, as well as for the evaluation of the functional state of the fetus during pregnancy examination. The theoretical basis accompanied by the examples of implementation of the discussed algorithms with the help of MATLAB. The teaching book is intended for students, graduate students, as well as specialists in the field of biomedical electronics and medical professionals. ISBN _____________________. © K.Ivanko, N.Ivanushkina, 2017. This manual has been published within the framework of the BME-ENA project under the responsibility of National Technical University of Ukraine. The BME-ENA “Biomedical Engineering Education Tempus Initiative in Eastern Neighbouring Area”, Project Number: 543904-TEMPUS-1-2013-1-GR-TEMPUSJPCR is a Joint Project within the TEMPUS IV program. This project has been funded with support from the European Commission. This publication reflects the views only of the author, and the Commission cannot be held responsible for any use which may be made of the information contained therein.. - 2-.

(3) CONTENT INTRODUCTION…………………………………………………………..…6 SECTION 1. BIOELECTRIC PROCESSES OF CELL MEMBRANE…10 1.1. Structure and functions of cell membrane………………………..10 1.2. Nature and modeling of membrane resting potential……………12 1.2.1. Nernst-Planck equation and Nernst potential………………13 1.2.2. Goldman equation……………………………………………15 1.2.3. Parallel-conductance model…………………………………17 1.3. Nature and modeling of action potential…………………………..18 1.3.1. Simulation of conductivity of ion channels …………………18 1.3.2. Simulation of action potential in cardiomyocytes………….23 1.3.3. Simulation of action potential in nerve fibers………………40 1.3.4. Modeling of the circulation of exciting impulse in the myocardium…................................................................................45 Conclusions………………………………………………………………...54 Control questions and tasks………………………………………………..54 References………………………………………………………………….56 SECTION 2. POTENTIALS. IDENTIFICATION AS. MARKERS. AND OF. ASSESSMENT CARDIAC. OF. LATE. ELECTRICAL. INSTABILITY……………………………………………………………….59 2.1. Electrical conduction system of the heart………………………….59 2.2. Arrhythmias and mechanisms of their occurrence….....................63 2.3. Circulation of excitation wave in the myocardium by the reentry mechanism…………………………………………………………64 - 3-.

(4) 2.4. Biophysical bases of genesis of atrial and ventricular late potentials…………………………………………………………………….68 2.5. Peculiarities of late potentials detection in the systems of high resolution electrography……………………………………………………71 2.6. Simulation of late potentials in ECG………………………………….80 2.7. Analysis of the subtle structure of the electrocardiogram…………83 2.7.1. Features of wavelet analysis of electrocardiosignals………84 2.7.2. Combined method for analysis of the subtle structure of electrocardiosignal for late potentials identification………………...90 Conclusions………………………………………………………………..109 Control questions and tasks………………………………………………110 References…………………………………………………………………112. SECTION 3. IDENTIFICATION AND ASSESSMENT OF T WAVE ALTERNANS. AS. A. PREDICTOR. OF. LIFE-THREATENING. VENTRICULAR TACHYARRHYTHMIAS………………………………120 3.1. T wave alternans phenomenon ………………………………….120 3.2. Detection and assessment of the level of T wave alternans by the scattergram method………………………………………………..123 3.3. Detection of T wave alternans using spectral analysis……..…132 3.4. Identification and evaluation of T wave alternans using principal component analysis………………………………………………..133 Conclusions ……………………………………………………………….139 Control questions and tasks……………………………………………..139 References…………………………………………………………………140. - 4-.

(5) SECTION 4.. ATRIAL ELECTRICAL ACTIVITY EXTRACTION. FOR. ATRIAL FIBRILLATION ASSESSMENT……………………………….142 4.1. Method of Blind Source Separation for atrial activity extraction…144 4.2. Average beat subtraction approach……………………………….151 4.3. Fibrillatory frequency tracking………………………………………157 Conclusions ……………………………………………………………….161 Control questions and tasks………………………………………………161 References…………………………………………………………………161. SECTION 5. FETAL ELECTROCARDIOGRAM EXTRACTION FROM MATERNAL ABDOMINAL SIGNALS…………………………………….164 5.1. Peculiarities of manifestation of fetal cardiac electrical activity…164 5.2. Fetal cardiocycle modeling…………………………………………170 5.3. Fetal ECG detection and analysis based on combined use of ICA and WT methods ..................................................................................174 Conclusions ………………………………………………………………189 Control questions and tasks…………………………………………….190 References………………………………………………………………..190 CONCLUSIONS………………………………………………………….195. - 5-.

(6) INTRODUCTION The development of digital processing of biomedical signals and implementation of information technology in medical practice create conditions for improvement of diagnostic methods in cardiology. The use of special equipment for electrocardiogram registration and digital signal processing techniques provide information that can not be obtained by analysis of standard electrocardiogram. The classic use of high resolution electrocardiography is noninvasive detection of markers of myocardial electrical instability to predict cardiac arrhythmias. It is believed that the first signs of electrical instability of the heart reflect the exhaustion of regulatory systems at the level of individual cells of the myocardium. Changes in energy and metabolic processes that occur in the deterioration of electrical stability of the myocardium at the cellular and subcellular levels and microstructural violations can not be detected by conventional clinical and physiological methods of functional diagnostics. These changes in the first stage of the disease often does not manifest clinically. The first part of the teaching book focuses on simulation of action potential and membrane ionic currents in cardiomyocytes and nerve fibers. It presents a mathematical model based on 6 ordinary differential equations that describe the ionic currents at the level of cardiomyocyte, transmembrane potential with action potential generation, as well as activation and inactivation functions for potassium, sodium and calcium channels of cardiomyocytes. The second part of the teaching book is devoted to development of diagnostic ECG systems of new generation, which are characterized by the use of complex mathematical transformations for estimation of - 6-.

(7) parameters. of. electrical. potentials. considering. biophysics. and. electrophysiology of myocardium. The systems of high resolution electrocardiography. (HR. ECG). allow. investigation. of. subtle. manifestations of cardiac electrical activity. Early diagnosis of cardiovascular system state, which developed on the basis of HR ECG methods, has essentially new opportunities and begins to be widely applied in clinical practice. To the tasks of investigation of subtle structure of ECG belongs detection of electrical violations of myocardial homogeneity based on registration and analysis of atrial and ventricular late potentials by the systems of high resolution electrocardiography. The presence of late potentials is associated with an increased probability of dangerous cardiac arrhythmias. Particularly difficult is the study of atrial late potentials, which occur as high frequency microbursts in terminal part of P wave and have nearly the same amplitude as noise components of cardiosignal. The amplitude of atrial late potentials is much lower than the amplitude of ventricular late potentials. Atrial late potentials reflect the presence of fragmented delayed depolarization of the atria and are markers of myocardial electrophysiological disorders that lead to supraventricular arrhythmias such as atrial fibrillation and paroxysmal atrial tachycardia. Detection of late potentials may be based on the use of a broad class of algorithms from the conventional analysis in the time and frequency domains to the complex procedures of pattern recognition. Combined method of analysis of low-amplitude components of electrocardiosignals is developed. It is intended for detection of late potentials representing pathological cardiac electrical activity and based on the creation of eigensubspaces of signals and noise at different. - 7-.

(8) combinations. of. wavelet. decomposition. and. decomposition. in. eigenvectors basis. The third part of the teaching book is dedicated to detection and evaluation of beat-to-beat T wave alternans in electrocardiogram. The task of developing of theoretically justified methods and algorithms for early prediction of potentially dangerous atrial arrhythmias is very important because of the wide spread of cardiovascular diseases and their significant contribution to the structure of mortality and disability of socially active groups. To reveal diagnostic features related to T wave alternans, variations in amplitude and shape of T waves from synthetic and clinical recordings were explored by means of time-domain analysis, Poincare plots, and principal component analysis together with clustering analysis. The proposed methods make it possible to eliminate noises and interferences in recorded electrocardiograms, and also to detect the morphology changes of T waves. The fourth part of the teaching book is focused on noninvasive research of atrial electrical activity by its extraction from the surface electrocardiogram recordings. Spectral and time-frequency analysis is performed for fibrillatory frequency detection and tracking. It is shown that atrial electrical activity and ventricular activity can be separated from real multichannel recordings with symptoms of atrial fibrillation by means of independent component analysis. To identify the source that corresponds to atrial electrical activity, the frequency spectrums and coefficients of skewness and kurtosis are determined. Also atrial activity extraction is performed by average beat subtraction method. The comparison of obtained results is carried out. Spectral analysis by Fast Fourier Transform and Multiple Signal Classification methods as well as. - 8-.

(9) time-frequency analysis are used to find and track dominant atrial fibrillation frequency. The fifth part of the teaching book demonstrates another example of the investigation of subtle structure of ECG, which consists in separation of fetal and maternal electrocardiosignals obtained by the abdominal. recording.. Despite. significant. advances. in. clinical. electrocardiography, this direction is still in early stage of development. It is caused by lack of necessary databases of fetal ECG in norm and with the presence of defects, low amplitude of fetal ECG compared to the high-amplitude maternal signal, as well as the presence of abdominal noise components in the abdominal ECG signal. The proposed method of fetal ECG separation from the abdominal maternal electrocardiosignals makes it possible to separate from the noisy mixture of signals the low-amplitude component of fetal ECG that is the independent determination. source of. of. fetal. cardiac heart. electrical. rate. activity.. variability. and. This. allows. morphological. parameters of fetal cardiac cycles for assessment of fetal condition during pregnancy. The teaching book presents the research of simulated and real electrocardiosignals. It can be useful for students, graduate students, as well as for developers of new methods and means of high resolution electrocardiography.. - 9-.

(10) SECTION 1 BIOELECTRIC PROCESSES OF CELL MEMBRANE. 1.1. Structure and functions of cell membrane Vital activity of any organism is accompanied by bioelectric processes. The appearance of electrical potentials in living cells is associated with the physicochemical properties of cell membranes. All living organisms consist of cells, which surrounded by biological membranes with thickness about 6-10 nm. Biological membranes are the limited supramolecular structures, which form the cells and their intracellular environment and provide all functions of organism. The most important function of biological membranes is regulation of ions and substances transport [1]. The structure of the membrane composed of lipids, carbohydrates and protein components (Fig. 1.1). One of the first membrane models is bilayer model of lipid membrane ("double sandwich"). In this model lipid molecules are arranged in two layers in such a way that the polar heads, which have hydrophilic properties, are directed outward, and nonpolar tails, which have hydrophobic properties, are directed into the membrane. Protein components penetrate the bilayer of the membrane. Lipid layer effectively prevents the process of free passage of ions through the membrane.. - 10-.

(11) 4 1 2 1 3. Fig. 1.1. Structure of the cell membrane: 1 ─ hydrophilic lipid head, 2 ─ of hydrophobic lipid molecules, 3 ─ protein molecule, 4 ─ carbohydrate molecules Hydrophobic tails form the internal structure of the membrane, which behaves like a dielectric with a relative dielectric constant and thickness d = 30Å. The capacity Cm of following structure is about 1µF/cm2 [1]. The functions of membrane proteins include formation of the channels, which enable exchange of ions between intra- and extracellular environments. Functional model of membrane channels (Fig.1.2) consists of the aqueous pores, through which ions have the ability to penetrate the membrane, selectivity filter, which displays property of selective permeability channels for different ions and sensor with control gates. On the base of functional model ion channels can be classified by selectivity according to the ions that pass through the channels: sodium, potassium, calcium, chloride. Also, functional model allows classification of ion channels depending on the control facilities: - independent channels; - potential-dependent channels; - 11-.

(12) - ligand-dependent channels; - jointly-controlled channels (activated by ligands and membrane potential); - stimulus-driven. channels. (mehano-sensitive,. temperature-. sensitive).. Fig. 1.2. Functional model of membrane channels [1] 1.2. Nature and modeling of membrane resting potential. One of the most characteristic features of living cells at rest is a potential difference between the inside part of cell and environment, which is known as resting potential. Registration of the resting potential is carried out by inserting the electrode into the cell (registered voltage is approximately -70 mV).. - 12-.

(13) 1.2.1. Nernst-Planck equation and Nernst potential In all excitable cells intracellular potassium concentration is much higher than the potassium concentration outside cell and extracellular concentrations of sodium and chlorine are much higher comparing to their intracellular concentrations [1]. Inequality of ions concentration in the intra- and extracellular space cause the diffusion of ions from the region with their high concentration to the region with low concentration. Due to the fact that ions have electric charges, membrane has a certain electric capacity. The charges are collected on the membrane, leading to the emergence of a potential difference. This potential difference creates in the thickness of the membrane the electric field, which generates forces acting on all the charged particles inside the membrane. Thus, the total current of ions through the membrane has diffusive and drift components:      J  J dif  J dr  ZDN  ZN ,. where Z is ion charge, D is diffusion coefficient, N is ion concentration as a function of coordinates, µ is mobility; φ is potential of electric field inside the membrane. Then, the full density of the ion current across the membrane can be found using Nernst-Planck equation:   Z N   J  ZD  N    , q T  . where φT is thermal potential, q is electron charge. Nernst-Planck equation is used to obtain quantitative estimation of the membrane potential. Under the equilibrium conditions the force of the - 13-.

(14) electric field exactly compensates diffusion force, and full density of the ion current across the membrane equals zero:   Z Nc   J  Dc Z  Nc     0 , q  T  . where Dc is diffusion coefficient of cations, Nc is concentration of cations. From this equation we can find potential difference with regard to assumption of integration on the membrane thickness from the intracellular space i to the extracellular space e. For biological cell membrane potential is determined as the difference between internal and external potentials. Thus, the potential difference on the membrane is Vm: Vm  i  e  T. q NCe ln Z NCi. The resulting value of the membrane potential is called the Nernst potential. It can be regarded as an electrical measure, which counterbalances the diffusion force that arises from the difference in concentrations on opposite sides of the permeable membrane. Nernst potential for cations can be defined as: VC  T. where. NCe. q NCe ln Z NCi. is extracellular concentration of cations, NCi. intracellular concentration of cations. Nernst potential for anions can be found as:. - 14-. is.

(15) U A  T. where. q N Ai , ln Z N Ae. N Ae is extracellular concentration of anions, N Ai. is. intracellular concentration of anions. 1.2.2. Goldman equation. Membrane potential is a weighted sum of the Nernst potentials for different ions. The weights in this amount depend on on the ability of membrane to pass through itself a variety of presented ions. To take into account for the influence of membrane permeability for different ions on the membrane potential, we have to solve Nernst-Planck equation, based on membrane processes. In general, biological membrane can not be in the state of equilibrium for all ions. Nernst potentials for K+, Na+, Cl- ions under their usual concentrations have different values. Thus, none of the membrane Nernst potentials can simultaneously balance all ions. Resting potential can be found from the condition of equilibrium J = 0.. pK [K ]e  pNa [Na]e  pCl [Cl ]i   pK [K ]i  pNa [K ]i  pCl [Cl ]e   e. Vm T. 0. where [K]i and [K]e are intra- and extracellular concentrations of potassium ions, [Na]i and [Na]e. are. intra- and extracellular. concentrations of sodium ions, [Cl]I and [Cl]e are intra- and extracellular concentrations of chlorine ions, pK is permeability of potassium ions, pNa is permeability of sodium ions, pCl is permeability of chlorine ions.. - 15-.

(16) From this formula we can obtain the resting potential:.  p [K ]  pNa [Na]e  pCl [Cl ]i  Vm  T ln  K e .  pK [K ]i  pNa [Na]i  pCl [Cl ]e . This equation is called the equation of Goldman. Goldman equation shows that the membrane potential is independent from the absolute values of membrane permeability for different ions and depends only on their relationship pK : pNa : pCl. If the membrane is permeable only to one type of ions (for example, pК>>pNa and pK>>pCl), then the membrane potential equals Nernst potential of that type of ions, for which membrane is permeable (such as K+ ions). The impact of chlorine ions to creation of resting potential in comparison with the impact of potassium ions is secondary. This is due to the fact that the intracellular concentration of chlorine ions is very small and under a small influx or outflow of ions is subjected to large relative changes. Consequently, the movement of chloride ions will lead its concentration on both sides of the membrane to practically fixed ratio of the concentrations of potassium ions. Therefore, in many practical problems, it can be assumed, that the chlorine ions are in equilibrium or absent. In this case the resting potential is defined as: Vm  T ln. where. . [K ]e  [Na]e , [K ]i  [Na]i. pNa . pK. - 16-.

(17) 1.2.3. Parallel-conductance model The assumption of the presence of independent conducting ion channels allows construction of a circuitry model of membrane, known as the parallel-conductance model [1]. Each of the branches of the electrical circuit reflects the contribution of certain type of ions to the total membrane current (Fig. 1.3).. Fig. 1.3. The parallel-conductance model of membrane (IN = intracellular space, OUT = extracellular space) [1] Independent conductance channels are presented for K+, Na+ and Cl− ions. The voltage sources in the equivalent circuit simulate Nernst potentials for each type of ions: EK  T ln ENa  T ln. [K ]e ; [K ]i. [Na]e ; [Na]i. ECl  T ln. [Cl ]i . [Cl ]e. Complete membrane current can be defined as: - 17-.

(18) J  gK Vm  EK   gNa Vm - ENa   gCl Vm  ECl   Cm. At equilibrium J = 0 and. Vm . dVm . dt. dVm  0 . Then dt gNaENa  g K EK  gCl ECl . g Na  g K  gCl. This equation is known as the equation of parallel conductances. According to it, Vm is a weighted average of values EK, ENa, EСl, depending on the relative conductivities. 1.3. Nature and modeling of action potential 1.3.1. Simulation of conductivity of ion channels. A typical current through a single channel consists of discrete random jumps. According to experimental data, current fluctuations can be explained by opening and closing of the channel, which has a fixed conductance. The current switches between two discrete levels and durations of open and closed positions are described by random variables. The current, which is created by a large number of channels, corresponds to the average current of the whole cell. Let's consider relationships that describe macroscopic properties of the membrane by a total effect of a large number (N) of individual channels, among which NO channels are open and NC channels are closed (N = NO + NC). It is believed that over time some channels randomly switch from the closed position to the open and vice versa [1-3].. - 18-.

(19) . NO. Nc . If  is the speed of channels switching from closed position to open and  is speed of switching from open position into closed, then the full speed of changing the open channels is determined as: dNO  NC  NO    N  NO   NO dt. Suppose that p . Nо is the probability that the channel is open, then N. dp  (1  p)  p dt The last equation can be rewritten as. dp p  p  dt  where  . 1 is relaxation period, p   is the steady-state . value. The decision of the equation is an expression . t . p  p  ( p  p0 )e ,. where p0 is the probability that the channel is in the open position at the moment t=0. Character of the solution depends on the ratio between the values of p0 and p (Fig. 1.4).. - 19-.

(20) p p p> p0 p = p 0. p0. p < p 0. p. t. Fig. 1.4. Probability function p(t) depending on the ratio between the values of p0 and p. If p= p0, then the probability that the channel is in the open state is independent from time and corresponds to the steady state. When p>p0, function p(t) increases due to the opening of the channels over time. This process is called activation of the channels, and the function p(t) is known as activation function. When p<p0, the inactivation of the channels occurs and the channels change their state to the closed position. In this case, the function p(t) is called inactivation function. Described dependences are fair for so-called individual channels. The real ion channels are composed of several subunits, each of which can be described as a single channel. The channel is considered to be open, if all the subunits are also open. If pi is the probability that i subunit is in the open state, then for the channel, which consists of k subunits, function p(t), is defined as k. p(t )   pi (t ) . i 1. - 20-.

(21) Biological membranes contain a large number of several types of ion channels. With stimulation such electrically active tissues can generate an electrical impulse - action potential (AP). Action potential depends on the voltage applied to the stimulating electrodes and its polarity. If the stimulus amplitude is small, then stimulus artifact can be registered, which coincides in time with the stimulus and occurs through direct capacitive connection of the stimulator and recording tool. Increasing the amplitude of the stimulus we can reach a level at which an answer, potential of the cell occurs, indicating the threshold nature of the phenomenon. With increasing of the stimulus amplitude above the threshold, amplitude of action potential remains unchanged. By placing a microelectrode inside the cell, a change in intracellular potential relative to extracellular potential can be measured. At rest state, the membrane has a negative potential of 60-100 mV. The cycle of the action potential consists of fast depolarization (1) that implies increasing of the transmembrane potential, then follows the slow repolarization (2) of the membrane and return to the resting state of the membrane (3) (Fig. 1.5).. Fig.1.5. The main phases of action potential [4] - 21-.

(22) The answer of the excitable membranes does not occur until the stimulus reaches a certain level, called the threshold potential. For all stimuli, which applied to the cell with amplitude higher than the threshold, action potentials have the same amplitude. To understand the origin of the action potential we must turn to the Goldman equation. In deriving the Goldman equation a key assumption implies the existence of a steady state, in which the full density of ion flow is equal to zero J=0. At the top of the action potential. dVm  0 , and dt. therefore total current density through the membrane must also be zero (quasi-stationary state). However, at the peak value of action potential membrane permeability is different from the permeability at rest. Experimentally determined, that the peak of action potential is approaching the Nernst potential for sodium, but never exceeds it. This result is consistent with a possible increase in sodium permeability. It has been shown that good match between theory and experiment is achieved if accept pK : pNa : pCl = 1.0 : 0.01 : 0.45 at rest and pK : pNa : pCl = 1.0 : 10.0 : 0.45 at the peak of action potential. The influence of chlorine on behavior of membrane at rest is negligible. As a result, given the difference of permeabilities pK and pNa , as a first approximation we can assume that at rest V m  E K   T ln. [ K ]e ; [K ]i  [N a ] . e at the peak of action potential V m  E N a   T ln  . [ N a ] i  . - 22-.

(23) 1.3.2. Simulation of action potential in cardiomyocytes Classification of models of heart electrical activity is based on the description of bioelectric processes, which take place at different levels of organization: cell, tissue and organ. Cellular level models reflect the occurrence of action potential in cardiomyocytes. Simulation at the tissue level describes AP distribution processes in myocardium, and organlevel models allow investigation of dynamics of cardiac conduction system. Cellular-level models are based on a fundamental approach, which was proposed by Hodgkin and Huxley in the simulation of ionic currents and membrane AP for nerve fiber [2]. For simulation of action potential in cardiomyocyte, it is necessary to take into account that AP of cardiomyocyte differs from the AP of nerve fiber due to the presence of weighty contribution of calcium current. Duration of AP in axon is 1 ms, whereas duration of AP in myocardial contractile cells is 250-300 ms, that allows implementation of synchronous excitation and contraction of the heart structures for blood ejection. There are important physiological differences between nodal cells and atrial or ventricular cells of heart due to the specific differences in ion channels and mechanisms of repolarization. Therefore, action potentials in cardiomyocytes of various heart structures (sinus and atrioventricular node, atrial and ventricular cardiomyocytes) are different in shape and duration. One of the first modern models of myocardium is a model by McAllister R.E., Noble D., and Tsien R.W. [5], which includes all membrane currents that were confirmed experimentally, as well as intracellular calcium currents, which are closely linked to muscle contraction. One of the latest models is the model of the myocardium by - 23-.

(24) Luo and Rudy [6, 7]. Model of Luo and Rudy is a modern model of heart membrane, which extends the previous models by using single-channel parameter measurements. The fundamental approach of modeling, which was proposed more than fifty years ago by Hodgkin and Huxley in [2], remains unchanged. Improvement of the model was performed taking into account experimental data accumulated by many authors by introducing new features and functions for membrane currents due to the complexity of cardiac membrane and bioelectric processes in cardiomyocytes. Arrhythmias, being the most common heart diseases, are caused by violations in the formation and/or propagation of the action potential. To investigate such violations of the heart electrical activity, experimental methods and methods of mathematical simulation are performed, complementing each other. Different researchers have created a wide variety of models that simulate the appearance of the action potential in cardiomyocytes [3]. There are complex detailed models for the simulation of action potential with using a large number of gate variables. For example, the Luo–Rudy model consists of 14 ODEs [6], the Courtemanche et al. model consists of 21 ODEs [8], the Winslow et al. model uses 32 ODEs [9], the Puglisi–Bersmode model contains 17 ODEs [10]. However, the use of the relatively simple model allows us to reproduce a variety of forms of action potentials in cardiomyocytes without unnecessary complications. For simulation of action potential in cardiomyocytes, modification of the generalized model of Hodgkin-Huxley was performed by adding a branch corresponding to the calcium ion channel [11]. Model of cell membrane at occurrence of action potential in cardiomyocyte is described by the parallel-conductance model (Fig. 1.6), - 24-.

(25) based on the assumption of the presence of independent ion channels for K+, Na+ and Ca2+ ions, as well as leakage through the membrane. Each branch of the electrical circuit reflects the contribution of one type of ions to total membrane current. The membrane of cardiomyocytes has selective permeability to ions of various types. The voltage sources in the equivalent circuit simulate Nernst potentials for each type of ion.. intracellular Im. gK. gNa. gCa. gl. Cm Vm EK. ENa. El. ECa. extracellular. Fig. 1.6. The parallel-conductance model of cardimyocyte’s cell membrane. Independent conductance channels are presented for K+, Na+, Ca2+ and leakage together with other factors. The battery polarity is chosen to show, that usually Nernst potentials of EK and El are negative (inside potential is more negative than outside), whereas ENa and ECa are positive (inside potential is more positive than outside). - 25-.

(26) Transmembrane potential at any moment of time can be represented as Vm (t )  Vm 0  v m (t ). (1.1). where v m (t ) is the alternating component of membrane potential, Vm 0 is resting potential determined from the parallel-conductance model. of membrane as Vm 0 .  g K 0 E K  g Na 0 E Na  g Ca 0 ECa  g l E l , g K 0  g Na 0  gCa 0  g l. (1.2). where g K 0 , g Na 0 and gCa 0 are the conductances of potassium, sodium and calcium ions at rest; E K , E Na , and ECa are electromotive forces of sources, which simulate Nernst potential of potassium, sodium and calcium respectively; g l is the leakage conductivity through the membrane; E l is electromotive force of source, which simulates Nernst potential for chlorine ions, leakage and other factors that affect the membrane potential in rest. Value of E l in practice is selected so, that for a given conductivity g l resting potential is equal to a predetermined value Vm 0 . As follows from the parallel-conductance model, the variable component of the membrane potential during the flow of the depolarizing current Id through the membrane must satisfy the differential equation dv m 1   IK (v m , t )  INa (v m , t )  ICa (v m , t )  Il  Id  , dt Cm. where IK (v m ,t )=g K (v m ,t )(Vm 0 +v m +EK ) ,. - 26-. (1.3).

(27) I Na (v m , t ) = g Na (v m , t )(Vm 0 + v m - E Na ) , ICa (v m , t ) = gCa (v m , t )(Vm 0 + v m - ECa ). are time dependences of potassium, sodium and calcium currents respectively; I l = g l (Vm 0 + v m  E l ) is the leakage current through the membrane;. I , 0  t  Td Id   d 0 t  Td  0,. is depolarizing current; Id0 is the amplitude of. depolarizing current; Td is duration of depolarizing pulse. During the membrane polarization conductance of potassium, sodium and calcium channels can be described by the equations:. g K (v m , t )  g K max n 4 (v m , t ) ,. (1.4). g Na (v m , t )  g Namax m3 (v m , t )h(v m , t ) ,. (1.5). gCa (v m , t )  g Ca max d (v m , t )f (v m , t ) ,. (1.6). where g K max , g Na max and gCamax are the membrane conductances for potassium, sodium and calcium ions respectively in the case all the channels of the membrane for this type of ions are in the open state; n is activation function of potassium channels; m is activation function of sodium channels; h is inactivation function for sodium channels; d is activation function of calcium channels and f is inactivation function for calcium channels.. - 27-.

(28) Activation and inactivation functions have the meaning of probability, that activation or inactivation ion channel subunit is in the open state at the moment. These gating variables change from zero to one and they are the solutions of the set of nonlinear ordinary differential equations (ODE). The proposed model contains ionic currents, which are determined by five gating variables n, m, h, d, and f. dn n  n  , dt n dh h  h  , dt h. dm m  m  , dt m. dd d   d  , dt d. (1.7). df f  f  , dt f. where. n . m . d . 1 ,  Vhn  v m   1  exp  VSn . 1 ,  Vhm  v m   1  exp V   Sm. 1 ,  Vhd  v m   1  exp V   Sd. h . 1 hmin.   v V    1 exp m hh    VSh   . f . Ph. 1fmin.   v V  1 exp m hf  VSf .   . Pf.  hmin ,.  fmin ,. n , m , and d  are the steady-state values of activation function for potassium, sodium and calcium channels respectively; h and f are the steady-state values of inactivation function for sodium and calcium channels; hmin and fmin are minimum values of inactivation for sodium and - 28-.

(29) calcium channels; Vhn , Vhm and Vhd are voltages of semi-activation for potassium, sodium and calcium channels; VSn , VSm and VSd are voltages of activation shape for potassium, sodium and calcium channels; Vhh and Vhf are voltages of semi-inactivation for sodium and calcium. channels; VSh and VSf are voltages of inactivation shape for sodium and calcium channels; Ph and Pf are degrees of relaxation of inactivation for sodium and calcium channels. Relaxation periods of activation n , m , and d for potassium, sodium and calcium channels were determined as follows.       1  rn n  n max   rn  , Pn  v m  Vhn     1  exp        Vsn           1  rm m  m max   r , m Pm     v  V   m hm    1  exp  V  s m        1  rd d  d max   v m  Vhd    1  exp  V s d  .   . Pd.     r d ,   . where  n max ,  m max , and  d max are maximum relaxation periods of activation for potassium, sodium and calcium channels; rn , rm , and rd are ratios of minimum relaxation period to maximum relaxation period of - 29-.

(30) activation for potassium, sodium and calcium channels; Vhn , Vhm , and Vhd are voltages of semi-relaxation of activation; VSn , VSm , and VSd are voltages of relaxation shape of activation; Pn , Pm , and P d are degrees of relaxation of activation for potassium, sodium and calcium channels, respectively. Relaxation periods of inactivation h and f for sodium and calcium channels were determined as follows       1 rh h  h max   rh  , Ph   1 exp  v m Vhh         VS h           1 rf f  f max   rf  , Pf   1 exp v m Vhf         VSf    . h max and f max are maximum relaxation periods of inactivation for sodium and calcium channels; rh and rf are ratios of minimum relaxation period to maximum relaxation period of inactivation for sodium and calcium channels; Vhh and Vhf. are voltages of semi-relaxation of. inactivation; VSh and VSf are voltages of relaxation shape of inactivation;. Ph and Pf are degrees of relaxation of inactivation for sodium and calcium channels, respectively. The initial conditions are. - 30-.

(31) v m (0)  0, n(0)  n0 , m(0)  m0 , d (0)  d0 ,. (1.8). h(0)  h0 , f (0)  f0 ,. where n0 , m0 , d0 , h0 and f0 can be found as n , m , d  , h and f for v m  0.. Maximal conductances of the membrane for K+, Na+ and Ca2+ channels are defined as g K max . gK 0 ; n04. g Na max . g Na 0 ; m03 h0. (1.10). g Ca max . g Na 0 . d 0 f0. (1.11). (1.9). Equation (1.3), the set of equations for all gating variables (1.7), as well as initial conditions (1.8) make Cauchy problem relatively functions v m , n, m, d , h and f , which is a set of stiff differential equations. Therefore,. to solve it, the implicit methods of integration should be used [12]. Using the predictor-corrector method, in order to avoid overflow of the computer grid, step t in the initial segment of integration should not be too large. The initial value should be selected as t  Td , where Td is duration of the depolarizing pulse. As a simplified, but a rough solution of the Cauchy problem, the following method can be used. Suppose that speed of transitions of - 31-.

(32) activation and inactivation subunits of ion channels from closed state into open and backwards changes slightly within the integration step. Then relaxation periods n , m , h , d , f and values n , m , d  , h , f can be considered as slightly changing. Accepting that these values are the constants depending only on voltage and not on time, differential equations can be solved explicitly. Then n(v m , t )  n  (n  n0 ) exp(  t. );. h(v m , t )  h  (h  h0 ) exp(  t. m );. n m(v m , t )  m  (m  m0 ) exp(  t. );. h t d (v m , t )  d   (d   d 0 ) exp(  ); d f (v m , t )  f  (f  f0 ) exp(  t ). f. Such approach makes it possible to analyze the parameters of the model and find their primary approximation in accordance with action potential shape. The. value. of. El. is. determined. from. (1.2).. Values. g K max , g Na max , and g Ca max are found from the equations (1.9-1.11) given. the fact, that n0 , m0 , d0 , h0 and f0 are the values of activation and inactivation functions n , m , d  , h and f for the steady state, when v m ( 0 )  0.. Together with transmembrane potential dependence on the time for action potential generation, changes of membrane currents for potassium ions, sodium, calcium and influence of leakage, chlorine ions, and other factors that affect the membrane potential are considered in the model using the Hodgkin-Huxley formalism. The model of electrical activity in cardiomyocyte simulates action potential dynamics and - 32-.

(33) dynamics of transmembrane currents.. Fig. 1.7. Activation function for potassium channel versus alternating component of membrane potential for the voltages of semi-activation Vhn = 50 mV and various voltages of activation shape Vsn To obtain the correct results with the considered model, it is necessary to identify the model parameters and the relationships between their values. For the occurrence of action potential, the condition I Na  I K. is necessary. In this condition the overall flow of. cations into the cell causes an increase in potential Vm and starts the process, which characterizes the growth phase of action potential. Growth of Vm causes growth of gNa , as channel conductivity is a function that depends continuously on the transmembrane potential, which in turn leads to a further increase of Vm , etc. The relaxation period of activation of sodium channels m is very small. The action potential can arise, when - 33-.

(34) m  n and h , in case of a state with a low value of n, normal value of m and increased value of h. This combination of three parameters leads. Fig.1.8. Gating variables, which determine ionic currents in the model: а) activation function (n) of potassium channels; b) activation (m) and inactivation (h) functions of sodium channels; c) activation (d) and inactivation (f) functions of calcium channels to the fulfillment of the condition INa  IK and excitation may occur in a cell. Resting membrane potential is mostly determined by high potassium permeability (conductance) and its value is close to the Nernst potential of potassium EK. Taking this into account, it can be seen from the - 34-.

(35) equation (1.2), that conductance of calcium ions at rest g Ca 0 should be small and should not practically affect the resting potential value. Meanwhile, in order to make calcium current to contribute to the stage of cardiomyocyte membrane repolarization, especially in the plateau phase of the action potential, the maximum conductance of the membrane for calcium channels g Ca max should be quite large and comparable to g K max . This can be achieved by a very small value of the parameter d 0 , which in its turn is provided by large value of Vhd and small value of Vsd . Maximum relaxation period of inactivation of calcium channels f max defines the duration of repolarization, which is long enough. in. cardiomyocytes.. specify dmaх  fmaх .. Therefore,. it. is. necessary. to. The equations, listed above, were solved for. simulation of action potential in ventricular cardiomyocytes using parameter values presented in TABLE I. Gating variables for potassium, sodium and calcium channels, obtained by solving the set of differential equations, are presented in Fig.1.8, a-c. Obtained simulated action potentials for atrial and ventricular cardiomyocytes (Fig. 1.9) in accordance with the theoretical information [1-4] have three characteristic phases: depolarization, plateau and repolarization. Depolarization phase is determined by a sharp increase of AP amplitude due to the growth of membrane permeability for sodium ions. Plateau phase describes the slow decline of action potential, which is explained by the simultaneous work of two types of channels - slow calcium channels and potassium channels. Repolarization phase is determined by the faster decline in AP amplitude, which is caused by closing of calcium channels and strengthening of the outward potassium current. - 35-.

(36) TABLE I.. PARAMETERS OF THE MODEL. Parameter. Value. Parameter Value. Parameter. Value. Cm. 1 F/cm2. Vhm. 40 mV. d max. 70 ms. Vm0. -85 mV. V hh. 7 mV.  f max. 260 ms. EK. 80 mV. Vhd. 34 mV. Vhn. 40 mV. ENa. 50 mV. Vhf. 5 mV. Vhm. 40 mV. ECa. 120 mV. Vsn. 50 mV. Vhh. 50 mV. gK0. 7.10-7 S/cm2. Vsm. 10 mV. Vhd. 20 mV. gNa0. 5.10-8 S/cm2. Vsh. -5 mV. Vhf. 75 mV. gCa0. 8.10-10 S/cm2. Vsd. 2.5 mV. Vsn. 20 mV. gl. 15.10-6 S/cm2. Vsf. -1 mV. Vsm. 10 mV. Id0. 7.10-5 A/cm2.  n max. 50 ms. Vsh. -50 mV. Td. 0.5 ms.  m max. 0.5 ms. Vsd. 20 mV. Vhn. 200 mV.  h max. 2 ms. Vsf. 10 mV. The model, based on the formalism of Hodgkin-Huxley, should be considered as a model describing the average values of membrane currents. It assumes, that membrane contains a large number of channels, which simultaneously operate in a cardiomyocyte, so the description of membrane currents results in averaging over the ensemble of channels.. - 36-.

(37) Fig. 1.9. Simulated action potentials for atrial and ventricular cardiomyocytes. Pacemaker signals arise in specialized cardiomyocytes, which able to generate periodical oscillations of their membrane potential [1, 13]. Pacemaker cells are located in the sinoatrial node (SAN). The pacemaker impulses propagate from the SAN through cardiac conduction system: the atrioventricular node (AVN), Bundle of His and Purkinje fibers network, and cause the contraction of the myocardium. The cycle generation of SAN action potential determines the heart rate. The spontaneous depolarization (Fig. 1.10) is unique property of SAN cells that necessary for SAN pacemaker activity. There is important physiological difference between nodal cells and ventricular cells due to the specific mechanisms in ion channels. The potential for cardiomyocyte action of the working myocardium (atria and ventricles) does not normally develop spontaneously. Generation of the action potential is caused by a wave of depolarization of neighboring cells of the - 37-.

(38) cardiomyocyte, and its membrane potential becomes less negative. Pacemaker heart cells (or rhythm drivers) do not need external stimuls to generate action potential, since they have their own automatism, that is, the ability to spontaneously depolarize and to develop the action potential. The properties of rhythm drivers are possessed by cells of the sinoatrial,. atrioventricular. nodes. and. Purkinje. fibers. (natural. pacemakers). The cells of the ventricular myocardium do not normally have automaticity, but they can acquire this ability in myocardial ischemia [13].. Fig. 1.10. Membrane potential Vm for ventricular cells (a) and SAN cells (b) [1].. The action potential of pacemaker cells (in contrast to atrial and ventricular myocardium cells) has three features: - 38-.

(39) 1. The negative potential in pacemaker cells reaches a value of -60 mV, which is much less than the resting potential of ventricular myocardial cells (-90 mV). The fast sodium channels of pacemaker cells remain inactive with such a potential. 2. The phase 4 of the action potential for pacemaker cells is not horizontal (Fig. 1.10, b), but has an upward character, which reflects a gradual, spontaneous depolarization due to pacemaker current, which is carried by Na+ ions. The ion channels, through which the pacemaker flow passes, differ from the fast sodium channels responsible for the depolarization phase of myocardial cells. Penetration of positively charged Na+ ions through the pacemaker channels helps the membrane potential during phase 4 to become less negative until it reaches the threshold potential. In this case, a slow inactivation of the pacemaker channels occurs. 3. The phase 0 of the action potential for pacemaker cells has a longer duration and a smaller amplitude than in myocardial cells. This can be explained by the fact that the fast sodium channels of pacemaker cells are inactive, and the excitation of the action potential is due to the flow of Ca2+ currents into the cell through slow calcium channels. The repolarization of pacemaker cells occurs in the same way as in ventricular myocardial cells. It is due to the inactivation of calcium channels, as well as the activation of potassium channels and the release of K+ ions from the cell. Different currents of the SAN ion channels contribute to the generation of action potential, but their functional roles and interactions are still not understood [13-15]. In the recent papers three hypotheses were proposed to explain the mechanism underling the cardiac automaticity. Authors [14] discribed the - 39-.

(40) “membrane clock” model of pacemaking, which includes the “funny” current. (If),. an. inward. Na+/K+ current,. activated. by. membrane. hyperpolarization at negative voltages. In the “calcium clock” model and “coupled-clock” model of pacemaking [15] the key mechanism in the diastolic depolarization was a spontaneous rhythmic phenomenon of Ca2+ release from the sarcoplasmic reticulum (SR). And authors [13] focused on functional role of voltage gated Ca2+ channels of two families: L-type and T-type . Ca2+ currents of L- and T-type were recorded in cells of SAN, AVN, and Purkinje Fibers.. The current of L-type (ICa,L ) was defined as a. “high”-threshold Ca2+ current, activated from about -30 mV and Ca2+ current of. T-type. the. (ICa,T) was defined as a “low” threshold. Ca2+ current, activated at −50 mV, moreover both currents participated at the latter half of the slow diastolic depolarization. 1.3.3. Simulation of action potential in nerve fibers. Simulation of action potential in the membrane of nerve fiber is based on parallel-conductance model (Fig. 1.11), in which independent conductance channels are presented for K+, Na+, and Cl− ions (channel for Cl− can be replaced on channel for leakage through the membrane). Total membrane current is the sum of currents for each branch of the electrical circuit. To solve the Cauchy problem a system of differential equations is constructed. This system contains four equations, that describe membrane potential vm, activation function for potassium n, activation function for sodium m, inactivation function for sodium h.. - 40-.

(41) dv m 1   IK (v m , t )  INa (v m , t )  Il  Id , dt Cm. (1.12). dn n  n ,  dt n. (1.13). dm m  m ,  dt m. (1.14). dh h  h ,  dt h. (1.15). where IK(vm, t) is dependent on time potassium current, INa(vm, t) is dependent on time sodium current, Il is leakage current through the membrane, n . 1 is the relaxation period of activation for potassium channels,  n  n. n . 10(0.01  v m ) is rate for switching of activation subunits for  0.01  v m  exp   1  0.01 . potassium channels from closed state to open,  v  n  125 exp   m   0.08 . is rate for switching of activation subunits for. potassium channels from open position to closed,. n   n n is the steady-state value of activation function for potassium m  m . 1 is the relaxation period of activation for sodium channels,  m  m 100(0.025  v m ) is rate for switching of activation subunits for  0.025  v m  exp   1 0.01  . sodium channels from closed position to open, - 41-.

(42) v   m  4  103 exp   m  is rate for switching of activation subunits for  0.018 . sodium channels from open position to closed,. m   m m is the steady-state value of activation function for sodium, h . 1 is the relaxation period of inactivation for sodium channels,  h  h. v    h  70 exp   m  is rate for switching of inactivation subunits for  0.002 . sodium channels from closed position to open, 103 is rate for switching of inactivation subunits for h   0.03  v m  exp   1  0.01 . sodium channels from open position to closed,. h   h h is the steady-state value of inactivation function for sodium. To simulate of action potential in the membrane of nerve a system of differential equations is solved with initial conditions: vm(0) = 0, n(0) = n0, m(0) = m0, h(0) = h0.. (1.16). Thus, equations (1.12) - (1.16) are Cauchy problem relatively functions vm(t), n, m and h. This Cauchy problem can be presented as: dv m  f1(v m , n, m, h, t ); dt dn  f2 (v m , n, m, h, t ); dt dm  f3 (v m , n, m, h, t ); dt dh  f4 (v m , n, m, h, t ), dt. - 42-.

(43) where the functions f1, f2, f3, f4 are determined by the equations (1.12) - (1.15) and the initial conditions are determined by the expressions (1.16). There is the period after appearance of the action potential, during which the membrane cannot be excited. This period is called refractory period (Fig. 1.11, 1.12). 120. 25 u. m. 100. I. d. 20. V, mV. 15. 60 40. 10. I d , 10 -6A. 80. 20 5 0 -20. 0. 20. 40. 60. 80. 100. 0 120. t , 10-3sec Fig. 1.11. Membrane potential under the excitation by two depolarizing pulses with a period less than the refractory period. The presence of this condition can be explained on the basis of the behavior of the inactivation function for sodium channels h. After appearance of action potential the function h reduces to a very low value, which prohibits occurrence of excitation. To return to a normal level of function h some time should pass. Also the activation function of. - 43-.

(44) potassium channels n have to decrease, since the implementation of condition: |INa|> |IK| is necessary for appearance of excitation. 120. 25 um. 100. Id. 20. V, mV. 15. 60 40. 10. I d , 10 -6A. 80. 20 5 0 -20. 0. 20. 40. 60. 80. 100. 0 120. t , 10-3sec Fig. 1.12. Membrane potential under the excitation by two depolarizing pulses with a period more than the refractory period. Experimental studies have shown, that excitation of the membrane may occur after prolonged hyperpolarization. Since the electrode in the external environment is an anode, this phenomenon is called anode break excitation of the cell. The presence of such effect can be explained as follows. Immediately before the removal of hyperpolarization the level of h goes down, while the levels of m and n rise up. However, after restoring of normal Vm, the activation function of sodium channels m quickly restores its value, since the relaxation period of activation for sodium channels  m is very small.. Consequently, as a result of.  m    n , h , there is a state, which is characterized by reduced value of n, normal value of m and increased value of h. This combination of - 44-.

(45) three parameters creates a condition |INa| > |IK|, resulting in a possible excitation occurrence (Fig. 1.13). 100. 0. -0.5 u. 50. V, mV. d. -1. -1.5. I d , 10 -6A. m. I. 0 -2. -50. 0. 50. 100. 150. -2.5 200. t , 10-3sec Fig. 1.13. Membrane potential under anode break excitation of the cell. 1.3.4. Modeling of the circulation of exciting impulse in the myocardium. The purpose of modeling of excitation wave circulation in the heart structures is the research at the level of myocardial cells of mechanisms that cause arrhythmia. Cardiac arrhythmias are basically produced by one of three mechanisms: enhanced automaticity, triggered activity, or reentry mechanism. Reentry mechanism occurs when a propagating impulse does not fade out after normal activation of the heart and persists to re-excite the heart after the end of refractory period. This electrophysiologic mechanism is responsible for the majority of - 45-.

(46) arrhythmias. For one-dimensional modeling of the excitation impulse circulation a closed loop for reentry mechanism simulation in the myocardium is obtained by connecting the ends of the fiber. The membranes of cardiomyocytes have selective permeability for different types of ions. Membrane model (MM) for the case of action potential emergence is presented by parallel-conductance model, based on the assumption of existence of independent channels for K+, Na+ and Ca2+ ions. In modeling of excitation impulse circulation by of reentry mechanism a modification of the model with leading core and a generalized model of Hodgkin-Huxley was performed by creating a onedimensional ring model in the form of a closed chain of N cardiomyocytes (Fig. 1.14) [16-21], N = 100-120.. reΔz. reΔz riΔz. riΔz riΔz. riΔz. reΔz. reΔz. Fig. 1.14. Ring model of reentry circle with inhomogeneity of electrophysiological parameters of myocardial cells - 46-.

(47) Membrane potential Vm (z, t) is a function of time t and spatial coordinate z taken along the direction of propagation of the action potential. In this case, Vm ( z, t )  Vm 0  v m ( z, t ), where Vm 0 is the resting membrane potential, vm (z, t) is the variable component of the membrane potential, which is a deviation from the membrane potential. Model of the leading core connects the second derivative on the membrane potential with current flowing through the membrane and the external excitation current:  2v m  (ri  re )i m  i p re , z 2. where im is membrane current per unit length, measured in А/m; re and ri are axial resistances of extracellular and intracellular environment per unit length, measured in Ohm/m; ip. is. external. stimulus,. measured. in. A/m.. Taking into account the fact that re . 4ρe ,  D2. ri . 4ρi ,  D2. i m  J m D ,. we can obtain that 4ρe  2v m 4  ρi  ρe   Jm  ip , 2 z D πD 2. where e and i are resistivities of extracellular and intracellular environment, measured in Ohm ∙ m; Jm is membrane current density, measured in А/m2; - 47-.

(48) D is diameter of the cylindrical cell. The spread of cardiac action potential along the fiber is represented as a differential equation: 4ρ  2v m 4  ρi  ρe  v m 4  ρi  ρe   C   JK  JNa  JCa  Jl   e2 i p . m 2 z D t D D. Pathway in myocardium [0, L] is divided by uniform grid with a step Δz . L , where L is length of reentry circle, N is number of mesh nodes. N. Cardiomyocytes. length. ranges. from. 20. microns. for. atria. to. 60-140 microns for ventricles and 150-200 microns for cells of Purkinje fibers.. Simulating. reentry. loop. with. length. in. range. of. 8 – 12 mm, number of mesh nodes zi must be chosen equal to the number of cells in the investigated area, for instance N = 100 – 120. The second spatial derivative of the membrane potential in the grid of nodes can be approximated by using an interpolation polynomial of Lagrange and calculated explicitly:  2v m ( zi , t ) v i -1(t ) - 2v i (t )  v i 1(t )  , i =2, N -1. z 2 z 2. Given the fact, that pathway is the reentry circle, the second spatial derivative of the membrane potential in the endpoints of the fiber is defined as:  2v m ( z1, t ) v N (t ) - 2v1(t )  v 2 (t )  ; z 2 z 2  2v m ( zN , t ) v1(t ) - 2v N (t )  v N -1(t ) .  z 2 z 2. - 48-.

(49) Determination of membrane potential vm and activation function of potassium channels n, activation function of sodium channels m, inactivation function for sodium channels h, activation function of calcium channels d and inactivation function for calcium channels f was carried out by solving the system of differential equations. Numerous parameters of. the. defined. functions. are. vectors. that. determine. the. electrophysiological properties of each individual cell in the reentry loop. The speed of action potential propagation through the cardiac tissue significantly varies in different areas of the heart. This speed directly depends on diameter of muscle fibers involved in the process. Conduction through cells of AV node, which have small diameter, is much slower than conduction through cells of Purkinje fibers with large diameter. The velocity of conduction also depends directly on the intensity of local depolarizing ion current, which in turn directly depends on intensity of action potential increase. Rapid depolarization facilitates rapid conduction. Differences in capacitance and (or) resistance of cell membranes are also the factors, which cause the difference in values of speed of action potential propagation at certain areas of the heart. The anatomical reentry mechanism is the simplest case for modeling the circular motion of the excitation pulse. When excitation pulse circulates in anatomically isolated enclosed chain, the path length is fixed and determined by perimeter of the anatomical structure, which forms the unexcited central part of the closed path. The reentry loop in this case is usually characterized by a large size of general perimeter of the circle (macro reentry). Modeling this case, we do not need to simulate a region with altered electrophysiological parameters. The same parameters can be used for all cells in the enclosed chain. The. functional. reentry. mechanism - 49-. is. determined. by. the.

(50) electrophysiological properties of the myocardial tissue. When these properties change, the length of the closed path changes too. The size of reentry circuit is not fixed, and its position is not strictly localized. This mechanism is usually characterized by a small loop and referred as micro reentry (6 - 10 mm in atrial fibrillation). Considering that differences in refractoriness period and other electrophysiological parameters of neighboring cell groups can cause circulation in the tissues of atria, ventricles and Purkinje fibers, the mathematical model was constructed in such a way, that the enclosed loop was inhomogeneous in its parameters. By varying the range of values of the variables that affect appearance and maintenance of reentry mechanism, and also changing the size of the loop and the ratio of the lengths of sections with different electrophysiological parameters, different patterns of propagation of excitation pulse along the investigated chain were obtained. The results of the simulation showed, that at different ratios of the parameters of cardiomyocytes in two sections of enclosed path, the existence of such variants is possible:  a single passage of the excitation pulse along the loop without the occurrence of circulation (Fig. 1.15);  occurrence of continuous circulation of the excitation pulse, for interruption of which the external intervention is required (Fig. 1.16);  several turns of the excitation pulse along enclosed chain with following attenuation (Fig. 1.17);  propagation of the excitation pulse only along the region with normal electrophysiological parameters and instantaneous attenuation in the affected region of the loop without the occurrence of circulation.. - 50-.

(51) Vm, mV. The area with the changed parameters. z, mm. t, ms. Vm, mV. Fig. 1.15. Single passage of the excitation pulse without circulation. z,mm t, ms. Fig. 1.16. Spatio-temporal representation of the continuous circulation of the excitation pulse along the loop with delayed propagation of action potentials - 51-.

(52) Ke1=4 mM/l. Vm, mV. Ke2=12 mM/l. z,mm. м t, мс. Fig. 1.17. Spatio-temporal propagation of action potential along the loop, which contains a region with an increased extracellular concentration of potassium ions Ke2 Action potential conduction, i.e. velocity of excitation propagation in cardiac structures, depends on the anatomical and physiological factors. The anatomical factors are diameter of the muscle fibers (10 – 15 μm for cardiomyocytes) and geometric arrangement of the muscle fibers. The velocity of conduction along the muscle fiber is greater than in the transverse direction. The physiological factors are amplitude of action potential, depolarization speed, as well as value of rest potential. Normally, the excitation pulse spreads along the tissues of the heart at speed of 0.5 - 5 m/s. One of the reasons leading to the development of continuous circulation of the excitation pulse is the slow propagation of action potentials. Fig. 1.16 shows the simulated results of the excitation pulse circulation by reentry mechanism, obtained by decrease in the conducting speed to 0.068 m/s, caused by a change in the - 52-.

(53) electrophysiological parameters of cardiomyocytes in the atrial affected region. The time of the excitation pulse turnover t t around the loop can be defined as: tt . where vc. is. L the. is speed. the of. the. L , vc. length excitation. of pulse. the. circle;. conduction.. Propagating at velocity of 0.068 m/s along a loop with L=10 mm, the excitation pulse returns to its initial point of motion through 147 ms, i.e. at the end of an effective refractory period of 145 ms. The cells, which recovered their excitability, reactivate and maintain the circulation of the excitation pulse according to reentry mechanism. The increase in the action potential on the affected areas occurs unevenly. In regions with a notable decrease in the action potential, the slow wave excitation can occur. Decrease in the negative value of the resting potential as a result of partial depolarization also slows down conduction velocity. At rest state, the intracellular concentration of potassium ions for myocardial cells is approximately 30 times higher than the extracellular concentration. Resting potential becomes less negative, when the extracellular concentration of potassium ions Ke increases. Changes in intracellular potassium concentration Ki are limited and do not significantly affect the value of Vm0. Fig. 1.17 depicts the simulation results for Ki = 150 mM/l throughout the entire length of the conductive path, Ke1 = 4 mM/l in the main region and Ke2 = 12 mM/l in the anomalous region. The excitation pulse in this case made 1 complete - 53-.

(54) turnover along the loop with L=12 mm, then again passed through the main loop section, equal in this case to 2/3 of the total circumference. The attenuation of the excitation pulse occurred at the boundary of the region with an increased extracellular concentration of potassium ions. CONCLUSIONS The changes of shape, duration, and pathways of action potentials propagation in the heart are manifestations of the mechanisms underlying the development of arrhythmias. In some cases, spontaneous attenuation of the excitation wave and termination of its circulation after several turnovers along a closed path are observed, in other cases a short strong single electric pulse is required to interrupt the circulation, which is used in clinical practice. The proposed models can be used to determine the boundary conditions, under which the circulation of the excitation pulse leads to a breakdown of heart rhythm. CONTROL QUESTIONS AND TASKS 1. Explain. the. phenomenon. of. selective. permeability. of. biomembranes. 2. What concentration is bigger for potassium ions in muscle and nerve cells: intracellular concentration Ki or extracellular concentration Ke? 3. What concentration is bigger for sodium ions in muscle and nerve cells: intracellular concentration Nai or extracellular concentration Nae? - 54-.

(55) 4. Explain the Nernst-Planck equation. 5. Explain how to identify the Nernst potential. 6. Explain the Goldman equation. 7. Explain the nature of the resting potential. Which factors does it depend on? 8. How does change the membrane resting potential, in case if - intracellular concentration of sodium ions Nai increases (decreases); - extracellular concentration of sodium ions Nae increases (decreases); - intracellular concentration of potassium ions Ki increases (decreases); - extracellular concentration of potassium ions Ke increases (decreases); - permeability of K+ ions increases (decreases); - permeability of Na+ ions increases (decreases)? 9. Draw the parallel-conductance model. 10. Explain the phenomena of action potential genesis. 11. What condition should be fulfilled between sodium and potassium currents (INa and IK ) for the emergence of membrane action potential? 12. How does the amplitude of action potential of the cell membrane change, if the intracellular concentration of sodium ions Nai increases? 13. What is the refractory period and what are the mechanisms of this phenomenon? 14. What is the mechanism of anode break excitation of the cell? 15. Explain the system of equations for solving the Cauchy problem for the process of propagation of the action potential along cellular fiber.. - 55-.

(56) REFERENCES [1] R. Plonsey and R. Barr, Bioelectricity. A Quantitative Approach, third ed. New York: Springer Science+Business Media, LLC, 2007. [2] A.Hodgkin and A. Huxley, "A quantitative description of membrane current and its application to conduction and excitation in nerve" J. Physiol., 117, pp. 500–544, 1952. [3] R.R. Aliev, "Computer modelling of heart electrical activity," Advances of Physiological Sciences, vol. 41. no. 3, pp. 44-63, 2010. [4] http://www.all-fizika.com/article/index.php?id_article=1977 [5] McAllister R.E., Noble D., Tsien R.W. "Reconstruction of the electrical activitity of cardiac purkinje ﬁbres," J. Physiol., vol. 251, pp. 1–59, 1975. [6] C. Luo and Y. Rudy, "A model of ventricular cardiac action potential," Circ. Res., 68, 1501–1526, 1991. [7] C. Luo and Y. Rudy, “A dynamic model of the cardiac ventricular action potential. i. Simulations of ionic currents and concentration changes,” Circ. Res., vol. 74, no. 6, pp. 1071–1096, 1994. [8] M. Courtemanche, R. J. Ramirez and S. Nattel, “Ionic mechanisms underlying human atrial action potential properties: Insights from a mathematical model,” Amer. J. Physiol., vol. 275, no. 1, pp. H301– H321, 1998. [9] R. L. Winslow et al., “Mechanisms of altered excitation-contraction coupling in canine tachycardia induced heart failure, ii model studies,” Circ.Res., vol.84, no.5, pp.571– 586, 1999. [10] J. L. Puglisi and D. M. Bers, “Labheart: An interactive computer model of rabbit ventricular myocyte ion channels and ca transport,” - 56-.

No results found