Forward Problem Of Electrocardiology

The forward problem of electrocardiology is a computational and mathematical approach to study the electrical activity of the heart through the body surface. 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.

The three main parts of a forward model for the ECG are:
* a model for the cardiac electrical activity;
* a model for the diffusion of the electrical potential inside the torso, which represents the extracardiac region;
* some specific heart-torso coupling conditions.
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.

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.

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

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.

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.
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.

On the other hand, microscopic models try to represent the behaviour of single cells and to connect them considering their electrical properties. 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.

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.

Bidomain model

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.

The standard formulation of the bidomain model, including a dynamical model for the ionic current, is the following
$\begin \nabla \cdot (\boldsymbol \sigma_i \nabla V_m ) + \nabla \cdot (\boldsymbol \sigma_i \nabla u_e ) = A_m \left( C_m \dfrac + I_\text(V_m,w)\right) + I_\text &\text \Omega_H\\ \nabla \cdot \left( (\boldsymbol \sigma_i + \boldsymbol \sigma_e) \nabla u_e \right) +\nabla \cdot(\boldsymbol \sigma_i \nabla V_m) = 0 & \text \Omega_H\\ \frac + g(V_m,w) = 0 &\text \Omega_H \end$
where $V_m$ and $u_e$ are the transmembrane and extracellular potentials respectively, $I_\text$ is the ionic current, which depends also from a so-called gating variable $w$ (accounting for cellular-level ionic behavior), and $I_\text$ is an external current applied to the domain. Moreover, $\boldsymbol \sigma_i$ and $\boldsymbol \sigma_e$ are the intracellular and extracellular conductivity tensors, $A_m$ is the surface to volume ratio of the cell membrane and $C_m$ is the membrane capacitance per unit area. Here the domain $\Omega_H$ represents the heart muscle.

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
$\boldsymbol \sigma_i \nabla V_m \cdot \mathbf n + \boldsymbol \sigma_i \nabla u_e \cdot \mathbf n = 0 \quad \quad \text\Sigma$
where $\Sigma = \partial\Omega_H$ denotes the boundary of the heart domain and $\mathbf n$ is the outward unit normal to $\Sigma$.

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 $V_m$.

The standard formulation is the following partial differential equation, whose only unknown $V_m$ is the transmembrane potential:
$\chi C_m \frac - \nabla \cdot \left(\boldsymbol \sigma_i \frac\nabla V_m\right) + \chi I_\text = I_\text \quad \quad \text\Omega_H$
where $\lambda$ is a parameter that relates the intracellular and extracellular conductivity tensors.

The boundary condition used for this model is
$\left(\boldsymbol \sigma_i \nabla V_m \right) \cdot \mathbf n = 0 \quad \quad \text\partial \Omega_H.$

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.

The standard formulation consists of a partial differential equation with one unknown scalar field, the torso potential $u_T$. Basically, the torso model is the following generalized Laplace equation
$\nabla \cdot (\boldsymbol \sigma_T \nabla u_T ) = 0 \quad \quad \text \Omega_T,$
where $\boldsymbol \sigma_T$ is the conductivity tensor and $\Omega_T$ is the domain surrounding the heart, ''i.e.'' the human torso.

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.

Under quasi-static assumption, the Maxwell's equations are
$\begin \nabla \cdot \mathbf E &= \frac \\ \nabla \times \mathbf E &= 0 \\ \nabla \cdot \mathbf B &= 0\\ \nabla \times \mathbf B &= \mu_0\mathbf J \end$
and the continuity equation is
$\nabla \cdot \mathbf J = 0 .$

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

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

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

Then, substituting () in ()
$\nabla \cdot (\boldsymbol \sigma_T \nabla u_T ) = \nabla \cdot \mathbf J_\text = I_v$

in which $I_v$ is the current per unit volume.

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
$\nabla \cdot (\boldsymbol \sigma_T \nabla u_T ) = 0 \quad \quad \text \Omega_T.$

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
$\boldsymbol \sigma_T \nabla u_T \cdot \mathbf n_T = 0 \quad \quad \text\Gamma_T$

where $\mathbf n_T$ is the unit outward normal to the torso and $\Gamma_T$ is the torso boundary, which means the torso surface.

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.

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.

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.

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
$\begin u_e &= u_T &\text\Sigma\\ (\boldsymbol \sigma_e \nabla u_e)\cdot \mathbf n_e &= - (\boldsymbol \sigma_T \nabla u_T) \cdot \mathbf n_T \quad &\text\Sigma. \end$
This equations ensure the continuity of both the potential and the current across the epicardium.

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.

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 $C_p$ and the resistive$R_p$ effect that the pericardium has, alternative boundary conditions that take into account this effect can be formulated as follows
$\begin R_p \boldsymbol \sigma_e \nabla u_e \cdot \mathbf n &= R_p C_p \frac + (u_T - u_e) \quad \quad &\text\Sigma \\ \boldsymbol \sigma_e \nabla u_e \cdot \mathbf n &= \boldsymbol \sigma_T \nabla u_T \cdot \mathbf n \quad \quad &\text\Sigma \end$

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

\begin
\nabla \cdot (\boldsymbol \sigma_i \nabla V_m ) + \nabla \cdot (\boldsymbol \sigma_i \nabla u_e ) = A_m \left( C_m \dfrac + I_\text(V_m,w)\right) + I_\text &\text \Omega_H\\ ex\nabla \cdot \left( (\boldsymbol \sigma_i + \boldsymbol \sigma_e) \nabla u_e \right) +\nabla \cdot(\boldsymbol \sigma_i \nabla V_m) = 0 & \text \Omega_H\\ ex\dfrac + g(V_m,w) = 0 &\text \Omega_H \\ ex\nabla \cdot (\boldsymbol \sigma_T \nabla u_T ) = 0 & \text \Omega_T\\ ex(\boldsymbol \sigma_i \nabla V_m) \cdot \mathbf n + (\boldsymbol \sigma_i \nabla u_e) \cdot \mathbf n = 0 & \text\Sigma\\ exu_e = u_T &\text\Sigma\\ ex(\boldsymbol \sigma_e \nabla u_e)\cdot \mathbf n = (\boldsymbol \sigma_T \nabla u_T) \cdot \mathbf n &\text\Sigma\\ ex(\boldsymbol \sigma_T \nabla u_T) \cdot \mathbf n = 0 &\text\Gamma_T
\end

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.

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:
$\nabla \cdot \left( (\boldsymbol \sigma_i + \boldsymbol \sigma_e) \nabla u_e \right) +\nabla \cdot(\boldsymbol \sigma_i \nabla V_m) = 0 \quad \quad \text \Omega_H.$

This way, the coupling conditions do not need to be changed, and the complete heart-torso model is composed of two different blocks:
* First the monodomain model with its usual boundary condition must be solved: \begin
\chi C_m \dfrac - \nabla \cdot \left(\boldsymbol \sigma_i \dfrac\nabla V_m\right) + \chi I_\text = I_\text &\text\Omega_H \\ ex\left(\boldsymbol \sigma_i \nabla V_m \right) \cdot \mathbf n = 0 &\text\partial \Omega_H.
\end

* Then, the coupled model that includes the computation of the extracellular potential, the torso model and the coupling conditions must be solved: $\begin \nabla \cdot \left( (\boldsymbol \sigma_i + \boldsymbol \sigma_e) \nabla u_e \right) +\nabla \cdot(\boldsymbol \sigma_i \nabla V_m) = 0 & \text \Omega_H \\ \nabla \cdot (\boldsymbol \sigma_T \nabla u_T ) = 0 & \text \Omega_T \\ \boldsymbol \sigma_T \nabla u_T \cdot \mathbf n = 0 &\text\Gamma_T \\ u_e = u_T &\text\Sigma\\ (\boldsymbol \sigma_e \nabla u_e)\cdot \mathbf n = (\boldsymbol \sigma_T \nabla u_T) \cdot \mathbf n &\text\Sigma \end$

Uncoupled heart-torso models

The fully coupled heart-torso models are very detailed models but they are also computationally expensive to solve. A possible simplification is provided by the so-called ''uncoupled assumption'' in which the heart is considered completely electrically isolated from the heart. Mathematically, this is done imposing that the current cannot flow across the epicardium, from the heart to the torso, namely
$\boldsymbol \sigma_e \nabla u_e \cdot \mathbf n = 0 \quad \quad \text\Sigma.$

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.

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:
* The bidomain model in its isolated form $\begin \nabla \cdot (\boldsymbol \sigma_i \nabla V_m ) + \nabla \cdot (\boldsymbol \sigma_i \nabla u_e ) = A_m \left( C_m \dfrac + I_\text(V_m,w)\right) + I_\text &\text \Omega_H\\ \nabla \cdot \left( (\boldsymbol \sigma_i + \boldsymbol \sigma_e) \nabla u_e \right) +\nabla \cdot(\boldsymbol \sigma_i \nabla V_m) = 0 & \text \Omega_H\\ \dfrac + g(V_m,w) = 0 &\text \Omega_H \\ (\boldsymbol \sigma_i \nabla V_m) \cdot \mathbf n + (\boldsymbol \sigma_i \nabla u_e) \cdot \mathbf n = 0 & \text\Sigma\\ \boldsymbol \sigma_e \nabla u_e \cdot \mathbf n = 0 & \text\Sigma. \end$
* The torso diffusive model in its standard formulation, with the potential continuity condition $\begin \nabla \cdot (\boldsymbol \sigma_T \nabla u_T ) = 0 & \text \Omega_T \\ \boldsymbol \sigma_T \nabla u_T \cdot \mathbf n = 0 &\text\Gamma_T \\ u_T = u_e &\text\Sigma\\ \end$

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:
* The monodomain model with its usual boundary condition: \begin
\chi C_m \dfrac - \nabla \cdot \left(\boldsymbol \sigma_i \dfrac\nabla V_m\right) + \chi I_\text = I_\text &\text\Omega_H \\ ex\left(\boldsymbol \sigma_i \nabla V_m \right) \cdot \mathbf n = 0 & \text\partial \Omega_H.
\end

* The problem to compute the extracellular potential with a boundary condition on the epicardium prescribing no intracellular current flow: \begin
\nabla \cdot \left( (\boldsymbol \sigma_i + \boldsymbol \sigma_e) \nabla u_e \right) +\nabla \cdot(\boldsymbol \sigma_i \nabla V_m) = 0 & \text \Omega_H \\ ex(\boldsymbol \sigma_i + \boldsymbol \sigma_e) \nabla u_e \cdot \mathbf n = - (\boldsymbol \sigma_i \nabla V_m ) \cdot \mathbf n &\text\Sigma
\end

* The torso diffusive model with the potential continuity boundary condition at the epicardium: $\begin \nabla \cdot (\boldsymbol \sigma_T \nabla u_T ) = 0 & \text \Omega_T \\ \boldsymbol \sigma_T \nabla u_T \cdot \mathbf n = 0 &\text\Gamma_T \\ u_T = u_e &\text\Sigma\\ \end$

Electrocardiogram computation

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
$\begin \mathrm & = u_T(L) - u_T(R) \\ \mathrm &= u_T(F) -u_T(R) \\ \mathrm &= u_T(F) - u_T(L) \\ \mathrm &= \frac(u_T(R) - u_w) \\ \mathrm &= \frac(u_T(L) - u_w) \\ \mathrm &= \frac(u_T(F) - u_w) \\ \mathrm_i &= u_T(V_i) - u_w \quad i = 1,\dots,6 \end$
where $u_W = (u_T(L) + u_T(R) + u_T(F))/3$ and $L, R, F, \_^6$ are the standard locations of the electrodes.

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.

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.

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.

Geometric torso model

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 is a useful dataset to create a three-dimensional torso model detailed with internal organs including the skeletal structure and muscles.

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.

The quasi-periodicity of the heart beat is reproduced by a three-dimensional trajectory around an attracting limit cycle in the $(x,y)$ plane. The principal peaks of the ECG, which are the P,Q,R,S and T, are described at fixed angles $\theta_P,\theta_Q,\theta_R,\theta_S \text \theta_T$, which give the following three ODEs

\begin
x' &= \alpha z - \omega y \\ exy' &= \alpha y + \omega x \\ exz' &= \textstyle - \sum_ a_i \Delta \theta_i \text\left( - /\right) - (z-z_0)
\end
with $\alpha = 1-\sqrt$, $\Delta \theta_i = (\theta -\theta_i) \bmod (2\pi)$, $\theta = \text2(y,x)$

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

See also

* Monodomain model
* Bidomain model
* Electrocardiography

References

{{reflist

Cardiac electrophysiology
Cardiac procedures
Electrodiagnosis
Electrophysiology
Mathematics in medicine
Medical tests
Partial differential equations
Mathematical modeling
Numerical analysis