Forward problem of electrocardiology

Last updated
Simulation of a realistic heart beat. CG Heart.gif
Simulation of a realistic heart beat.
Schematic diagram of normal sinus rhythm for a human heart as seen on ECG. SinusRhythmLabels.svg
Schematic diagram of normal sinus rhythm for a human heart as seen on ECG.

The forward problem of electrocardiology is a computational and mathematical approach to study the electrical activity of the heart through the body surface. [1] The principal aim of this study is to computationally reproduce an electrocardiogram (ECG), which has important clinical relevance to define cardiac pathologies such as ischemia and infarction, or to test pharmaceutical intervention. Given their important functionalities and the relative small invasiveness, the electrocardiography techniques are used quite often as clinical diagnostic tests. Thus, it is natural to proceed to computationally reproduce an ECG, which means to mathematically model the cardiac behaviour inside the body. [1]

Contents

The three main parts of a forward model for the ECG are:

Thus, to obtain an ECG, a mathematical electrical cardiac model must be considered, coupled with a diffusive model in a passive conductor that describes the electrical propagation inside the torso. [1]

The coupled model is usually a three-dimensional model expressed in terms of partial differential equations. Such model is typically solved by means of finite element method for the solution's space evolution and semi-implicit numerical schemes involving finite differences for the solution's time evolution. However, the computational costs of such techniques, especially with three dimensional simulations, are quite high. Thus, simplified models are often considered, solving for example the heart electrical activity independently from the problem on the torso. To provide realistic results, three dimensional anatomically realistic models of the heart and the torso must be used. [1]

Another possible simplification is a dynamical model made of three ordinary differential equations. [3]

Heart tissue models

The electrical activity of the heart is caused by the flow of ions across the cell membrane, between the intracellular and extracellular spaces, which determines a wave of excitation along the heart muscle that coordinates the cardiac contraction and, thus, the pumping action of the heart that enables it to push blood through the circulatory system. The modeling of cardiac electrical activity is thus related to the modelling of the flow of ions on a microscopic level, and on the propagation of the excitation wave along the muscle fibers on a macroscopic level. [1] [4]

Between the mathematical model on the macroscopic level, Willem Einthoven and Augustus Waller defined the ECG through the conceptual model of a dipole rotating around a fixed point, whose projection on the lead axis determined the lead recordings. Then, a two-dimensional reconstruction of the heart activity in the frontal plane was possible using the Einthoven's limbs leads I, II and III as theoretical basis. [5] Later on, the rotating cardiac dipole was considered inadequate and was substituted by multipolar sources moving inside a bounded torso domain. The main shortcoming of the methods used to quantify these sources is their lack of details, which are however very relevant to realistically simulate cardiac phenomena. [4]

On the other hand, microscopic models try to represent the behaviour of single cells and to connect them considering their electrical properties. [6] [7] [8] These models present some challenges related to the different scales that need to be captured, in particular considering that, especially for large scale phenomena such as re-entry or body surface potential, the collective behaviour of the cells is more important than that of every single cell. [4]

The third option to model the electrical activity of the heart is to consider a so-called "middle-out approach", where the model incorporates both lower and higher level of details. This option considers the behaviour of a block of cells, called a continuum cell, thus avoiding scale and detail problems. The model obtained is called bidomain model, which is often replaced by its simplification, the monodomain model. [4]

Bidomain model

Stylized representation of a human torso that describes the domain and the notation considered for the forward problem of electrocardiography. Two separated areas and their boundary are considered, which represent the heart and human torso around it. Torso domain for the forward problem of electrocardiography.png
Stylized representation of a human torso that describes the domain and the notation considered for the forward problem of electrocardiography. Two separated areas and their boundary are considered, which represent the heart and human torso around it.

The basic assumption of the bidomain model is that the heart tissue can be divided in two ohmic conducting continuous media, connected but separated through the cell membrane. This media are called intracellular and extracellular regions, the former representing the cellular tissues, and the latter representing the space between cells. [2] [1]

The standard formulation of the bidomain model, including a dynamical model for the ionic current, is the following [2]

where and are the transmembrane and extracellular potentials respectively, is the ionic current, which depends also from a so-called gating variable (accounting for cellular-level ionic behavior), and is an external current applied to the domain. Moreover, and are the intracellular and extracellular conductivity tensors, is the surface to volume ratio of the cell membrane and is the membrane capacitance per unit area. Here the domain represents the heart muscle. [2]

The boundary conditions for this version of the bidomain model are obtained through the assumption that there is no flow of intracellular potential outside of the heart, which means that

where denotes the boundary of the heart domain and is the outward unit normal to . [2]

Monodomain model

The monodomain model is a simplification of the bidomain model that, in spite of some unphysiological assumptions, is able to represent realistic electrophysiological phenomena at least for what concerns the transmembrane potential . [2] [1]

The standard formulation is the following partial differential equation, whose only unknown is the transmembrane potential:

where is a parameter that relates the intracellular and extracellular conductivity tensors. [2]

The boundary condition used for this model is [9]

Torso tissue model

In the forward problem of electrocardiography, the torso is seen as a passive conductor and its model can be derived starting from the Maxwell's equations under quasi-static assumption. [1] [2]

The standard formulation consists of a partial differential equation with one unknown scalar field, the torso potential . Basically, the torso model is the following generalized Laplace equation

where is the conductivity tensor and is the domain surrounding the heart, i.e. the human torso. [2]

Derivation

As for the bidomain model, the torso model can be derived from the Maxwell's equations and the continuity equation after some assumptions. First of all, since the electrical and magnetic activity inside the body is generated at low level, a quasi-static assumption can be considered. Thus, the body can be viewed as a passive conductor, which means that its capacitive, inductive and propagative effect can be ignored. [1]

Under quasi-static assumption, the Maxwell's equations are [1]

and the continuity equation is [1]

Since its curl is zero, the electrical field can be represented by the gradient of a scalar potential field, the torso potential

 

 

 

 

(1)

where the negative sign means that the current flows from higher to lower potential regions. [1]

Then, the total current density can be expressed in terms of the conduction current and other different applied currents so that, from continuity equation, [1]

 

 

 

 

(2)

Then, substituting ( 1 ) in ( 2 )

in which is the current per unit volume. [1]

Finally, since aside from the heart there is no current source inside the torso, the current per unit volume can be set to zero, giving the generalized Laplace equation, which represents the standard formulation of the diffusive problem inside the torso [1]

Boundary condition

The boundary conditions accounts for the properties of the media surrounding the torso, i.e. of the air around the body. Generally, air has null conductivity which means that the current cannot flow outside the torso. This is translated in the following equation [1]

where is the unit outward normal to the torso and is the torso boundary, which means the torso surface. [1] [2]

Torso conductivity

Usually, the torso is considered to have isotropic conductivity, which means that the current flows in the same way in all directions. However, the torso is not an empty or homogeneous envelope, but contains different organs characterized by different conductivity coefficients, which can be experimentally obtained. A simple example of conductivity parameters in a torso that considers the bones and the lungs is reported in the following table. [2]

Values of the torso conductivity. [2]
(S/cm)(S/cm)(S/cm)

Heart-torso models

The coupling between the electrical activity model and the torso model is achieved by means of suitable boundary conditions at the epicardium, i.e. at the interface surface between the heart and the torso. [1] [2]

The heart-torso model can be fully coupled, if a perfect electrical transmission between the two domains is considered, or can be uncoupled, if the heart electrical model and the torso model are solved separately with a limited or imperfect exchange of information between them. [2]

Fully coupled heart-torso models

The complete coupling between the heart and the torso is obtained imposing a perfect electrical transmission condition between the heart and the torso. This is done considering the following two equations, that establish a relationship between the extracellular potential and the torso potential [2]

This equations ensure the continuity of both the potential and the current across the epicardium. [2]

Using these boundary conditions, it is possible to obtain two different fully coupled heart-torso models, considering either the bidomain or the monodomain model for the heart electrical activity. From the numerical viewpoint, the two models are computationally very expensive and have similar computational costs. [2]

Alternative boundary conditions

Boundary conditions that represent a perfect electrical coupling between the heart and the torso are the most used and the classical ones. However, between the heart and the torso there is the pericardium, a sac with a double wall that contains a serous fluid which has a specific effect on the electrical transmission. Considering the capacitance and the resistive effect that the pericardium has, alternative boundary conditions that take into account this effect can be formulated as follows [10]

Formulation with the bidomain model

The fully coupled heart-torso model, considering the bidomain model for the heart electrical activity, in its complete form is [2]

where the first four equations are the partial differential equations representing the bidomain model, the ionic model and the torso model, while the remaining ones represent the boundary conditions for the bidomain and torso models and the coupling conditions between them. [2]

Formulation with the monodomain model

The fully coupled heart-torso model considering the monodomain model for the electrical activity of the heart is more complicated that the bidomain problem. Indeed, the coupling conditions relate the torso potential with the extracellular potential, which is not computed by the monodomain model. Thus, it is necessary to use also the second equation of the bidomain model (under the same assumptions under which the monodomain model is derived), yielding: [2]

This way, the coupling conditions do not need to be changed, and the complete heart-torso model is composed of two different blocks: [2]

  • First the monodomain model with its usual boundary condition must be solved:
  • Then, the coupled model that includes the computation of the extracellular potential, the torso model and the coupling conditions must be solved:

Uncoupled heart-torso models

The fully coupled heart-torso models are very detailed models but they are also computationally expensive to solve. [2] A possible simplification is provided by the so-called uncoupled assumption in which the heart is considered completely electrically isolated from the heart. [2] Mathematically, this is done imposing that the current cannot flow across the epicardium, from the heart to the torso, namely [2]

Applying this equation to the boundary conditions of the fully coupled models, it is possible to obtained two uncoupled heart-torso models, in which the electrical models can be solved separately from the torso model reducing the computational costs. [2]

Uncoupled heart-torso model with the bidomain model

The uncoupled version of the fully coupled heart-torso model that uses the bidomain to represent the electrical activity of the heart is composed of two separated parts: [2]

  • The bidomain model in its isolated form
  • The torso diffusive model in its standard formulation, with the potential continuity condition

Uncoupled heart-torso model with the modomain model

As in the case of the fully coupled heart-torso model which uses the monodomain model, also in the corresponding uncoupled model extracellular potential needs to be computed. In this case, three different and independent problems must be solved: [2]

  • The monodomain model with its usual boundary condition:
  • The problem to compute the extracellular potential with a boundary condition on the epicardium prescribing no intracellular current flow:
  • The torso diffusive model with the potential continuity boundary condition at the epicardium:

Electrocardiogram computation

Precordial leads in ECG Precordial leads in ECG.png
Precordial leads in ECG

Solving the fully coupled or the uncoupled heart-torso models allows to obtain the electrical potential generated by the heart in every point of the human torso, and in particular on the whole surface of the torso. Defining the electrodes positions on the torso, it is possible to find the time evolution of the potential on such points. Then, the electrocardiograms can be computed, for example according to the 12 standard leads, considering the following formulas [2]

where and are the standard locations of the electrodes. [2]

Numerical methods

The heart-torso models are expressed in terms of partial differential equations whose unknowns are function of both space and time. They are in turn coupled with an ionic model which is usually expressed in terms of a system of ordinary differential equations. A variety of numerical schemes can be used for the solution of those problems. Usually, the finite element method is applied for the space discretization and semi-implicit finite-difference schemes are used for the time discretization. [1] [2]

Uncoupled heart-torso model are the simplest to treat numerically because the heart electrical model can be solved separately from the torso one, so that classic numerical methods to solve each of them can be applied. This means that the bidomain and monodomain models can be solved for example with a backward differentiation formula for the time discretization, while the problems to compute the extracellular potential and torso potential can be easily solved by applying only the finite element method because they are time independent. [1] [2]

The fully coupled heart-torso models, instead, are more complex and need more sophisticated numerical models. For example, the fully heart-torso model that uses the bidomain model for the electrical simulation of the cardiac behaviour can be solved considering domain decomposition techniques, such as a Dirichlet-Neumann domain decomposition. [2] [11]

Geometric torso model

Three-dimensional torso model including the most organs. Journal.pone.0141573.g004.png
Three-dimensional torso model including the most organs.

To simulate and electrocardiogram using the fully coupled or uncoupled models, a three-dimensional reconstruction of the human torso is needed. Today, diagnostic imaging techniques such as MRI and CT can provide a sufficiently accurate images that allow to reconstruct in detail anatomical human parts and, thus, obtain a suitable torso geometry. For example, the Visible Human Data [13] is a useful dataset to create a three-dimensional torso model detailed with internal organs including the skeletal structure and muscles. [1]

Dynamical model for the electrocardiogram

Even if the results are quite detailed, solving a three-dimensional model is usually quite expensive. A possible simplification is a dynamical model based on three coupled ordinary differential equations. [3]

The quasi-periodicity of the heart beat is reproduced by a three-dimensional trajectory around an attracting limit cycle in the plane. The principal peaks of the ECG, which are the P,Q,R,S and T, are described at fixed angles , which give the following three ODEs [3]

with , ,

The equations can be easily solved with classical numerical algorithms like Runge-Kutta methods for ODEs. [3]

See also

Related Research Articles

Continuum mechanics is a branch of mechanics that deals with the deformation of and transmission of forces through materials modeled as a continuous medium rather than as discrete particles. The French mathematician Augustin-Louis Cauchy was the first to formulate such models in the 19th century.

In particle physics, the Dirac equation is a relativistic wave equation derived by British physicist Paul Dirac in 1928. In its free form, or including electromagnetic interactions, it describes all spin-12 massive particles, called "Dirac particles", such as electrons and quarks for which parity is a symmetry. It is consistent with both the principles of quantum mechanics and the theory of special relativity, and was the first theory to account fully for special relativity in the context of quantum mechanics. It was validated by accounting for the fine structure of the hydrogen spectrum in a completely rigorous way.

<span class="mw-page-title-main">Navier–Stokes equations</span> Equations describing the motion of viscous fluid substances

The Navier–Stokes equations are partial differential equations which describe the motion of viscous fluid substances. They were named after French engineer and physicist Claude-Louis Navier and the Irish physicist and mathematician George Gabriel Stokes. They were developed over several decades of progressively building the theories, from 1822 (Navier) to 1842–1850 (Stokes).

The vorticity equation of fluid dynamics describes the evolution of the vorticity ω of a particle of a fluid as it moves with its flow; that is, the local rotation of the fluid. The governing equation is:

In fluid dynamics, Stokes' law is an empirical law for the frictional force – also called drag force – exerted on spherical objects with very small Reynolds numbers in a viscous fluid. It was derived by George Gabriel Stokes in 1851 by solving the Stokes flow limit for small Reynolds numbers of the Navier–Stokes equations.

In mechanics, virtual work arises in the application of the principle of least action to the study of forces and movement of a mechanical system. The work of a force acting on a particle as it moves along a displacement is different for different displacements. Among all the possible displacements that a particle may follow, called virtual displacements, one will minimize the action. This displacement is therefore the displacement followed by the particle according to the principle of least action.

The work of a force on a particle along a virtual displacement is known as the virtual work.

In differential geometry, the four-gradient is the four-vector analogue of the gradient from vector calculus.

<span class="mw-page-title-main">Cartesian tensor</span>

In geometry and linear algebra, a Cartesian tensor uses an orthonormal basis to represent a tensor in a Euclidean space in the form of components. Converting a tensor's components from one such basis to another is done through an orthogonal transformation.

<span class="mw-page-title-main">Stokes flow</span> Type of fluid flow

Stokes flow, also named creeping flow or creeping motion, is a type of fluid flow where advective inertial forces are small compared with viscous forces. The Reynolds number is low, i.e. . This is a typical situation in flows where the fluid velocities are very slow, the viscosities are very large, or the length-scales of the flow are very small. Creeping flow was first studied to understand lubrication. In nature, this type of flow occurs in the swimming of microorganisms and sperm. In technology, it occurs in paint, MEMS devices, and in the flow of viscous polymers generally.

The intent of this article is to highlight the important points of the derivation of the Navier–Stokes equations as well as its application and formulation for different families of fluids.

The Cauchy momentum equation is a vector partial differential equation put forth by Cauchy that describes the non-relativistic momentum transport in any continuum.

The Clausius–Duhem inequality is a way of expressing the second law of thermodynamics that is used in continuum mechanics. This inequality is particularly useful in determining whether the constitutive relation of a material is thermodynamically allowable.

The bidomain model is a mathematical model to define the electrical activity of the heart. It consists in a continuum (volume-average) approach in which the cardiac microstructure is defined in terms of muscle fibers grouped in sheets, creating a complex three-dimensional structure with anisotropical properties. Then, to define the electrical activity, two interpenetrating domains are considered, which are the intracellular and extracellular domains, representing respectively the space inside the cells and the region between them.

<span class="mw-page-title-main">Stokes' theorem</span> Theorem in vector calculus

Stokes' theorem, also known as the Kelvin–Stokes theorem after Lord Kelvin and George Stokes, the fundamental theorem for curls or simply the curl theorem, is a theorem in vector calculus on . Given a vector field, the theorem relates the integral of the curl of the vector field over some surface, to the line integral of the vector field around the boundary of the surface. The classical theorem of Stokes can be stated in one sentence: The line integral of a vector field over a loop is equal to its curl through the enclosed surface. It is illustrated in the figure, where the direction of positive circulation of the bounding contour ∂Σ, and the direction n of positive flux through the surface Σ, are related by a right-hand-rule. For the right hand the fingers circulate along ∂Σand the thumb is directed along n.

In continuum mechanics, a compatible deformation tensor field in a body is that unique tensor field that is obtained when the body is subjected to a continuous, single-valued, displacement field. Compatibility is the study of the conditions under which such a displacement field can be guaranteed. Compatibility conditions are particular cases of integrability conditions and were first derived for linear elasticity by Barré de Saint-Venant in 1864 and proved rigorously by Beltrami in 1886.

The monodomain model is a reduction of the bidomain model of the electrical propagation in myocardial tissue. The reduction comes from assuming that the intra- and extracellular domains have equal anisotropy ratios. Although not as physiologically accurate as the bidomain model, it is still adequate in some cases, and has reduced complexity.

<span class="mw-page-title-main">Rock mass plasticity</span>

Plasticity theory for rocks is concerned with the response of rocks to loads beyond the elastic limit. Historically, conventional wisdom has it that rock is brittle and fails by fracture while plasticity is identified with ductile materials. In field scale rock masses, structural discontinuities exist in the rock indicating that failure has taken place. Since the rock has not fallen apart, contrary to expectation of brittle behavior, clearly elasticity theory is not the last word.

<span class="mw-page-title-main">Objective stress rate</span>

In continuum mechanics, objective stress rates are time derivatives of stress that do not depend on the frame of reference. Many constitutive equations are designed in the form of a relation between a stress-rate and a strain-rate. The mechanical response of a material should not depend on the frame of reference. In other words, material constitutive equations should be frame-indifferent (objective). If the stress and strain measures are material quantities then objectivity is automatically satisfied. However, if the quantities are spatial, then the objectivity of the stress-rate is not guaranteed even if the strain-rate is objective.

The variational multiscale method (VMS) is a technique used for deriving models and numerical methods for multiscale phenomena. The VMS framework has been mainly applied to design stabilized finite element methods in which stability of the standard Galerkin method is not ensured both in terms of singular perturbation and of compatibility conditions with the finite element spaces.

The streamline upwind Petrov–Galerkin pressure-stabilizing Petrov–Galerkin formulation for incompressible Navier–Stokes equations can be used for finite element computations of high Reynolds number incompressible flow using equal order of finite element space by introducing additional stabilization terms in the Navier–Stokes Galerkin formulation.

References

  1. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 Pullan, Andrew J.; Buist, Martin L.; Cheng, Leo K. (2005). Mathematically modelling the electrical activity of the heart : from cell to body surface and back again. World Scientific. ISBN   978-9812563736.
  2. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 Boulakia, Muriel; Cazeau, Serge; Fernández, Miguel A.; Gerbeau, Jean-Frédéric; Zemzemi, Nejib (24 December 2009). "Mathematical Modeling of Electrocardiograms: A Numerical Study" (PDF). Annals of Biomedical Engineering. 38 (3): 1071–1097. doi:10.1007/s10439-009-9873-0. PMID   20033779. S2CID   10114284.
  3. 1 2 3 4 McSharry, P.E.; Clifford, G.D.; Tarassenko, L.; Smith, L.A. (March 2003). "A dynamical model for generating synthetic electrocardiogram signals". IEEE Transactions on Biomedical Engineering. 50 (3): 289–294. doi:10.1109/TBME.2003.808805. PMID   12669985. S2CID   544816.
  4. 1 2 3 4 Lines, G.T.; Buist, M.L.; Grottum, P.; Pullan, A.J.; Sundnes, J.; Tveito, A. (1 July 2002). "Mathematical models and numerical methods for the forward problem in cardiac electrophysiology". Computing and Visualization in Science. 5 (4): 215–239. doi:10.1007/s00791-003-0101-4. S2CID   123211416.
  5. Einthoven, Waller (1903). "Die galvanometrische Registrierung des menschlichen Elektrokardiogramms, zugleich eine Beurteilung der Anwendung des Capillarelektrometers in der Physiologie". Pflügers Archiv. 99 (9–10): 472–480. doi:10.1007/BF01811855. S2CID   10400303.
  6. Henriquez, C.S.; Plonsey, R. (1987). "Effect of resistive discontinuities on waveshape and velocity in a single cardiac fibre". Med. Biol. Eng. Comput. 25 (4): 428–438. doi:10.1007/BF02443364. PMID   3450994. S2CID   3038844.
  7. Muller-Borer, B.J.; Erdman, D.J.; Buchanan, J.W. (1994). "Electrical Coupling and Impulse Propagation in Anatomically Modeled Ventricular Tissue". IEEE Trans. Biomed. Eng. 41 (5): 445–454. doi:10.1109/10.293219. PMID   8070804. S2CID   14407776.
  8. Hren, R.; Nenonen, J.; Horacek, B.M. (1998). "Simulated Epicardial Potential Maps During Paced Activation Reflect Myocardial Fibrous Structure". Annals of Biomedical Engineering. 26 (6): 1022–1035. doi:10.1114/1.73. PMID   9846940. S2CID   9978399.
  9. Keener, James; Sneyd, James (2009). Mathematical physiology 2009 : systems physiology ii (2nd revisedition ed.). Springer. ISBN   978-1-4939-3709-7.
  10. Boulakia, Muriel; Fernández, Miguel A.; Gerbeau, Jean-Frédéric; Zemzemi, Nejib (2007). "Towards the Numerical Simulation of Electrocardiograms". Functional Imaging and Modeling of the Heart. Lecture Notes in Computer Science. Vol. 4466. Springer. pp. 240–249. doi:10.1007/978-3-540-72907-5_25. ISBN   978-3-540-72906-8.
  11. Fernández, Miguel A.; Zemzemi, Nejib (1 July 2010). "Decoupled time-marching schemes in computational cardiac electrophysiology and ECG numerical simulation" (PDF). Mathematical Biosciences. 226 (1): 58–75. doi:10.1016/j.mbs.2010.04.003. ISSN   0025-5564. PMID   20416327. S2CID   8792966.
  12. Ferrer, Ana; Sebastián, Rafael; Sánchez-Quintana, Damián; Rodríguez, José F.; Godoy, Eduardo J.; Martínez, Laura; Saiz, Javier (2 November 2015). "Detailed Anatomical and Electrophysiological Models of Human Atria and Torso for the Simulation of Atrial Activation". PLOS ONE. 10 (11): e0141573. doi: 10.1371/journal.pone.0141573 . ISSN   1932-6203. PMC   4629897 . PMID   26523732.
  13. Spitzer, V.; Ackerman, M. J.; Scherzinger, A. L.; Whitlock, D. (1 March 1996). "The Visible Human Male: A Technical Report". Journal of the American Medical Informatics Association. 3 (2): 118–130. doi: 10.1136/jamia.1996.96236280 . PMC   116294 . PMID   8653448.