Unsupervised stochastic learning and reduced order modeling for global sensitivity analysis in cardiac electrophysiology models

BACKGROUND AND OBJECTIVE

Numerical simulations in electrocardiology are often affected by various uncertainties inherited from the lack of precise knowledge regarding input values including those related to the cardiac cell model, domain geometry, and boundary or initial conditions used in the mathematical modeling. Conventional techniques for uncertainty quantification in modeling electrical activities of the heart encounter significant challenges, primarily due to the high computational costs associated with fine temporal and spatial scales. Additionally, the need for numerous model evaluations to quantify ubiquitous uncertainties increases the computational challenges even further.

METHODS

In the present study, we propose a non-intrusive surrogate model to perform uncertainty quantification and global sensitivity analysis in cardiac electrophysiology models. The proposed method combines an unsupervised machine learning technique with the polynomial chaos expansion to reconstruct a surrogate model for the propagation and quantification of uncertainties in the electrical activity of the heart. The proposed methodology not only accurately quantifies uncertainties at a very low computational cost but more importantly, it captures the targeted quantity of interest as either the whole spatial field or the whole temporal period. In order to perform sensitivity analysis, aggregated Sobol indices are estimated directly from the spectral mode of the polynomial chaos expansion.

RESULTS

We conduct Uncertainty Quantification (UQ) and global Sensitivity Analysis (SA) considering both spatial and temporal variations, rather than limiting the analysis to specific Quantities of Interest (QoIs). To assess the comprehensive performance of our methodology in simulating cardiac electrical activity, we utilize the monodomain model. Additionally, sensitivity analysis is performed on the parameters of the Mitchell-Schaeffer cell model.

CONCLUSIONS

Unlike conventional techniques for uncertainty quantification in modeling electrical activities, the proposed methodology performs at a low computational cost the sensitivity analysis on the cardiac electrical activity parameters. The results are fully reproducible and easily accessible, while the proposed reduced-order model represents a significant contribution to enhancing global sensitivity analysis in cardiac electrophysiology.

Balancing conduction velocity error in cardiac electrophysiology using a modified quadrature approach

Conduction velocity error is often the main culprit behind the need for very fine spatial discretizations and high computational effort in cardiac electrophysiology problems. In light of this, a novel approach for simulating an accurate conduction velocity in coarse meshes with linear elements is suggested based on a modified quadrature approach. In this approach, the quadrature points are placed at arbitrary offsets of the isoparametric coordinates. A numerical study illustrates the dependence of the conduction velocity on the spatial discretization and the conductivity when using different quadrature rules and calculation approaches. Additionally, examples using the modified quadrature in coarse meshes for wave propagation demonstrate the improved accuracy of the conduction velocity with this method. This novel approach possesses great potential in reducing the computational effort required but remains limited to specific linear elements and experiences a reduction in accuracy for irregular meshes and heterogeneous conductivities. Further research can focus on developing an adaptive quadrature and extending the approach to other element formulations in order to make the approach more generally applicable.

A numerical study on the effects of spatial and temporal discretization in cardiac electrophysiology

Millions of degrees of freedom are often required to accurately represent the electrophysiology of the myocardium due to the presence of discretization effects. This study seeks to explore the influence of temporal and spatial discretization on the simulation of cardiac electrophysiology in conjunction with changes in modeling choices. Several finite element analyses are performed to examine how discretization affects solution time, conduction velocity and electrical excitation. Discretization effects are considered along with changes in the electrophysiology model and solution approach. Two action potential models are considered: the Aliev-Panfilov model and the ten Tusscher-Noble-Noble-Panfilov model. The solution approaches consist of two time integration schemes and different treatments for solving the local system of ordinary differential equations. The efficiency and stability of the calculation approaches are demonstrated to be dependent on the action potential model. The dependency of the conduction velocity on the element size and time step is shown to be different for changes in material parameters. Finally, the discrepancies between the wave propagation in coarse and fine meshes are analyzed based on the temporal evolution of the transmembrane potential at a node and its neighboring Gauss points. Insight obtained from this study can be used to suggest new methods to improve the efficiency of simulations in cardiac electrophysiology.

Stabilized hybrid discontinuous Galerkin finite element method for the cardiac monodomain equation

Numerical methods for solving the cardiac electrophysiology model, which describes the electrical activity in the heart, are proposed. The model problem consists of a nonlinear reaction-diffusion partial differential equation coupled to systems of ordinary differential equations that describes electrochemical reactions in cardiac cells. The proposed methods combine an operator splitting technique for the reaction-diffusion equation with primal hybrid methods for spatial discretization considering continuous or discontinuous approximations for the Lagrange multiplier. A static condensation is adopted to form a reduced global system in terms of the multiplier only. Convergence studies exhibit optimal rates of convergence and numerical experiments show that the proposed schemes can be more efficient than standard numerical techniques commonly used in this context when preconditioned iterative methods are used for the solution of linear systems.

Modeling and simulation of hypothermia effects on cardiac electrical dynamics

Previous experimental evidence has shown the effect of temperature on the action potential duration (APD). It has also been demonstrated that regional cooling of the heart can prolong the APD and promote the termination of ventricular tachycardia. The aim of this study is to demonstrate the effect of hypothermia in suppressing cardiac arrhythmias using numerical modeling. For this purpose, we developed a mathematical model that couples Pennes’ bioheat equation and the bidomain model to simulate the effect of heat on the cardiac action potential. The simplification of the proposed heat–bidomain model to the heat–monodomain model is provided. A suitable numerical scheme for this coupling, based on a time adaptive mesh finite element method, is also presented. First, we performed two-dimensional numerical simulations to study the effect of heat on a regular electrophysiological wave, with the comparison of the calculated and experimental values of Q₁₀. Then, we demonstrated the effect of global hypothermia in suppressing single and multiple spiral waves.

An adaptive hybridizable discontinuous Galerkin approach for cardiac electrophysiology

Cardiac electrophysiology simulations are numerically challenging because of the propagation of a steep electrochemical wave front and thus require discretizations with small mesh sizes to obtain accurate results. In this work, we present an approach based on the hybridizable discontinuous Galerkin method (HDG), which allows an efficient implementation of high-order discretizations into a computational framework. In particular, using the advantage of the discontinuous function space, we present an efficient p-adaptive strategy for accurately tracking the wave front. The HDG allows to reduce the overall degrees of freedom in the final linear system to those only on the element interfaces. Additionally, we propose a rule for a suitable integration accuracy for the ionic current term depending on the polynomial order and the cell model to handle high-order polynomials. Our results show that for the same number of degrees of freedom, coarse high-order elements provide more accurate results than fine low-order elements. Introducing p-adaptivity further reduces computational costs while maintaining accuracy by restricting the use of high-order elements to resolve the wave front. For a patient-specific simulation of a cardiac cycle, p-adaptivity reduces the average number of degrees of freedom by 95% compared to the nonadaptive model. In addition to reducing computational costs, using coarse meshes with our p-adaptive high-order HDG method also simplifies practical aspects of mesh generation and postprocessing.

Performance evaluation of GPU parallelization, space-time adaptive algorithms, and their combination for simulating cardiac electrophysiology

The use of computer models as a tool for the study and understanding of the complex phenomena of cardiac electrophysiology has attained increased importance nowadays. At the same time, the increased complexity of the biophysical processes translates into complex computational and mathematical models. To speed up cardiac simulations and to allow more precise and realistic uses, 2 different techniques have been traditionally exploited: parallel computing and sophisticated numerical methods. In this work, we combine a modern parallel computing technique based on multicore and graphics processing units (GPUs) and a sophisticated numerical method based on a new space-time adaptive algorithm. We evaluate each technique alone and in different combinations: multicore and GPU, multicore and GPU and space adaptivity, multicore and GPU and space adaptivity and time adaptivity. All the techniques and combinations were evaluated under different scenarios: 3D simulations on slabs, 3D simulations on a ventricular mouse mesh, ie, complex geometry, sinus-rhythm, and arrhythmic conditions. Our results suggest that multicore and GPU accelerate the simulations by an approximate factor of 33×, whereas the speedups attained by the space-time adaptive algorithms were approximately 48. Nevertheless, by combining all the techniques, we obtained speedups that ranged between 165 and 498. The tested methods were able to reduce the execution time of a simulation by more than 498× for a complex cellular model in a slab geometry and by 165× in a realistic heart geometry simulating spiral waves. The proposed methods will allow faster and more realistic simulations in a feasible time with no significant loss of accuracy.

Simulation of cardiac electrophysiology on next-generation high-performance computers

Models of cardiac electrophysiology consist of a system of partial differential equations (PDEs) coupled with a system of ordinary differential equations representing cell membrane dynamics. Current software to solve such models does not provide the required computational speed for practical applications. One reason for this is that little use is made of recent developments in adaptive numerical algorithms for solving systems of PDEs. Studies have suggested that a speedup of up to two orders of magnitude is possible by using adaptive methods. The challenge lies in the efficient implementation of adaptive algorithms on massively parallel computers. The finite-element (FE) method is often used in heart simulators as it can encapsulate the complex geometry and small-scale details of the human heart. An alternative is the spectral element (SE) method, a high-order technique that provides the flexibility and accuracy of FE, but with a reduced number of degrees of freedom. The feasibility of implementing a parallel SE algorithm based on fully unstructured all-hexahedra meshes is discussed. A major computational task is solution of the large algebraic system resulting from FE or SE discretization. Choice of linear solver and preconditioner has a substantial effect on efficiency. A fully parallel implementation based on dynamic partitioning that accounts for load balance, communication and data movement costs is required. Each of these methods must be implemented on next-generation supercomputers in order to realize the necessary speedup. The problems that this may cause, and some of the techniques that are beginning to be developed to overcome these issues, are described.