Solving the Inverse Problem of Electrocardiography on the Endocardium Using a Single Layer Source

The inverse problem of electrocardiography consists in reconstructing cardiac electrical activity from given body surface electrocardiographic measurements. Despite tremendous progress in the field over the last decades, the solution of this problem in terms of electrical potentials on both epi- and the endocardial heart surfaces with acceptable accuracy remains challenging. This paper presents a novel numerical approach aimed at improving the solution quality on the endocardium. Our method exploits the solution representation in the form of electrical single layer densities on the myocardial surface. We demonstrate that this representation brings twofold benefits: first, the inverse problem can be solved for the physiologically meaningful single layer densities. Secondly, a conventional transfer matrix for electrical potentials can be split into two parts, one of which turned out to posess regularizing properties leading to improved endocardial reconstructions. The method was tested in-silico for ventricular pacings utilizing realistic CT-based heart and torso geometries. The proposed approach provided more accurate solution on the ventricular endocardium compared to the conventional potential-based solutions with Tikhonov regularization of the 0th, 1st, and 2nd orders. Furthermore, we show a uniform spatio-temporal behavior of the single layer densities over the heart surface, which could be conveniently employed in the regularization procedure.

1. Introduction

Non-invasive electrocardiographic imaging (ECGI) is a novel imaging modality which is based on numerical reconstruction of cardiac electrical activity using the so-called body surface potential maps (BSPM) and patient-specific heart and torso geometries (Ramanathan et al., 2004). The ultimate goal of the ECGI is to allow non-invasive panoramic cardiac mapping in a beat-to-beat mode, thus facilitating diagnostics and treatment planning for non-sustained, aperiodic or non-tolerable cardiac arrhythmia.

A mathematical problem underlaying ECGI is known as the inverse problem of ECG. This problem can be formulated in several ways with respect to the unknown physical values that are sought after. Its early formulation concerned pericardial (also called epicardial) potentials, which by definition neglect electrophysiological processes ongoing on the inner heart walls (Franzone et al., 1978; Rudy and Messinger-Rapport, 1988). Mathematically, this is a Cauchy problem for the Laplace equation, a classical example of ill-posed problem: even small amount of noise in the Cauchy boundary data can lead to arbitrary high errors in the solution (Kubo, 1994; Takeuchi and Yamamoto, 2008). Therefore, special regularizing algorithms must be used in order to obtain a stable solution approximation.

With this respect, significant progress has been recently made in the development of numerical algorithms for solving the inverse problem of ECG in terms of epicardial potentials, which is reflected in a constantly increasing number of clinical applications of the ECGI. This methodology was successfully used for optimizing cardiac resynchronization therapy, guiding catheter ablation of origins of focal atrial and ventricular tachycardia, detecting macro-reentrant circuits and electrical rotors in patients with reentrant ventricular tachycardia, atrial flutter and atrial fibrillation (Guillem et al., 2013; Erkapic et al., 2014; Shah et al., 2014; Dubois et al., 2015; Varma, 2015; Rodrigo et al., 2017). Moreover, Cuculich et al. demonstrated the way ECGI can provide a support for guiding non-invasive ablation of cardiac arrhythmia (Cuculich et al., 2017).

Unfortunately, epicardial imaging bears one significant limitation of potentially losing electrophysiologically relevant information about cardiac electrical activity on the cardiac endocardium and, especially, on the interventricular and interatrial septum. To overcome this drawback, one can consider reconstructing electrical potentials on both epicardial and endocardial heart surfaces. Formally, this statement results in the same Cauchy problem for the Laplace equation. Due to a more complex non-convex geometry of the epi-endocardial heart surface compared to its relatively simple “convex hull” (epicardium/pericardium), the inverse problem for endo-epicardial reconstruction becomes even worse conditioned.

Nevertheless, solution of the endo-epicadial inverse problem was employed for detection of origins of focal ventricular tachycardia (Revishvili et al., 2015; Wissner et al., 2016), determination of electrical rotors in atrial fibrillation (Metzner et al., 2017), exploring morphology of unipolar epicardial and endocardial electrograms in the right ventricular outflow tract in patients with Brugada syndrome (Rudic et al., 2016), analysis of excitation patterns in reentrant ventricular tachycardia (Tsyganov et al., 2017) and atrial flutter (Wissner et al., 2018). In these studies, a numerical algorithm based on a combination of Tikhonov and iterative regularization was used (Bokeriya et al., 2008; Kalinin, 2011).

Alternative to the potential-based statement, the problem can also be formulated in terms of surface electrical layer source models. The most prominent example of such statement is the equivalent double layer (EDL) defined on both epi- and endocardial surfaces of the heart (van Oosterom, 2014). According to the bidomain model (Tung, 1978), the EDL is proportional to the transmembrane potential when the body electrical conductivity as well as the extracellular and intracellular myocardial conductivities are considered to be isotropic and the sum of the extracellular and intracellular conductivities is equal to those of the body (Geselowitz, 1989; Kalinin et al., 2017). This electrophysiological meaning was shown to be highly beneficial for construction of ECGI-specific regularization techniques (Berger et al., 2006, 2011; van Dam et al., 2009).

In contrast to the surface EDL source model, the electrical sources inside the myocardium cannot be reconstructed unambiguously (Geselowitz, 1989; Kalinin et al., 2017). However, employement of proper regularization schemes targeting intramural transmembrane potentials or current densities was reported to overcome this rather theoretical limitation delivering promising results (see for example He et al., 2003; Skipa, 2004; Schulze et al., 2013; Wang et al., 2013; Xu et al., 2014; Zhou et al., 2016).

Overall, despite the efforts and progress made in the ECGI field, non-invasive reconstruction of the local cardiac activity on both epi- and endocardium of the heart remains a challenging task for clinical, mathematical and engineering research. Furthermore, it is evident that, irrespectively of the source model under consideration, effectiveness of Tikhonov regularization method strongly depends on the choice of a regularization operator R. In addition to the simplest option, i.e., using an identity matrix, the surface Laplacian L as well as an operator D mapping the electrical potential on the cardiac surface to its normal derivative or the transmural gradient were used for the ECGI applications (Horácek and Clements, 1997; Erem et al., 2014; Wang et al., 2016). However, the problem of an optimal choice of the regularization operator in Tikhonov regularization is still open.

In this article, we describe a novel numerical approach for treating the epi-endocardial reconstruction problem by introducing an alternative source model formulation, the single layer density. We consider this problem from three interrelated perspectives. From the numerical algebraic point of view, we introduce the involved transfer matrices associated to the boundary elements method. Furthermore, we investigate regularizing properties of the inverse single layer operator for Tikhonov regularization. Finally, we investigate spatio-temporal behavior of the single layer density source model, which can be employed in the regularization procedures.

2. Methods

In this paper we use the geometry notation reported in Figure 1. Let Γ₀ be a body surface and Γ₁ be a surface of ventricles (or atria) circumpassing both epi- and endocardial parts. Let Ω ⊂ ℝ³ be a body domain bounded from the outside by the surface Γ₀ and from the inside by surface Γ₁, with outward unit normal vectors. Let Ω_M ⊂ Ω represent the myocardial domain bounded by Γ₁ with inward unit normal vectors. Surfaces Γ₀ and Γ₁ are supposed to be sufficiently smooth.

FIGURE 1

Figure 1. Physical model underlying the inverse ECG problem.

In the physical model considered here, the electrical field is originated by the electrical sources situated in the myocardium domain Ω_M only. We neglect electrical sources in the human body domain Ω and consider the body domain as a passive volume conductor. This body domain Ω includes extracardiac organs as well as the cardiac chambers filled with the blood. Moreover, for those time moments when the ventricles/atria are in the resting state the atrial/ventricular myocardium can be also considered as a passive volume conductor.

The commonly used approach for the mathematical description of myocardial electrical activity is based on so the called bidomain model. Within this model, myocardial tissue is considered to consist of homogenized intracellular and extracellular spaces. The set of accompanying differential equations establishes the relationship between the intra- and extracellular potentials taking into account cell membrane properties and ionic concentrations (see e.g., Tung, 1978; Bourgault et al., 2009). The electrical conductivity in the intra- and extracellular media are provided in tensor form reflecting faster excitation propagation along the myocardial fibers than across them. In this work, we employed the bidomain model with isotropic cardiac electrical conductivity values and a homogeneous torso model to simulate the electrical potentials throughout the geometry volume for ectopic ventricular stimuli. The potentials obtained on the body surface, the body surface potential maps (BSPM), were then used as the input for validating the proposed approach of non-invasive reconstruction of cardiac electrical activity. This, so-called forward problem of ECG, was solved with the CHASTE software (Mirams et al., 2013), which provided the reference transmembrane potentials in the heart and electrical potentials in the whole geometry volume.

For biological tissues frequencies under the conventional assumptions the Maxwell equations, describing propagation of electromagnetic fields in the body volume conductor, can be simplified to the quasi-static form (Gulrajani, 1998). It allows temporal separation of the cardiac sources, meaning that one can solve the inverse problem of ECG for distinct time instants t₀, t₁, …, t_M independently from each other.

Furthermore, the inverse problem of ECG can be treated using direct and indirect regularization approaches. The direct way is the computation of the harmonic function value on Γ₁ without considering myocardial electrical sources. The indirect way consists of presenting the electrical potential on Γ₁ as a potential of the myocardial sources. It is well known that an endless number of the sources distributions in the myocardium domain can generate the same potential in the passive volume conductor domain. Therefore, they use an “effective” unique representation of the electrical sources in form of sources on the myocardial surface. In this article we consider the direct and the indirect ways for numerical solving the inverse electrocardiography problem.

2.1. Computational Method for the Inverse Potential Problem: A Conventional Approach

The inverse problem of ECG in terms of electrical potentials for the geometry depicted in Figure 2 reads to find a function u(x) in $\overline{\Omega }$ such that

$\begin{array}{cc}\begin{array}{r}\Delta u(x)=0,x\in \Omega ,\end{array}& (1)\end{array}$

$\begin{array}{cc}\begin{array}{r}u(x)=\phi (x),x\in {\text{Γ}}_0,\end{array}& (2)\end{array}$

$\begin{array}{cc}\begin{array}{r}\frac{\partial u(x)}{\partial n}=0,x\in {\text{Γ}}_0,\end{array}& (3)\end{array}$

where φ(x) is the measured BSPM. Problem (1)–(3) is known as the Cauchy problem for the Laplace equation. Its solution is unique, however, the problem is ill-posed: even a small amount of noise in the boundary conditions can lead to an arbitrary high error in the solution.

FIGURE 2

Figure 2. Schematic geometric relationships of the inverse potential problem in the internal statement. Ω is the passive volume conductor domain, Ω_M is the myocardial domain, Γ₀ is the body surface, Γ₁ is the myocardial surface (endo- and epicardial surface), $\overrightarrow{\text{n}}$ its unit normal vector directed inward, P_i, i = 1..N₀+N₁ are collocation points used in direct boundary element method, N₀ is the number of collocation points on the Γ₀, N₁ is the number of collocation points on Γ₁.

To solve the problem a direct boundary element method (BEM) can be used. Accordingly, the problem (1)–(3) is reformulated as a boundary integral equation based on the third Green's identity: for any given point P ∈ Γ₀ ∪ Γ₁ and harmonic function u in domain Ω

$\begin{array}{cc}\begin{array}{r}c(P)u(P)+{\displaystyle \underset{{\text{Γ}}_0\cup {\text{Γ}}_1}{\int }}u(Q)\frac{\partial G(P,Q)}{\partial n}d\text{Γ}={\displaystyle \underset{{\text{Γ}}_0\cup {\text{Γ}}_1}{\int }}\frac{\partial u(Q)}{\partial n}G(P,Q)d\text{Γ},\end{array}& (4)\end{array}$

where Q ∈ dΓ is the integration variable and

$\begin{array}{cc}\begin{array}{r}G(P,Q)=\frac{1}{|P-Q|}\end{array}& (5)\end{array}$

is the inverse Euclidean distance between P and Q, c(P) is the solid angle at the point P.

Next step is to establish a system of linear algebra equations suitable for numerical calculations from the continues statement (4). This step is called discretization. We use the following discretization scheme: (a) approximation of the surfaces Γ₀ and Γ₁ by the triangular meshes, (b) approximation of the functions u(x) and $\frac{\partial u(x)}{\partial n}$ by series of linear basis functions, and (c) computation of the single and double layer type integrals over basis functions. Computation of such integrals over basis function is most important step. Fortunately, it is well studied, see e.g., Dunavant (1985), Davey and Hinduja (1989), and van Oosterom (2012). Full details of the discretization process are described in the Appendix A.

After the discretization we get the following system of linear equations:

$\begin{array}{cc}\begin{array}{r}\begin{array}{c}{H}_{00}{u}_0+H_{01}{u}_1=G_{01}{q}_1\\ H_{10}{u}_0+H_{11}{u}_1=G_{11}{q}_1\end{array}\end{array}& (6)\end{array}$

where u₀, u₁ are electrical potentials on the surfaces Γ₀ and Γ₁ respectively, q₁ is the normal derivative of the electrical potential on the surface Γ₁, matrices G_ij arise from the discretization of the surface integrals corresponding to the single layer

$\begin{array}{cc}\begin{array}{r}{\displaystyle \underset{{\text{Γ}}_j}{\int }}\frac{\partial u(Q)}{\partial n}G(P_{{\text{Γ}}_i},Q)d{\text{Γ}}_Q,\end{array}& (7)\end{array}$

while matrices H_ij arise from the discretization of the surface integrals corresponding to the double layer

$\begin{array}{cc}\begin{array}{r}{\displaystyle \underset{{\text{Γ}}_j}{\int }}u(Q)\frac{\partial G(P_{{\text{Γ}}_i},Q)}{\partial n}d{\text{Γ}}_Q,\end{array}& (8)\end{array}$

Finally, i is the index of the surface containing the point P, j is that of the surface containing Q.

In Figure 3 we provide the plots of the singular values' decay of the matrices above. In agreement with the boundary element theory, matrices G₁₁, H₀₀ and H₁₁ are well-conditioned and can be inverted without regularization.

FIGURE 3

Figure 3. SVD plot of the considered matrices, (A) is ill-conditioned matrices, (B) is well-conditioned matrices.

From the system (6) we can get the transfer matrix u₁→u₀ (relating EP to BSPM, EP stands for endo- and epicardial potentials):

$\begin{array}{cc}\begin{array}{r}(-H_{01}+G_{01}{G}_{11}^{-1}{H}_{11})u_1=(H_{00}-G_{01}{G}_{11}^{-1}{H}_{10})u_0,\end{array}& (9)\end{array}$

in short form

$\begin{array}{cc}\begin{array}{r}A{u}_1=f_0,\end{array}& (10)\end{array}$

where f₀ is known right-hand side of Equation (9). In the following, we call the inverse ECG problem statement (1)–(3) the internal statement and Equation (9) the internal equation.

In order to find the normal derivative of the potential on Γ₁ let us introduce the Dirichlet-Neumann mapping matrix u₁ → q₁ relating EP to its normal derivative on the heart. This matrix can be derived from the system (6) in the form

$\begin{array}{cc}\begin{array}{r}{(G_{11}-H_{10}{H}_{00}^{-1}{G}_{01})}^{-1}(H_{11}-H_{10}{H}_{00}^{-1}{H}_{01})u_1=q_1,\end{array}& (11)\end{array}$

or

$\begin{array}{cc}\begin{array}{r}D{u}_1=q_1\end{array}& (12)\end{array}$

Note that representations (9) and (11) require inversion of well-conditioned matrices only. Matrices A and D are well known in the literature (e.g., see Yun et al., 1997; Gulrajani, 1998).

Matrix A is ill-conditioned, therefore the numerical solution of Equation (10) requires suitable regularization techniques. The commonly used approach is the Tikhonov regularization method:

$\begin{array}{cc}\begin{array}{r}{u}_1^{\lambda }=\text{arg min}(\parallel A{u}_1-f_0\parallel {}_2^2+{\lambda }^2\parallel R{u}_1\parallel {}_2^2),\end{array}& (13)\end{array}$

where $u_1^{\lambda }$ is the regularized solution, λ² is the regularization parameter and R is the regularization operator. Minimization problem (13) has the closed-form solution:

$\begin{array}{cc}\begin{array}{r}{u}_1^{\lambda }={(A^{T}A+{\lambda }^2{R}^{T}R)}^{-1}{A}^T{f}_0.\end{array}& (14)\end{array}$

Regularization operator R can be taken, for example, as R = I₁₁ (identity matrix) for 0 order, R = D for 1th order and as a Laplace-Beltrami (“surface Laplacian”) operator Lu₁ = Δ_Γ1u₁ (see for example Huiskamp, 1991) for the 2nd order Tikhonov regularization.

2.2. Computational Method for the Inverse Potential Problem: A Single Layer Approach

In this section we will formulate an alternative representation of the u₁ → u₀ transfer matrix and propose a new statement of the inverse problem in terms of the equivalent single layer (ESL). Geometry notations for this statement are depicted on Figure 4.

FIGURE 4

Figure 4. Schematic geometric relationships of the inverse potential problem in the external statement. Ω_M is the myocardial domain, Γ₀ is the body surface, Γ₁ is the myocardial surface (endo- and epicardial surface), $\overrightarrow{\text{n}}$ is unit normal vector, P_i, i = 1..N₀+N₁ are collocation points used in direct boundary element method, N₀ is the number of collocation points on the Γ₀, N₁ is the number of collocation points on Γ₁.

Although the cardiac electrical potential u(x) is not a harmonic function in the domain Ω_M, we can define another function û harmonic in Ω_M whose boundary values on Γ₁ are equal to those of u(x), the solution of problem (1)–(3). The function û can be defined as the unique solution of the following Dirichlet problem for the Laplace equation:

$\begin{array}{cc}\begin{array}{r}\Delta û(x)=0,x\in {\Omega }_M,\end{array}& (15)\end{array}$

$\begin{array}{cc}\begin{array}{r}û(x)=u(x),x\in {\text{Γ}}_1.\end{array}& (16)\end{array}$

The idea under this definition is following. With the function û(x) harmonic in the domain Ω_M we force the electric sources to be only on the endo- and epicardial surface. And with the condition û(x) = u(x) on Γ₁ we can link such sources densities with the actual electrical potential u(x). To derive such relationships we use boundary element method.

For any given point P ∈ Γ₁ and harmonic function û in domain Ω_M the third Green's identity gives us the following boundary integral equation

$\begin{array}{cc}\begin{array}{l}ĉ(P)û(P)+{\displaystyle \underset{{\text{Γ}}_1}{\int }}û(Q)\frac{\partial G(P,Q)}{\partial n}d{\text{Γ}}_Q\\ ={\displaystyle \underset{{\text{Γ}}_1}{\int }}\frac{\partial û(Q)}{\partial n}G(P,Q)d{\text{Γ}}_Q,P\in {\text{Γ}}_1.\end{array}& (17)\end{array}$

Moreover, for any point P ∈ Γ₀ and harmonic function û in domain Ω_M the third Green's identity give us the following boundary integral equation

$\begin{array}{cc}\begin{array}{r}{\displaystyle \underset{{\text{Γ}}_1}{\int }}û(Q)\frac{\partial G(P,Q)}{\partial n}d{\text{Γ}}_Q={\displaystyle \underset{{\text{Γ}}_1}{\int }}\frac{\partial û(Q)}{\partial n}G(P,Q)d{\text{Γ}}_Q,P\in {\text{Γ}}_0,\end{array}& (18)\end{array}$

Applying the same discretization as in section 2.1 we get the following algebraic system:

$\begin{array}{cc}\begin{array}{r}\begin{array}{c}{Ĥ}_{11}{û}_1={Ĝ}_{11}{\widehat{q}}_1\\ {Ĥ}_{01}{û}_1={Ĝ}_{01}{\widehat{q}}_1.\end{array}\end{array}& (19)\end{array}$

where û₁ is a vector containing values of function û(x) at the points on the surface Γ₁, ${\widehat{q}}_1$ is a vector containing values of function $\frac{\partial û(x)}{\partial n}$ on Γ₁, matrices Ĝ_ij arise from the discretization of the surface integrals corresponding to the single layer

$\begin{array}{cc}\begin{array}{r}{\displaystyle \underset{{\text{Γ}}_j}{\int }}\frac{\partial û(Q)}{\partial n}G(P_{{\text{Γ}}_i},Q)d{\text{Γ}}_Q,\end{array}& (20)\end{array}$

matrices Ĥ_ij arise from the discretization of the surface integrals corresponding to the double layer

$\begin{array}{cc}\begin{array}{r}{\displaystyle \underset{{\text{Γ}}_j}{\int }}û(Q)\frac{\partial G(P_{{\text{Γ}}_i},Q)}{\partial n}d{\text{Γ}}_Q,\end{array}& (21)\end{array}$

i is the index of the surface with the fixed point P, j is the index of the surface with points of integration Q.

We can express the unknown variable ${\widehat{q}}_1$ from the first equation and obtain a new matrix-vector identity for the variable û₁

$\begin{array}{cc}\begin{array}{r}{Ĥ}_{01}{û}_1={Ĝ}_{01}{Ĝ}_{11}^{-1}{Ĥ}_{11}{û}_1\end{array}& (22)\end{array}$

Using (16) we can write (22) as

$\begin{array}{cc}\begin{array}{r}{Ĥ}_{01}{u}_1={Ĝ}_{01}{Ĝ}_{11}^{-1}{Ĥ}_{11}{u}_1.\end{array}& (23)\end{array}$

Let us compare matrices Ĝ₀₁, Ĝ₁₁, Ĥ₀₁, and Ĥ₁₁ with matrices G₀₁, G₁₁, H₀₁, and H₁₁ defined in section 2.1. All matrices are determined only by the same surfaces Γ₀ and Γ₁. However, normal vectors to the surface Γ₁ were directed outwards related to the domain Ω and inwards related to the domain Ω_M. Taking into account these facts it is easy to demonstrate that

$\begin{array}{cc}\begin{array}{r}{Ĝ}_{01}=G_{01},\end{array}& (24)\end{array}$

$\begin{array}{cc}\begin{array}{r}{Ĝ}_{11}=G_{11},\end{array}& (25)\end{array}$

$\begin{array}{cc}\begin{array}{r}{Ĥ}_{01}=H_{01},\end{array}& (26)\end{array}$

$\begin{array}{cc}\begin{array}{r}{Ĥ}_{11}=H_{11}-4\pi I_{11},\end{array}& (27)\end{array}$

where I₁₁ is the identity matrix.

Now we can rewrite (23) using the matrix defined in section 2.1:

$\begin{array}{cc}\begin{array}{r}{H}_{01}{u}_1=G_{01}{G}_{11}^{-1}(H_{11}-4\pi I_{11})u_1\end{array}& (28)\end{array}$

or

$\begin{array}{cc}\begin{array}{r}4\pi G_{01}{G}_{11}^{-1}{u}_1=(-H_{01}+G_{01}{G}_{11}^{-1}{H}_{11})u_1\end{array}& (29)\end{array}$

Comparing (29) with (9) we can see that the right hand side of the Equation (29) is the same as the left hand side of the Equation (9). Therefore

$\begin{array}{cc}\begin{array}{r}4\pi G_{01}{G}_{11}^{-1}{u}_1=(H_{00}-G_{01}{G}_{11}^{-1}{H}_{10})u_0\end{array}& (30)\end{array}$

or in short form

$\begin{array}{cc}\begin{array}{r}4\pi G_{01}{G}_{11}^{-1}{u}_1=f_0.\end{array}& (31)\end{array}$

Equation (31) is a new representation of the u₁ → u₀ transfer matrix. In this article we call this approach the external statement of the inverse problem and transfer matrix (30) the external form of the transfer matrix.

Matrix G₁₁ is well-conditioned (see Figure 3), so we can define a new function

$\begin{array}{cc}\begin{array}{r}{w}_1\equiv G_{11}^{-1}{u}_1\end{array}& (32)\end{array}$

and new equation for the function w₁

$\begin{array}{cc}\begin{array}{r}4\pi G_{01}{w}_1=f_0.\end{array}& (33)\end{array}$

Matrices G₁₁ and G₀₁ are discretizations of the single layer integral operators. Therefore, the function w₁ has a physical meaning of electrical sources in form of single layer on the myocardial surface. In this paper we call function w₁ the equivalent single layer (ESL) density. Equation (33) allows us to solve the inverse ECG problem in terms of the ESL.

Furthermore, we propose two methods for regularizing the inverse potential problem. The first method consists of the ESL computation by solving the Equation (33) and reconstruction of the potential u₁ from the obtained ESL by formula u₁ = G₁₁w₁.

Matrix G₀₁ is ill-conditioned, therefore the numerical solution of Equation (33) requires suitable regularization algorithms. The Tikhonov regularization method of 0th order consists in solving

$\begin{array}{cc}\begin{array}{r}{w}_1^{\lambda }=\text{arg min}(4\pi \parallel G_{01}{w}_1-f_0\parallel {}_2^2+{\lambda }^2\parallel w_1\parallel {}_2^2),\end{array}& (34)\end{array}$

whose solution reads

$\begin{array}{cc}\begin{array}{r}{w}_1^{\lambda }={(G_{01}^T{G}_{01}+{\lambda }^2{I}_{11})}^{-1}{G}_{01}^T{f}_0,\end{array}& (35)\end{array}$

and next we compute

$\begin{array}{cc}\begin{array}{r}{u}_1^{\lambda }=G_{11}{w}_1^{\lambda },\end{array}& (36)\end{array}$

where $w_1^{\lambda }$, $u_1^{\lambda }$ are the regularized solutions in terms of the ESL and potentials respectively, λ² is the regularization parameter.

The second method relies on solving the transfer Equation (10) using Tikhonov regularization method with constraints on the ESL, i.e., using inverse ESL matrix $G_{11}^{-1}$ as regularization operator:

$\begin{array}{cc}\begin{array}{r}{u}_1^{\lambda }=\text{arg min}(\parallel A{u}_1-f_0\parallel {}_2^2+{\lambda }^2\parallel w_1\parallel {}_2^2),\end{array}& (37)\end{array}$

or

$\begin{array}{cc}\begin{array}{r}{u}_1^{\lambda }=\text{arg min}(\parallel A{u}_1-f_0\parallel {}_2^2+{\lambda }^2\parallel G_{11}^{-1}{u}_1\parallel {}_2^2),\end{array}& (38)\end{array}$

with solution:

$\begin{array}{cc}\begin{array}{r}{u}_1^{\lambda }={(A^{T}A+{\lambda }^2{G_{11}^{-1}}^T{G}_{11}^{-1})}^{-1}{A}^T{f}_0.\end{array}& (39)\end{array}$

2.3. Experimental Methods and Evaluation Protocols

Accuracy of the numerical algorithms for solving the inverse problem of ECG was tested on realistic in-silico data of cardiac electrical activity. Computer tomography (CT) scans of three patients were used for obtaining the personalized anatomy of the torso and heart surfaces.

These patients were examined in Bakulev National Medical Research Center of Cardiovascular Surgery (Moscow, Russia) in 2017 year. The study was performed in accordance with principles of the Declaration of Helsinki. The protocol of the study was approved by local ethics committee of Bakulev Research Center. All patients provided written informed consent to the CT scanning procedures and agreed to data retrieval, analysis and publication.

A patient-specific anatomical model was comprised of a biventricular cardiac model and a homogeneous torso model. Data processing steps included segmentation of the torso and the ventricles CT volume models, generation of the polygonal surface meshes (“Amycard 01 C” software, EP Solutions SA) and creation of tetrahedral final element meshes (“Gmsh” software Geuzaine and Remacle, 2009).

The bidomain model with the strong heart-torso coupling and TNNP cellular model (ten Tusscher et al., 2004) were used for simulation of a myocardium electrical activity. Transmural and apico-basal cellular heterogeneities were simulated using the approaches proposed in Keller et al. (2012) and ten Tusscher and Panfilov (2006), respectively. Simulations of cardiac electrical activity were performed with the methodology described in Ushenin et al. (2017) using the Cardiac CHASTE software (Mirams et al., 2013). In each anatomical patient model, these finite-element calculations resulted in the transmembrane potentials in the myocardial volume, while the electrical potentials were simulated at each node of the tetrahedral meshes. The torso surface potentials were used as the input for testing the proposed inverse routines.

We simulated three focal type electrical activation patterns. The focal origins were in the lateral wall of the left ventricle (LV) for the patient heart 1, in the right ventricular apex (RVA) for the patient heart 2 and in the right ventricular outflow tract (RVOT) for the patient heart 3. A “virtual” rectangular current pulse of −50μA with a duration of 3 ms applied to an area of 6 mm was used for initiation of a cardiac excitation. “Virtual” action potential signals and local unipolar electrograms on the ventricular surface as well as the body surface ECG were obtained as a result of the simulation.

The first part of the evaluation protocol included computation of the ESL (function w₁) from the simulated electrical potential u₁ on the myocardial surface according to the formula (32) and comparison of this function's morphology in space and time with the transmembrane action potentials.

The second part of the evaluation protocol included an actual testing of the proposed algorithms for the solution of the inverse problem. For this, the BSPM were distorted by an additive Gaussian noise of 50 dB SNR and used for reconstruction of the electrical potential on the endocardial and epicardial surface of the ventricles. The reconstructed local unipolar electrograms on the myocardial surface were compared to the references electrograms obtained from the simulations.

We used the following metrics to estimate solution quality:

$\begin{array}{cc}r{e}_x=\frac{1}{M}\sum _{i=1}^M\frac{\sqrt{{\sum }_{j=1}^{N_1}{(u_1(x_j,t_i)-u_1^{num}(x_j,t_i))}^2}}{\sqrt{{\sum }_{j=1}^{N_1}{u}_1{(x_j,t_i)}^2}}& (40)\end{array}$

$\begin{array}{cc}\begin{array}{r}c{c}_t=\frac{1}{N_1}{\displaystyle \sum _{i=1}^{N_1}}cc(u_1(x_i,t),u_1^{num}(x_i,t)),\end{array}& (41)\end{array}$

where u₁(x, t) is the given solution, $u_1^{num}(x,t)$ is the inverse reconstruction, cc(·, ·) – is the correlation coefficient, M is the number of time instances, N₁ number of nodes in the heart triangular mesh.

The quantity re_x provides the spatial error of the solution for each instant t_i, whereas cc_t shows the accuracy of the reconstructed electrogram at each node on the heart mesh. These errors are calculated separately for epicardial and endocardial surfaces of the heart ventricles.

We tested the following inverse numerical scheme:

1. Tikhonov solution of the equation Au₁ = f₀ with the conventional EP transfer matrix A and three types of constraints: 0th order $\parallel u_1\parallel {}_2^2\equiv \parallel I_{11}{u}_1\parallel {}_2^2$, 1st order $\parallel \frac{\partial u_1}{\partial n}\parallel {}_2^2\equiv \parallel D{u}_1\parallel {}_{L_2}$ and 2nd order $\parallel {\Delta }_{{\Gamma }_1}{u}_1\parallel {}_2^2\equiv \parallel L{u}_1\parallel {}_2^2$.

2. Tikhonov solution of the equation 4πG₀₁w₁ = f₀, u₁ = G₁₁w₁ with the ESL transfer matrix G₀₁ and the constraint $\parallel w_1\parallel {}_2^2\equiv \parallel I_{11}{w}_1\parallel {}_2^2$, see Equation (34).

3. Combination of these two approaches, i.e., Tikhonov solution of the equation Au₁ = f₀ with conventional EP matrix A but with the new type of constraint $\parallel w_1\parallel {}_2^2\equiv \parallel G_{11}^{-1}{u}_1\parallel {}_2^2$ for the ESL, see Equations (37), (38).

We found the value of the regularization parameter λ thanks to the L-Curve method (Hansen, 2000).

3. Results

3.1. Properties of the ESL (Function w₁)

Figure 5 shows the distribution of the ESL (function w₁(x, t)) on the heart surface for several time instants of the cardiocycle. Figure 6 shows ESL as time signals at several points of the ventricular surface. The ESL in space as well as in time domains looks like a sparsed function and well reflects the fronts of myocardial depolarization and repolarization.

FIGURE 5

Figure 5. Distribution of the ESL density (function w₁) on the heart surface for the fixed time moments of cardiocycle. Cardiac excitation was initiated in the apical area. LV, left ventricle; RV, right ventricle; RVOT, right ventricle outflow tract.

FIGURE 6

Figure 6. Example of ESL and transmembrane action potential signals (simulation data). Cardiac excitation was initiated in the apical area. (A) is locations of the points where the signal was computed, (B) is transmembrane action potentials (left panel) and the ESL signals (right panel), (C) is merged transmembrane action potentials and ESL signals. Notation mC/m² is the millicoulomb per square meter, the unit for an electrical charge density.

The ESL signals shape has a form of bipolar spike associated with time moments of myocardial depolarization and low-amplitude wave associated with myocardium repolarization. Depolarization part of the signal has similar morphology at all points of the myocardial surface with the first peak being positive and the second one negative (see Figure 6B). Duration of depolarization spike varied from 12 to 22 ms, its magnitude varied from 0.5 to 1.2 mC/m².

The repolarization waves were positive in 91% of ESL signals and negative in 2% of the signals. In 7% of the signals the repolarization waves were biphasic with the first negative half-wave. Duration of the repolarization waves varied from 18 to 27 ms, their magnitude varied from 0.05 to 0.2 mC/m².

The zero-crossing value between positive and negative peaks of the signals matched with depolarization time moment (see Figure 6C). The mean difference between that zero point of the single layer density signal and maximal slope of transmembrane action potential (TMP) up-stroke was 1.3±2.4 ms.

The time moment of maximum positive monophasic repolarization waves in ESL signals corresponded to 91±4% level of TMP repolarization (see Figure 6C). Biphasic and negative repolarization waves were also associated with repolarization phase of TMP, but reliable identification of connections between the moment of their appearance and the level of myocardial repolarization requires more data and further analysis. Therefore, ESL signals can be potentially used for detection of depolarization and repolarization of the myocardium. However, development of this method requires further investigations.

3.2. Accuracy of the Inverse Solutions

Table 1 shows results of the numerical experiments. The first column in the table shows the type of equation used, the second column shows the regularization constraint and the last columns show the re_x and cc_t metric values (see (40), (41)) separately for the LV, RVA and RVOT simulations and epicardial and endocardial surfaces.

TABLE 1

Table 1. Results of the inverse reconstruction.

Electrical signals from the several epicardial and endocardial sites of the ventricles were selected as representative examples in order to provide visual evaluation of their morphology (see Figure 7). The center of LV lateral wall (Epi LV lateral), the apical region (Epi LV apex), anterior-lateral zone of the RVOT (Epi RVOT) and the center of RV lateral wall (Epi RV lateral) were taken on the epicardial surface. The center of LV lateral wall (Endo LV lateral), the central zones of the left (Septum LV) and the right sides (Septum RV) of the ventrical septum and the center of RV lateral wall (Endo RV lateral) were taken on the endocardial surface.

FIGURE 7

Figure 7. Given (red curves) and inverse reconstructed electrograms (blue curves) in different point of epicardial and endocardial surface for the simulation of the pacing from the RVOT. (A) Is the reconstruction with the conventional EP transfer matrix A and Tikhonov regularization of 2^nd order, (B) is the reconstruction with the ESL transfer matrix G₀₁ and Tikhonov regularization of 0 order.

In Figure 7 we present the results of the inverse reconstruction of the electrograms in the defined above points on the heart surface. Figure 7A shows results of reconstruction with the conventional EP transfer matrix A and Tikhonov regularization of the 2nd order. Figure 7B shows results of reconstruction with the ESL transfer matrix G₀₁ and Tikhonov regularization of 0th order.

All algorithms demonstrated similar accuracy on the epicardial surface, but their accuracy on the endocardial surface was significantly different.

On the endocardial surface algorithm in terms of EP with conventional transfer matrix A and 0th order regularization demonstrated poor accuracy. The reconstructed electrograms have near-zero magnitude and poor correlation with the references electrograms. Algorithms with the 1st order regularization provided poor accuracy in terms of relative error, but the reconstructed electrograms correlated better with the reference signals (see Table 1). Algorithm with the 2nd order regularization showed slightly better results in comparison to 1^st order. However, this algorithm did not allow to reconstruct the electrograms morphology with acceptable quality. In particular, reconstructed electrograms at the Endo LV lateral, RV lateral and Septum RV sites has opposite polarity at QRS part and at the Endo LV lateral site has opposite polarity at ST part than the reference electrograms (see Figure 7A). Detailed metrics values are given in the Table 1.

Both algorithms in terms of ESL significantly improved the accuracy on the endocardial surface. These algorithms showed low relative error and high correlation coefficient. The morphology of endocardial electrograms on the LV and RV lateral walls of the ventricles as well as at RVOT were reconstructed with enough accuracy. The electrograms on the LV and RV ventricular septum were more smoothed, but the basic elements of their morphology (polarities of the electrogram waves) were reconstructed correctly.

Activation and recovery times are commonly used in the clinical practice as one of the important outputs of non-invasive cardiac imaging. Some numerical results in detection of activation and repolarization times from electrograms reconstructed by the proposed ESL algorithm are given in the Supplementary Material 1.

4. Discussion

Non-invasive cardiac electrical mapping on both epi- and edocardial surfaces of the heart can provide more detailed information about cardiac electrical activity. However, this methodology is more challenging compared to the non-invasive epicardial mapping. In previous works the problem of endo-epicardial mapping was attacked in two directions. The former was to extend the inverse electrocardiography problem in terms of epicardial potentials to a problem in terms of epicardial and endocardial potentials. The second one was to reconstruct cardiac electrical activity on the epicardium and endocardium in terms of the EDL or in terms of intramural “equivalent” electrical sources related to cardiac transmembrane potentials.

In this article, we introduced a novel representation of cardiac sources in terms of the ESL potential. This approach, in a sense, combines these two directions. Utilizing EDL for representation of cardiac electrical activity was motivated by the following reasons. First, it is well known that the electrical potential on the cardiac surface can be understood as a sum of two components: so-called “near field,” reflecting local myocardial electrical activity and so-called “far field,” which is generated by electrical sources at remote segments of the heart. This fact leads to certain difficulties in interpretation of local unipolar electrograms with respect to depolarization and repolarization times. In contrast, EDL allows detecting the local electrical activity of the myocardium with greater precision. Secondly, EDL signals have higher level of regularity in comparison to local unipolar electrograms. This trait of EDL may provide additional opportunities for regularization of the inverse problem.

Our results of the presented in-silico experiments showed that the ESL representation of cardiac electrical activity has also some attractive properties. ESL density correlated well with the local electrical activity of the myocardium. ESL density as time signals can be used for detection of activation and recovery times, calculating activation and recovery intervals and reconstruction of activation and recovery sequences. Note that detection of repolarization sequences by local electrograms meets some methodological difficulties (Cluitmans et al., 2017). However, repolarization abnormalities can be an important substrate of reentrant atrial and ventricular arrhythmias. Therefore, possible application of ESL for detection of repolarization abnormalities seems to be promising.

Moreover, ESL density is a temporally localized function exhibiting very similar morphology for all ventricular sites. These ESL features can be potentially used for narrowing down the set of admissible solutions in construction of regularization methods.

Note, these results were obtained for the myocardial model with isotropic electrical conductivities. The proposed approach for ESL computation does not require assumptions of myocardium anisotropy. Therefore, it can be translated directly to the more realistic anisotropic model of the myocardium. However, ESL been computed this way may slightly differ from the “physical” single layer density in case of the media with anisotropic electrical conductivity. To emphasize this fact, we used the term “equivalent” single layer density (ESL). However, we suppose that investigation of electrophysiological meaning of representation of cardiac electrical activity in form of electrical single layer for the more realistic anisotropic model requires more precise mathematical definition of electrical single layer and more complex algorithm for its computation. We address this task to further research.

The most common discretization method for the inverse potential problem, i.e., for reducing the boundary value problem for Laplace equation to a system of linear algebraic equations is BEM (Yun et al., 1997; Gulrajani, 1998 and description in section 2.1). It this work we proposed an alternative BEM scheme for assembling the transfer matrix, which is closely related to the ESL representation of the cardiac electrical field. Though we considered this method in the context of the endo-epicardial potential inverse problem, this approach can be used for the reconstruction of electrical potentials on the epicardial surface only, in contrast to the EDL.

Moreover, the presented derivation of the ESL transfer matrix allowed to identify the intrinsic structure of the conventional one, that can be split into two matrices: a well-conditioned G₁₁ and ill-conditioned G₀₁, whereby their elements are the inverse euclidean distances not depending on the solid angles involved by EDL computations.

A technical benefit of the splitting lies in this simple structure allowing greater numerical precision of the two matrices. In particular, the novel way for construction of the transfer matrix does not require calculation of normal vectors, thus eliminating possible mesh-related artifacts. Furthermore, a simple structure of ill-conditioned matrix G₀₁ provides an alternative basis for regularization approaches.

The above results were obtained for the simplified torso model with homogeneous electrical conductivities. We suppose that the translation of the reconstruction algorithm to the clinical practice requires more realistic human torso model with different electrical conductivities of the internal organs. In that cases the structure of the transfer matrix for the inverse problem in terms of ESL will be more complex. Identification of its structure is a task for further research.

We also presented a two-step method for solving the inverse potential problem including computation of the ESL density as an intermediate step. This method has some formal similarities with the method of fundamental solutions (MFS). The MFS was proposed for solving the inverse problem in the epicardial statement and showed promising results (Wang and Rudy, 2006). Briefly, MFS is based on computation of values of “virtual” point electrical sources placed outside of the domain of interest and subsequent computation of the cardiac electrical potential as a linear combination of these electrical sources. The MFS also allows usage of meshless construction of the transfer matrix. However, in contrast to the ESL on the myocardial surface, the electrical sources do not have a physiological meaning, i.e., they cannot be used for evaluating the local electrical activity of the myocardium. Next, in contrast to the ESL matrix G₁₁, the MFS matrix mapping electrical source values to the electrical potential on the cardiac surface is ill-conditioned. Thus, this matrix cannot be used as a regularization operator.

The most interesting and unexpected result was, in our opinion, obtained in the investigation of regularization effects of the ESL density on the heart. We first tested the 0th order Tikhonov regularization for the inverse problem in the ESL statement. Then, we used $G_{11}^{-1}$ as a regularization operator for the conventional EP statement of the inverse problem. We found that both schemes provided significantly more accurate solutions compared to the other Tikhonov regularizations.

Still, the regularizing properties of the inverse single layer operator require further investigation and theoretical explanation. It is worth noting that the ESL regularization increased the accuracy predominately on the endocardial surface of the heart. We can hypothesize that independence of the ESL from solid angles contributed to the accuracy increase for the complex “w-shape” geometry of the endocardial surface. Nevertheless, a detailed mathematical interpretation of the obtained results as well as a more general challenge of developing an optimal regularization operator for the inverse ECG problem is a subject for further investigations.

As our conceptual approach targeted an improved reconstruction of the electrogram itself, further efforts should be undertaken to quantify ESL advantages in estimation of derived parameters, such as activation times, frequency maps etc. The presented work focused solely on a mathematical description of the novel ESL formulation and its simulation-based proof-of-concept, making thorough in silico and clinical evaluation needed in order to translate our findings into practical benefits.

5. Conclusions

In this article, we proposed a novel statement of the inverse problem of ECG which is based on a representation of the electrical potential on the cardiac surface as ESL on the same surface. The results of in-silico experiments using personalized cardiac models demonstrated that the introduced ESL density well correlates with local electrical activity of the myocardium.

The reconstruction method was considered in two basic versions. The first version included assembling of a transfer matrix mapping the ESL density to the body surface potentials using the BEM and solving the matrix equation with Tikhonov regularization of 0th order. The second version used the conventional transfer matrix mapping EP on the cardiac surface to the BSPM and then applied Tikhonov regularization imposing constraints based on the single layer operator on the heart. The results demonstrated that both versions provided more accurate solution on the ventricular endocardial surface compared to the classical approach with Tikhonov regularization of 0th, 1st, and 2nd orders.

The proposed modifications in the solution scheme may improve non-invasive reconstruction of cardiac electrical activities on the endocardial part of the heart.

6. Limitations

We used only limited numerical simulations cases for testing feasibility of the proposed method. In our future work, we intend to extensively study performance of the presented approach in both in silico and clinical setups. Forward and inverse simulations were performed for a homogeneous torso model with an isotropic heart. While heart anisotropy affects both heart and torso surface potentials, it does not influence the relationship between the ESL source model and EP on the heart surface. The ESL variable obtained this way may slightly differ from the “physical” single layer density in case of the media with anisotropic electrical conductivity. To emphasize this fact, we used the term “equivalent” single layer density.

Only ventricular focal activation patterns were considered in the present work. Further studies should include more complex excitation propagation patterns as well as a comparative analysis between the derived (e.g., activation times) clinically relevant parameters.

Author Contributions

All authors contributed equally to this manuscript including development of algorithms, data processing and manuscript writing.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

The authors would like to thank Prof. Quarteroni for the thorough reading and his constructive comments and suggestions, Dr. Sopov, MD, PhD and Dr. Magomedova, MD for the help in patient's examinations.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fphys.2019.00058/full#supplementary-material

References

Berger, T., Fischer, G., Pfeifer, B., Modre, R., Hanser, F., Trieb, T., et al. (2006). Single-beat noninvasive imaging of cardiac electrophysiology of ventricular pre-excitation. J. Am. Coll. Cardiol. 48, 2045–2052. doi: 10.1016/j.jacc.2006.08.019

PubMed Abstract | CrossRef Full Text | Google Scholar

Berger, T., Pfeifer, B., Hanser, F. F., Hintringer, F., Fischer, G., Netzer, M., et al. (2011). Single-beat noninvasive imaging of ventricular endocardial and epicardial activation in patients undergoing crt. PLoS ONE 6:e16255. doi: 10.1371/journal.pone.0016255

PubMed Abstract | CrossRef Full Text | Google Scholar

Bokeriya, L., Revishvili, A., Kalinin, A., Kalinin, V., Lyadzhina, O., and Fetisova, E. (2008). Hardware–software system for noninvasive electrocardiographic heart examination based on inverse problem of electrocardiography. Biomed. Eng. 42, 273–279. doi: 10.1007/s10527-009-9065-0

CrossRef Full Text | Google Scholar

Bourgault, Y., Coudiere, Y., and Pierre, C. (2009). Existence and uniqueness of the solution for the bidomain model used in cardiac electrophysiology. Nonlin. Anal. 10, 458–482. doi: 10.1016/j.nonrwa.2007.10.007

CrossRef Full Text | Google Scholar

Cluitmans, M. J. M., Bonizzi, P., Karel, J. M. H., Das, M., Kietselaer, B. L. J. H., de Jong, M. M. J., et al. (2017). In vivo validation of electrocardiographic imaging. JACC 3, 232–242. doi: 10.1016/j.jacep.2016.11.012

PubMed Abstract | CrossRef Full Text | Google Scholar

Cuculich, P. S., Schill, M. R., Kashani, R., Mutic, S., Lang, A., Cooper, D., et al. (2017). Noninvasive cardiac radiation for ablation of ventricular tachycardia. N. Engl. J. Med. 377, 2325–2336. doi: 10.1056/NEJMoa1613773

PubMed Abstract | CrossRef Full Text | Google Scholar

Davey, K., and Hinduja, S. (1989). Analytical integration of linear three-dimensional triangular elements in bem. Appl. Math. Model. 13, 450–461. doi: 10.1016/0307-904X(89)90093-0

CrossRef Full Text | Google Scholar

Dubois, R., Shah, A. J., Hocini, M., Denis, A., Derval, N., Cochet, H., et al. (2015). Non-invasive cardiac mapping in clinical practice: application to the ablation of cardiac arrhythmias. J. Electrocardiol. 48, 966–974. doi: 10.1016/j.jelectrocard.2015.08.028

PubMed Abstract | CrossRef Full Text | Google Scholar

Dunavant, D. (1985). High degree efficient symmetrical gaussian quadrature rules for the triangle. Int. J. Numer. Methods Eng. 21, 1129–1148. doi: 10.1002/nme.1620210612

CrossRef Full Text | Google Scholar

Erem, B., Coll-Font, J., Orellana, R. M., Stovícek, P., and Brooks, D. H. (2014). Using transmural regularization and dynamic modeling for noninvasive cardiac potential imaging of endocardial pacing with imprecise thoracic geometry. IEEE Trans. Med. imaging 33, 726–738. doi: 10.1109/TMI.2013.2295220

PubMed Abstract | CrossRef Full Text | Google Scholar

Erkapic, D., Greiss, H., Pajitnev, D., Zaltsberg, S., Deubner, N., Berkowitsch, A., et al. (2014). Clinical impact of a novel three-dimensional electrocardiographic imaging for non-invasive mapping of ventricular arrhythmias—a prospective randomized trial. EP Europace 17, 591–597. doi: 10.1093/europace/euu282

PubMed Abstract | CrossRef Full Text | Google Scholar

Franzone, P. C., Taccardi, B., and Viganotti, C. (1978). “An approach to inverse calculation of epicardial potentials from body surface maps,” in Electrocardiology III/Vectorcardiography. 3rd International Congress / 17th International Symposium on Vectorcardiography Vol. 21. ed F. Kornreich (Brussels: Karger Publishers), 50–54. doi: 10.1159/000400421

CrossRef Full Text | Google Scholar

Geselowitz, D. B. (1989). On the theory of the electrocardiogram. Proc. IEEE 77, 857–876. doi: 10.1109/5.29327

CrossRef Full Text | Google Scholar

Geuzaine, C., and Remacle, J.-F. (2009). Gmsh: a 3-d finite element mesh generator with built-in pre-and post-processing facilities. Int. J. Numer. Methods Eng. 79, 1309–1331. doi: 10.1002/nme.2579

CrossRef Full Text | Google Scholar

Guillem, M. S., Climent, A. M., Millet, J., Arenal, Á., Fernández-Avilés, F., Jalife, J., et al. (2013). Non-invasive localization of maximal frequency sites of atrial fibrillation by body surface potential mapping. Circ. Arrhythm. Electrophysiol. 6, 294–301. doi: 10.1161/CIRCEP.112.000167

CrossRef Full Text | Google Scholar

Gulrajani, R. M. (1998). The forward and inverse problems of electrocardiography. IEEE Eng. Med. Biol. Mag. 17, 84–122. doi: 10.1109/51.715491

PubMed Abstract | CrossRef Full Text | Google Scholar

Hansen, P. C. (2000). “The l-curve and its use in the numerical treatment of inverse problems,” in in Computational Inverse Problems in Electrocardiology, Advances in Computational Bioengineering, ed P. Johnston (Southampton: WIT Press), 119–142.

Google Scholar

He, B., Li, G., and Zhang, X. (2003). Noninvasive imaging of cardiac transmembrane potentials within three-dimensional myocardium by means of a realistic geometry anisotropic heart model. IEEE Trans. Biomed. Eng. 50, 1190–1202. doi: 10.1109/TBME.2003.817637

PubMed Abstract | CrossRef Full Text | Google Scholar

Horácek, B. M., and Clements, J. C. (1997). The inverse problem of electrocardiography: a solution in terms of single- and double-layer sources on the epicardial surface. Math. Biosci. 144, 119–154. doi: 10.1016/S0025-5564(97)00024-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Huiskamp, G. (1991). Difference formulas for the surface laplacian on a triangulated surface. J. Comput. Phys. 95, 477–496. doi: 10.1016/0021-9991(91)90286-T

CrossRef Full Text | Google Scholar

Kalinin, A. (2011). Iterative algorithm for the inverse problem of electrocardiography in a medium with piecewise-constant electrical conductivity. Comput. Math. Model. 22, 30–34. doi: 10.1007/s10598-011-9086-x

CrossRef Full Text | Google Scholar

Kalinin, V., Kalinin, A., Schulze, W., Potyagaylo, D., and Shlapunov, A. (2017). On the correctness of the transmembrane potential based inverse problem of ECG. Comput. Cardiol. 44, 1–4. doi: 10.22489/CinC.2017.077-438

CrossRef Full Text | Google Scholar

Keller, D. U., Weiss, D. L., Dossel, O., and Seemann, G. (2012). Influence of I_Ks heterogeneities on the genesis of the t-wave: a computational evaluation. IEEE Trans. Biomed. Eng. 59, 311–322. doi: 10.1109/TBME.2011.2168397

CrossRef Full Text | Google Scholar

Kubo, M. (1994). L2-conditional stability estimate for the cauchy problem for the laplace equation. J. Inverse Ill-Posed Prob. 2, 253–262. doi: 10.1515/jiip.1994.2.3.253

CrossRef Full Text | Google Scholar

Metzner, A., Wissner, E., Tsyganov, A., Kalinin, V., Schlüter, M., Lemes, C., et al. (2017). Noninvasive phase mapping of persistent atrial fibrillation in humans: comparison with invasive catheter mapping. Ann. Noninvasive Electrocardiol. 23:e12527. doi: 10.1111/anec.12527

PubMed Abstract | CrossRef Full Text | Google Scholar

Mirams, G. R., Arthurs, C. J., Bernabeu, M. O., Bordas, R., Cooper, J., Corrias, A., et al. (2013). Chaste: an open source c++ library for computational physiology and biology. PLoS Comput. Biol. 9:e1002970. doi: 10.1371/journal.pcbi.1002970

PubMed Abstract | CrossRef Full Text | Google Scholar

Ramanathan, C., Ghanem, R. N., Jia, P., Ryu, K., and Rudy, Y. (2004). Noninvasive electrocardiographic imaging for cardiac electrophysiology and arrhythmia. Nat. Med. 10, 422–428. doi: 10.1038/nm1011

PubMed Abstract | CrossRef Full Text | Google Scholar

Revishvili, A. S., Wissner, E., Lebedev, D. S., Lemes, C., Deiss, S., Metzner, A., et al. (2015). Validation of the mapping accuracy of a novel non-invasive epicardial and endocardial electrophysiology system. EP Europace 17, 1282–1288. doi: 10.1093/europace/euu339

PubMed Abstract | CrossRef Full Text | Google Scholar

Rodrigo, M., Climent, A. M., Liberos, A., Fernández-Avilés, F., Berenfeld, O., Atienza, F., et al. (2017). Technical considerations on phase mapping for identification of atrial reentrant activity in direct-and inverse-computed electrograms. Circ. Arrhythm. Electrophysiol. 10:e005008. doi: 10.1161/CIRCEP.117.005008

PubMed Abstract | CrossRef Full Text | Google Scholar

Rudic, B., Chaykovskaya, M., Tsyganov, A., Kalinin, V., Tülümen, E., Papavassiliu, T., et al. (2016). Simultaneous non-invasive epicardial and endocardial mapping in patients with brugada syndrome: new insights into arrhythmia mechanisms. J. Am. Heart Assoc. 5:e004095. doi: 10.1161/JAHA.116.004095

PubMed Abstract | CrossRef Full Text | Google Scholar

Rudy, Y., and Messinger-Rapport, B. J. (1988). The inverse problem in electrocardiography: solutions in terms of epicardial potentials. Crit. Rev. Biomed. Eng. 16, 215–268.

PubMed Abstract | Google Scholar

Schulze, W. H., Henar, F. E., Potyagaylo, D., Loewe, A., Stenroos, M., and Dössel, O. (2013). “Kalman filter with augmented measurement model: an ECG imaging simulation study,” in Functional Imaging and Modeling of the Heart. FIMH 2013. Lecture Notes in Computer Science Vol. 7945, eds S. Ourselin, D. Rueckert, and N. Smith (Berlin; Heidelberg: Springer), 200–207. doi: 10.1007/978-3-642-38899-6_24

CrossRef Full Text | Google Scholar

Shah, A. J., Lim, H. S., Yamashita, S., Zellerhoff, S., Berte, B., Mahida, S., et al. (2014). Non invasive ECG mapping to guide catheter ablation. J. Atrial Fibrill. 7:1139. doi: 10.4022/jafib.1139

PubMed Abstract | CrossRef Full Text | Google Scholar

Skipa, O. (2004). Linear Inverse Problem of Electrocardiography: Epicardial Potentials and Transmembrane Voltages. Karlsruhe: Helmesverl.

Google Scholar

Takeuchi, T., and Yamamoto, M. (2008). Tikhonov regularization by a reproducing kernel hilbert space for the cauchy problem for an elliptic equation. SIAM J. Sci. Comput. 31, 112–142. doi: 10.1137/070684793

CrossRef Full Text | Google Scholar

ten Tusscher, K., Noble, D., Noble, P. J., and Panfilov, A. V. (2004). A model for human ventricular tissue. Am. J. Physiol. Heart Circul. Physiol. 286, H1573–H1589. doi: 10.1152/ajpheart.00794.2003

PubMed Abstract | CrossRef Full Text | Google Scholar

ten Tusscher, K. H., and Panfilov, A. V. (2006). Alternans and spiral breakup in a human ventricular tissue model. Am. J. Physiol. Heart Circul. Physiol. 291, H1088–H1100. doi: 10.1152/ajpheart.00109.2006

PubMed Abstract | CrossRef Full Text | Google Scholar

Tsyganov, A., Wissner, E., Metzner, A., Mironovich, S., Chaykovskaya, M., Kalinin, V., et al. (2017). Mapping of ventricular arrhythmias using a novel noninvasive epicardial and endocardial electrophysiology system. J. Electrocardiol. 51, 92–98. doi: 10.1016/j.jelectrocard.2017.07.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Tung, L. (1978). A Bi-domain Model for Describing Ischemic Myocardial d-c Potentials. PhD thesis, MIT, Cambridge, MA.

Ushenin, K. S., Dokuchaev, A., Magomedova, S., Sopov, O. V., Kalinin, V. V., and Solovyova, O. (2017). “Role of myocardial properties and pacing lead location on ECG in personalized paced heart models,” in Computing in Cardiology Conference (CinC) (Rennes: IEEE), 56. Available online at: https://ieeexplore.ieee.org/document/8331631

Google Scholar

van Dam, P. M., Oostendorp, T. F., Linnenbank, A. C., and Van Oosterom, A. (2009). Non-invasive imaging of cardiac activation and recovery. Ann. Biomed. Eng. 37, 1739–1756. doi: 10.1007/s10439-009-9747-5

PubMed Abstract | CrossRef Full Text | Google Scholar

van Oosterom, A. (2012). Closed-form analytical expressions for the potential fields generated by triangular monolayers with linearly distributed source strength. Med. Biol. Eng. Comput. 50, 1–9. doi: 10.1007/s11517-011-0837-9

PubMed Abstract | CrossRef Full Text | Google Scholar

van Oosterom, A. (2014). A comparison of electrocardiographic imaging based on two source types. Europace 16, iv120–iv128. doi: 10.1093/europace/euu268

PubMed Abstract | CrossRef Full Text | Google Scholar

Varma, N. (2015). Left ventricular electrical activation during right ventricular pacing in heart failure patients with lbbb: visualization by electrocardiographic imaging and implications for cardiac resynchronization therapy. J. Electrocardiol. 48, 53–61. doi: 10.1016/j.jelectrocard.2014.09.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, L., Dawoud, F., Yeung, S. K., Shi, P., Wong, K. C., Liu, H., et al. (2013). Transmural imaging of ventricular action potentials and post-infarction scars in swine hearts. IEEE Trans. Med. Imaging 32, 731–747. doi: 10.1109/TMI.2012.2236567

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, L., Gharbia, O. A., Horáček, B. M., and Sapp, J. L. (2016). Noninvasive epicardial and endocardial electrocardiographic imaging of scar-related ventricular tachycardia. J. Electrocardiol. 49, 887–893. doi: 10.1016/j.jelectrocard.2016.07.026

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Y., and Rudy, Y. (2006). Application of the method of fundamental solutions to potential-based inverse electrocardiography. Ann. Biomed. Eng. 34, 1272–1288. doi: 10.1007/s10439-006-9131-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Wissner, E., Kalinin, V., Kalinin, A., Sopov, O., Chmelevsky, M., Metzner, A., et al. (2018). Noninvasive phase mapping of atrial flutter in humans – comparison with invasive mapping. Eur. J. Arrhythm. Electrophysiol. 4, 15–20. doi: 10.17925/EJAE.2018.04.01.15

CrossRef Full Text | Google Scholar

Wissner, E., Revishvili, A., Metzner, A., Tsyganov, A., Kalinin, V., Lemes, C., et al. (2016). Noninvasive epicardial and endocardial mapping of premature ventricular contractions. Europace 19, 843–849. doi: 10.1093/europace/euw103

PubMed Abstract | CrossRef Full Text | Google Scholar

Xu, J., Dehaghani, A. R., Gao, F., and Wang, L. (2014). Noninvasive transmural electrophysiological imaging based on minimization of total-variation functional. IEEE Trans. Med. Imaging 33, 1860–1874. doi: 10.1109/TMI.2014.2324900

PubMed Abstract | CrossRef Full Text | Google Scholar

Yun, Z. Q., Tan, B. D., and Huang, J. (1997). A comparison of different transfer matrices for the calculation of an inverse problem using bem. IEEE Trans. Magn. 33, 1189–1191. doi: 10.1109/20.582465

CrossRef Full Text | Google Scholar

Zhou, Z., Jin, Q., Yu, L., Wu, L., and He, B. (2016). Noninvasive imaging of human atrial activation during atrial flutter and normal rhythm from body surface potential maps. PLoS ONE 11:e0163445. doi: 10.1371/journal.pone.0163445

PubMed Abstract | CrossRef Full Text | Google Scholar