Parametrization of activation based cardiac electrophysiology models using bidomain model simulations

Birgit Stender

doi:10.1515/cdbme-2016-0135

Article Open Access

Parametrization of activation based cardiac electrophysiology models using bidomain model simulations

Birgit Stender

Published/Copyright: September 30, 2016

Published by

Become an author with De Gruyter Brill

Author Information Explore this Subject

From the journal Current Directions in Biomedical Engineering Volume 2 Issue 1

Abstract

Eikonal models are useful to compute approximate solutions of cardiac excitation propagation in a computationally efficient way. In this work the underlying conduction velocities for different cell types were computed solving the classical bidomain model equations for planar wavefront propagation. It was further investigated how changes in the conductivity tensors within the bidomain model analytically correspond to changes in the conduction velocity. The error in the presence of local front curvature for the derived eikonal model parametrization were analyzed. The conduction velocity simulated based on the bidomain model was overestimated by a maximum of 10%.

Keywords: bidomain model; cardiac electrophysiology; conduction velocity; eikonal model; monodomain model

1 Introduction

One treatment option in case of abnormal cardiac excitation patterns is ablation therapy. It is commonly used in case of paroxysmal or persistent atrial fibrillation. Ablation therapy is further applied in presence of premature ventricular contractions arising from focal activity in the transition zone between an infarction scar and healthy myocardial tissue.

In clinical practice the intervention is currently not planned but conducted in accordance with the effective treatment guidelines. Intervention planning in contrast is typically an iterative process: It is required to predict the result of a certain treatment suggestion which is then varied based on the computed outcome. Therefore the simulation run needs to rely on computationally efficient forward model formulations. The bidomain model and its simplification, the monodomain model, are well established in computational cardiology. They are computationally demanding due to the fine discretization in space and time [9], [10]. To overcome the resulting requirements concerning memory and CPU usage different approximations by so-called eikonal models have been derived based on the classical monodomain model equations [4], [6].

2 Methods

Eikonal models describe the myocardial state transition in terms of the depolarization wavefront arrival time. The wavefront contour Γ(t) is formed by all positions x having a first arrival time τ equal to t:

(1)Γ⁢(t):={x|τ⁢(x)=t}

Neglecting the curvature of the front the simplest eikonal model formulation is

(2)〈∇ τ,ACV∗⁢CV∗⁢AT⁢∇ τ〉=1,

where A denotes the local fiber coordinate system. CV^* is the local conduction velocity tensor. It contains the scalar planar conduction velocities CV^l, CV^t and CVⁿ for the respective propagation in longitudinal, tangential and normal direction of the local fiber coordinate system. All variables given in local coordinates will be marked by an asterisk in the following.

While the distribution of local conduction velocities just represents the model parameters in terms of an eikonal model formulation it depends on many parameters within the bidomain model equations. Among these are the entries of the intra- and extracellular conductivity tensors M_i^* and M_e^*, parameters of the single cell model, the stimulation current density and the basic cycle length of the stimulation protocol as illustrated in Figure 1. The conduction velocity further depends on the local front curvature.

Figure 1

Left: Atrial action potential waveforms for different basic cycle lengths. Right: Change caused by diffusive coupling within the bidomain model in comparison to the uncoupled single cell response.

A known result from eikonal model derivations is

(3)CV∝n→⋅rm⁢M⁢n→,

where M is the conductivity tensor of the monodomain model. CV denotes the conduction velocity along the depolarization front normal n→. r_m is defined by the passive conductance of the cell membrane per unit volume at resting membrane potential v_{m, rest}. Assuming just coincident basis vectors for the intracellular and extracellular space within the bidomain model formulation M_e^* can be expressed as

(4)Me∗=Λ⁢Mi∗ ⇔

(5)=(σelσil000σetσit000σenσin)⁢Mi∗.

Using the notation

(6)λl:=σelσilλt:=σetσitλn:=σenσin

and neglecting the stimulation currents the monodomain model equation can then be restated as

(7)∇⋅(M⁢∇ vm)=χ⁢(cm⁢∂ vm∂ t+iion)

where χ is the surface-to-volume ratio. The local conductivity tensor M^* = A^TMA is given by

(8)M∗=(I+(I+Λ)−1)⁢Mi∗

(9)=(λl1+λl000λt1+λt000λn1+λn)⁢Mi∗.

This formulation was introduced as augmented monodomain model by Bishop and Plank [2], [3].

The trace of ventricular action potentials is known to vary in dependence on the heart chamber and the location of the cell in transmural and apico-basal direction [1], [12], [13]. These effects are called repolarization gradients or dispersion of action potential duration. The question arise, if repolarization gradients influence as well the conduction velocity distribution.

The study was performed in the following steps: Initially the entries of the conduction velocity tensor CV^* were derived based on planar wavefront propagation. The 1D bidomain model simulations comprised atrial and ventricular single cell models to facilitate the parametrization of whole heart models. In a second step the sensitivity of CV^* in response to changes in the intra- and extracellular conductivity was investigated. In a third step the error

(10)e=n→⋅cv⁢n→−CVexactCVexact

relying on CV^* in the presence of nonnegligible front curvature was determined. CV_exact denotes the conduction velocity of the bidomain model in direction of the respective front normal. The bidomain model equations were computed with Galerkin finite element methods (FEM) and operator splitting techniques as described by Sundness et al. [11]. Two different human cell models of ventricular myocardium were taken into account. The first ventricular cell model integrated is a model with minimum complexity presented by Bueno-Orovio et al. in 2008 [1]. It will be referred as BCF model in the following. The second ventricular single cell model is the ten Tusscher-Noble-Noble-Panfilov (TNNP) model [12], [13]. For action potential propagation in atrial myocardium the Courtemanche-Ramirez-Nattel (CRN) model [5] was applied.

Assumptions on the electrical parameters differ among published simulation studies. Parameter values for the capacitance per unit cell membrane and the amount of cell area per unit volume of ventricular myocardium were chosen according to Li et al. [7]. The entries of the conductivity tensors applied by Sundnes et al. were scaled to achieve an unchanged diffusion constant with respect to Bueno-Orovio et al. The resulting values for the entries of the conductivity tensors in longitudinal direction σil=7.12⁢ mS/cm and σel=4.74⁢ mS/cm were well in accordance with the values applied by Weiss et al. [15]. In the last-mentioned work the anisotropic ratios between longitudinal and tangential directions were considerably larger than those taken into account by Sundnes et al. The resulting entries of M^*_e and M^*_i are far from being linearly dependent.

Atrial myocardium was modeled as isotropic. Parameter values were set according to Mainardi et al. [8].

Spatial discretization for all simulations was Δ⁢x=25⁢μ⁢m in accordance with [10]. A discretization in time Δ⁢t=1⁢μ⁢s was chosen based on an evaluation of the single cell model convergence.

3 Results

3.1 Planar conduction velocities

For planar propagation the conduction velocities are assembled in Table 1. In contrast to the BCF model the differences in conduction velocity among the cell types in transmural direction were marginal for the TNNP model. The relative deviations with respect to the subendocardial cell accounted for less than 0.04%. The bidomain model simulations for the planar wave propagation in the local fiber directions were repeated with a scaling factor in the respective direction varied in the range [0.2, 1.8]. The scaling was applied to the entries of both M_i^* and M_e^*, M_i^* only, and M_e^* only. The relative changes in the conduction velocities with respect to the ones for the reference conditions were determined at steady state. These relative changes will be denoted in the following as CVI+ECV0 for the first case, CVICV0 for the second case and CVECV0 for the third case. The results for the last two cases are displayed in Figure 2. The relative changes were found to be independent from the choice of the single cell model and its parameters. This can be explained based on equation (3). Scaling the entries of M_i and M_e with the same scaling factor α_k, k∈{l,t,n} results in the corresponding scaling

Table 1

Planar conduction velocities in local fiber directions computed with different single cell models. Values are given in (mm/s).

	CV^l	CV^t	CVⁿ
BCF
subendo	612.0	440.9	282.6
m-cell	705.0	507.9	325.5
subepi	570.8	411.2	263.5
BCF
subendo	612.0	440.9	282.6
m-cell	705.0	507.9	325.5
subepi	570.8	411.2	263.5
TNNP
RV
subendo	555.5	403.5	262.9
m-cell	555.4	403.5	262.9
subepi	555.4	403.5	262.9
TNNP
LV
subendo	555.5	403.5	262.9
m-cell	554.1	402.5	262.3
subepi	554.3	402.7	262.4
CRN	742

(11)CVI+EkCV0k=αl

of the conduction velocity components with respect to the ones computed for reference parameters. The respective results for changing the entries of M_i or M_e are

(12)CVIkCV0k=αk⁢(1+λk)αk+λk

(13)CVEkCV0k=αk⁢(1+λk)1+αk⁢λk

Figure 2

Sensivity of the relative planar conduction velocity towards changes in the entries of the local intra- and extracellular conductivity tensors.

3.2 Ellipsoidal wavefront propagation

In the presence of front curvature the local conduction velocity is decreased in comparison to the planar case as illustrated in Figure 3 on the left based on simple spatial configurations. The conduction velocities were derived integrating the BCF cell model with the endocardial parameter set. The configuration with maximum deviation for the simplest assumption of just constant planar conduction velocity is displayed in Figure 3 on the right hand-side. The maximum relative error in the overshoot accounted for 10%. In a distance greater than 5 mm from the center of the ellipsoid a relative error of 5% was not exceeded.

Figure 3

Left: Comparison of conduction velocities computed for planar propagation, radial propagation as cylindrical wave and radial propagation as spherical depolarization wave. Right: Maximum relative errors introduced by the use of different eikonal model approximations.

4 Discussion

The conduction velocity is an important quantity in cardiac electrophysiology. I derived a parametrization of the conduction velocity suitable for whole human hearts. The results are based on the bidomain model in case of planar wavefront propagation. The parametrization was found to be scalable in a manner which does not depend on the particular single cell model. Relying on different intra- or extracellular conductivity tensors the results can be adapted without the need for repeated bidomain model computations. In addition the impact of the front curvature was investigated which can be high in the presence of focal activity. In case of only moderate spatial discretization bidomain model simulations have been reported to show a relative error comparable in magnitude [10]. Pezzuto et al. reported a relative error in conduction velocity exceeding 10% within the whole domain at a spatial discretization coarser than 0.3 mm to 0.7 mm, depending on the mass lumping scheme. Using finite elements the conduction velocity is commonly overestimated while finite difference schemes tend to underestimate the correct value.

5 Conclusion

The presented findings support the admissibility and transferability of eikonal model simulations in comparison to the bidomain model in case of moderate front curvature. This requirement is commonly fulfilled in healthy and ischemic ventricular myocardium. In case of atrial fibrillation showing spiral wavefront propagation patterns bidomain model equations seem to be more advantageous. The inconsistent results on the influence of dispersion might be further investigated in terms of the single cell models.

Acknowledgement

I would like to thank Gunnar Seemann for pointing out the work of Bishop and Plank.

Author’s Statement

Research funding: The author state no funding involved. Conflict of interest: Authors state no conflict of interest. Material and Methods: Informed consent: Informed consent is not applicable. Ethical approval: The conducted research is not related to either human or animal use.

References

[1] Bueno-Orovio A, Cherry EM, Fenton FH. Minimal model for human ventricular action potentials in tissue. J Theor Biol. 2008;253:544–60.10.1016/j.jtbi.2008.03.029Search in Google Scholar PubMed

[2] Bishop MJ, Plank G. Bidomain ECG simulations using an augmented monodomain model for the cardiac source. IEEE Trans Biomed Eng. 2011;58.10.1109/TBME.2011.2148718Search in Google Scholar PubMed PubMed Central

[3] Bishop MJ, Plank G. Representing cardiac bidomain bath-loading effects by an augmented monodomain approach: application to complex ventricular models. IEEE Trans Biomed Eng. 2011;58:1066–75.10.1109/TBME.2010.2096425Search in Google Scholar PubMed PubMed Central

[4] Colli Franzone P, Guerri L, Rovida S. Wavefront propagation in an activation model of the anisotropic cardiac tissue: asymptotic analysis and numerical simulations. J Math Biol. 1990;28:121–76.10.1007/BF00163143Search in Google Scholar PubMed

[5] Courtemanche M, Ramirez RJ, Nattel S. Ionic mechanisms underlying human atrial action potential properties: insights from a mathematical model. Am J Physiol. 1998;275:H301–21.10.1152/ajpheart.1998.275.1.H301Search in Google Scholar PubMed

[6] Keener JP. An eikonal-curvature equation for action potential propagation in myocardium. J Math Bio. 1991;29:629–51.10.1007/BF00163916Search in Google Scholar PubMed

[7] Li GR, Yang B, Feng J, Bosch RF, Carrier M, Nattel S. Transmembrane ica contributes to rate-dependent changes of action potentials in human ventricular myocytes. Am J Physiol. 1999;276:98–106.10.1152/ajpheart.1999.276.1.H98Search in Google Scholar

[8] Mainardi L, Sörnmo S, Cerutti S. Understanding atrial fibrillation: the signal processing contribution. Morgan Claypool. 2008.10.2200/S00152ED1V01Y200809BME024Search in Google Scholar

[9] Niederer SA, Kerfoot E, Benson AP, Bernabeu MO, Bernus O, Bradley C, et al. Verification of cardiac tissue electrophysiology simulators using an N-version benchmark. Philos Trans A Math Phys Eng Sci. 2011;369:4331–51.10.1098/rsta.2011.0139Search in Google Scholar PubMed PubMed Central

[10] Pezzuto S, Hake J, Sundnes J. Space-discretization error analysis and stabilization schemes for conduction velocity in cardiac electrophysiology. Int J Numer Method Biomed Eng. 2016;doi: 10.1002/cnm.2762.Search in Google Scholar PubMed

[11] Sundness J, Lines GT, Cai X, Nielsen BF, Mardal K, Tveito A. Computing the electrical activity in the heart. Berlin: Springer; 2006.Search in Google Scholar

[12] ten Tusscher KH, Noble D, Noble PJ, Panfilov AV. A model for human ventricular tissue. Am J Physiol. 2004;286:H1573–89.10.1152/ajpheart.00794.2003Search in Google Scholar PubMed

[13] ten Tusscher KH, Panfilov AV. Alternans and spiral breakup in a human ventricular tissue model. Am J Physiol Heart Circ Physiol. 2006;291:H1088–100.10.1152/ajpheart.00109.2006Search in Google Scholar PubMed

[14] Wallman M, Smith NP, Rodriguez B. A comparative study of graph-based, eikonal, and monodomain simulations for the estimation of cardiac activation times. IEEE Trans Biomed Eng. 2012;59:1739–48.10.1109/TBME.2012.2193398Search in Google Scholar PubMed

[15] Weiss DL, Ifland M, Sachse FB, Dössel O. Modeling of cardiac ischemia in human myocytes and tissue including spatiotemporal electrophysiological variations Biomed Tech. 2009;54:107–25.10.1515/BMT.2009.016Search in Google Scholar PubMed

Published Online: 2016-9-30

Published in Print: 2016-9-1

This work is licensed under the Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 License.

Articles in the same Issue

https://doi.org/10.1515/cdbme-2016-0135

Keywords for this article

bidomain model; cardiac electrophysiology; conduction velocity; eikonal model; monodomain model

Creative Commons

BY-NC-ND 4.0