EP-PINNs: Cardiac Electrophysiology Characterisation Using Physics-Informed Neural Networks

Accurately inferring underlying electrophysiological (EP) tissue properties from action potential recordings is expected to be clinically useful in the diagnosis and treatment of arrhythmias such as atrial fibrillation. It is, however, notoriously difficult to perform. We present EP-PINNs (Physics Informed Neural Networks), a novel tool for accurate action potential simulation and EP parameter estimation from sparse amounts of EP data. We demonstrate, using 1D and 2D in silico data, how EP-PINNs are able to reconstruct the spatio-temporal evolution of action potentials, whilst predicting parameters related to action potential duration (APD), excitability and diffusion coefficients. EP-PINNs are additionally able to identify heterogeneities in EP properties, making them potentially useful for the detection of fibrosis and other localised pathology linked to arrhythmias. Finally, we show EP-PINNs effectiveness on biological in vitro preparations, by characterising the effect of anti-arrhythmic drugs on APD using optical mapping data. EP-PINNs are a promising clinical tool for the characterisation and potential treatment guidance of arrhythmias.

1. Introduction

Cardiac arrhythmias are extremely common pathologies caused by disturbances in the generation or propagation of electrical signals across the heart. Atrial fibrillation (AF), the most common sustained arrhythmia, affects 0.5% of the world's population and accounts for 1% of the NHS's total budget through its large impact on patient mortality and morbidity, especially stroke (1). Catheter ablation of atrial myocardium believed to host the sources of the arrhythmia is the mainstay of AF treatment, but its long-term efficacy is disappointing (54%), especially in patients with persistent forms of the disease (43%) (2).

The mechanisms behind AF are very complex, involving the interplay of several factors at different scales, from changes in membrane proteins to alterations in cardiac tissue composition and organ shape (3). To characterise the arrhythmia, information about cardiac activity can be acquired by recording electrical potentials using electrodes placed on the chest (electrocardiogram, ECG) or, in a catheter lab, placed in direct contact with the myocardium (contact electrograms, EGMs). Expert analysis of these signals is extremely successful in the clinical diagnosis of arrhythmias and other types of cardiovascular disease (1). However, as sparse measurements of the combined electrical activity of large areas of the myocardium, electrical signals provide little direct information about local EP properties.

The ability to perform a detailed EP characterisation in the clinical setting could lead to improved treatments for arrhythmias. For example, evidence suggests that areas with abnormal EP properties (such as fibrotic or ischaemic regions) and their border-zone are often the sites of the abnormal electrical activity driving arrhythmias (3). Cardiac regions characterised by EP changes such as low conduction velocity, heightened excitability or shortened action potential duration (APD) could be prime targets for localised therapies such as catheter ablation, likely improving their efficacy. So far, ablation strategies that target these EP heterogeneities have not been successful (4), partly due to the difficulty in identifying suitable ablation sites.

In this study, we present EP-PINNs, a Physics-Informed Neural Network, as an artificial intelligence tool capable of inferring EP properties from sparse measurements of transmembrane potential, V, in cardiac tissue. We test EP-PINNs using in silico data from EP biophysical simulations in several conditions and also in vitro optical mapping data. EP-PINNs are deployed in forward mode, as high-resolution solvers of the biophysical equations that control EP systems, and also in inverse mode, as estimators of EP parameters. Tests are performed in 1D and 2D for single waves and spiral waves in homogeneous and heterogeneous conditions. We further demonstrate a pharmacological application of EP-PINNs, as a tool to characterise the effect of two different channel blockers in in vitro optical mapping data.

In the next section, we will introduce the EP biophysical model used and the PINNs technique, as applied to the EP problem. We will contextualise our work within the available techniques for parameter estimation in EP and other cardiovascular applications of PINNs.

2. Background

2.1. Biophysical Models of Cardiac Electrophysiology

Biophysical models of cardiac electrophysiology (5) are an important tool to understand how cardiac tissue properties affect the generation and propagation of cardiac electrical signals (action potentials, APs). They also offer an ideal means for the training and development of computational tools that may aim to infer EP properties from electrical and optical mapping measurements, such as EP-PINNs.

Several biophysical EP models have been proposed (see models.cellml.org/electrophysiology), each with varying degrees of detail aiming to reproduce different EP features, cardiac regions or animal/human experimental findings. Mathematically, these EP models usually take the form of a reaction-diffusion system where a diffusion term or equivalent (6) models the propagation of the electrical signal across the cardiac tissue. In the monodomain formulation, a partial differential equation (PDE) describes the spatio-temporal variations in the electrical potential across a myocyte cell membrane (V). This PDE is usually coupled to one or more ordinary differential equations (ODEs) describing how, at each point in time and space, V and other local state variables both determine and are determined by the flux of ions across the cell membrane (5).

The most parsimonious model of the action potential describes it as a travelling excitation wave followed by a non-excitable (refractory) region. This representation requires at least two state variables: V, which spreads (diffuses) across neighbouring regions, and a non-observable, non-diffusible recovery variable W which effectively controls the refractoriness and restitution properties of the model. One of the simplest models that captures these properties is the 6-parameter canine ventricular Aliev-Panfilov model (7), which models the transmembrane ionic currents (V − W relationship) using smooth, differentiable functions. Furthermore, the diffusion of V across the cardiac tissue can be described by the monodomain equation (5), which, when combined with the Aliev-Panfilov model gives:

$\begin{array}{ll}\frac{\partial V}{\partial t}=\overrightarrow{\nabla }.(D\overrightarrow{\nabla }V)-kV(V-a)(V-1)-VW& (1)\end{array}$

$\begin{array}{ll}\frac{\text{d}W}{\text{d}t}=(ϵ+\frac{{\mu }_{1}W}{V+{\mu }_2})(-W-kV(V-b-1))& (2)\end{array}$

The diffusion term $\overrightarrow{\nabla }.(D\overrightarrow{\nabla }V)$ reduces to D∇²V in the case of homogeneous and isotropic conduction, i.e., when the diffusion tensor D is approximated by the same scalar throughout. Intuitively, this term quantifies how fast V is able to spread to its immediate neighbourhood to become more spatially homogeneous. D is determined mostly by the electrical conductivity of the myocardium and is a strong determinant of the propagation velocity of the AP. To prevent a leakage of V to regions outside the heart domain, the system described by Equations (1, 2) usually obeys no flux Neumann boundary conditions: $\frac{\partial V}{\partial \overrightarrow{n}}=0$ in the boundary of the heart tissue.

The −kV(V − a)(V − 1) − VW term in Equations (1, 2) models the rapid changes in V caused by ionic fluxes across the cell membrane. a is related to the excitation threshold (i.e., the minimum V value that leads to the onset of an AP). The model's APD and refractoriness can, in turn, be controlled using b. The values for each of the model parameters are typically chosen empirically to reproduce observed electrical signals - we use the values listed in Supplementary Table 1. The Aliev-Panfilov model uses rescaled units: V is adimensional (typically in the [0, 1] AU interval) and time is measured in temporal units, referred to as TU throughout this study. 1 TU corresponds to approximately 13 ms (7).

Other more complex EP models exist, through which it is possible to model individual membrane ionic currents and other biological components relevant for the AP and its propagation. One example used in the current study is a 14-current 30-variable canine atrial model that incorporates different degrees of EP remodelling caused by atrial fibrillation (8).

By assigning different sets of parameters and/or initial conditions to EP mathematical models, they can represent the electrical behaviour of the heart in both healthy and arrhythmic conditions. Healthy conditions are usually represented as unidirectional smooth propagation from a single source (Figures 1a,b). Large sources produce wave fronts which are close to planar, whereas point-like wavefronts lead to centrifugal (convex) wavefronts. Arrhythmias are usually modelled as one or more re-entrant waves (called spiral waves or rotors in 2D–Figure 1c). Moreover, the wavefronts and/or wavebacks of spiral waves can, in some instances, fragment (break-up), leading to complex activation patterns (9). These models can also consider localised pathology such as fibrosis, scar or ischaemia as heterogeneities in one or more model parameters. Localised reductions in D, for instance, can be used to reproduce the slow-down of AP propagation in fibrotic lesions (10). These heterogeneities can lead to local changes in the curvature of the activation wavefront (Figures 1d,e).

FIGURE 1

Figure 1. Numerical solutions to the Aliev-Panfilov monodomain system for: (a) Planar wave. (b) Centrifugal wave. (c) Spiral wave. (d) Centrifugal wave in the presence of a square heterogeneity in D. (e) Spiral wave in the presence of a square heterogeneity in D.

2.2. Parameter Estimation in Cardiac Electrophysiology

Biophysical models are usually employed in forward mode, with the aim of reproducing the system's behaviour assuming that all model parameters are perfectly known. In many circumstances, such as the identification of pathology from EP measurements, it is more desirable to use these biophysical models in inverse mode, inferring the tissue parameters that underlie an observed system behaviour. This task, often called parameter estimation (or systems identification), is extremely challenging for several reasons.

The observed data are typically insufficient to identify a unique parameter value, since the observations are usually sparse, incomplete and polluted by noise. This can be handled by optimisation methods, such as least squares, that fit data and model in an optimal way, combined with problem-dependent regularisation terms to stabilise the estimate (11). Typically, the optimisation requires many forward runs of the underlying model. Unfortunately, in physical systems that are described by PDEs such as in EP, the forward runs are very computationally expensive in realistic 2D and 3D settings. These difficulties are exacerbated by the fact that the parameters in many EP models are heterogeneous. Hence, we often do not infer a scalar quantity but a (discrete) function in space and/or time. This increases the cost and complexity of the inverse mode even further.

Modern inverse problem solvers such as ensemble Kalman filters (12) sequential Monte Carlo (13) and parameter estimation based on Markov chain Monte Carlo (14) are often combined with reduced order models (15) or multi-fidelity approaches (16) to decrease the computational complexity of parameter estimation. These very complex methods typically make strong assumptions about the statistical distribution of model parameters and require dedicated problem-specific parameterisation. It is not clear how well they can generalise when applied to a different EP model or task.

2.3. Physics Informed Neural Networks (PINNs)

Physics Informed Neural Networks (PINNs) (17) are an exciting new tool for the study of physical systems modelled by PDEs and/or ODEs. PINNs have been shown to both efficiently find high-resolution solutions (forward problem) and perform parameter estimation (inverse problem) in a variety of systems (18). As opposed to most types of neural networks (NNs), whose inputs are exclusively empirical data, PINNs incorporate explicit knowledge about the physical laws that govern a system. This allows PINNs to compute solutions to initial and/or boundary value problems with comparatively less training data than conventional NNs.

Very briefly, PINNs make use of automatic differentiation (19) to rewrite the differential equations that a system obeys as the minimisation of a functional f. For example, for Equations (1, 2), this functional would be defined as:

$\begin{array}{ll}{f}_V=-\frac{\partial V}{\partial t}+\overrightarrow{\nabla }.(D\overrightarrow{\nabla }V)-kV(V-a)(V-1)-VW& (3)\end{array}$

$\begin{array}{ll}{f}_W=-\frac{\text{d}W}{\text{d}t}+(ϵ+\frac{{\mu }_{1}W}{V+{\mu }_2})(-W-kV(V-b-1))& (4)\end{array}$

PINNs are trained to minimise a hybrid loss function, L, which ensures the system obeys known physical laws (described by f_V and f_W), whilst simultaneously fitting known empirical system measurements. L thus includes terms to account for:

• agreement with the experimental measurements, L_data;

• consistency with the physical laws of the system, L_fV + L_fW;

• consistency with boundary L_VBC and initial value conditions L_VIC.

Mathematically:

$\begin{array}{ll}L=L_{data}+L_{f_V}+L_{f_W}+L_{V_{BC}}+L_{V_{IC}}& (5)\end{array}$

$\begin{array}{l}L=\frac{1}{N}{\displaystyle \sum _{i=1}^N}{(V(x_i,t_i)-{V_{GT}}_i)}^2+\frac{1}{N_f}{\displaystyle \sum _{j=1}^{N_f}}(f_V{(x_j,t_j)}^2+f_W{(x_j,t_j)}^2)+\\ \frac{1}{N_b}{\displaystyle \sum _{k=1}^{N_b}}{(\frac{\partial V}{\partial \overrightarrow{n}}(x_k,t_k))}^2+\frac{1}{N_0}{\displaystyle \sum _{l=1}^{N_0}}{(V(x_l,t_0)-V_0)}^2& (6)\end{array}$

Each of the terms of the loss function is typically computed in different domains:

• (x_i, t_i) are the N measurement points, where ground truth (GT) experimental measurements, V_GT, are known and should be reproduced as closely as possible.

• (x_j, t_j) are the N_f residual points, where the fulfilment of the biophysical equations is tested.

• (x_k, t_k) are the N_b boundary points, where the network aims to fulfil the Neumann boundary condition for V.

• (x_i, t₀) are the N₀ initial points, where the known initial condition, V₀, is replicated as closely as possible.

There is an asymmetry in L between measurable (V) and latent (W) variables: only experimental measurements (and initial conditions) for V are usually available. Moreover, as W does not diffuse across the tissue, it does not obey any boundary conditions.

PINNs have typically been used for two main purposes (18). In the so-called forward mode, the NN's parameters are optimised to provide a representation of the physical system of interest consistent with observations. Using this representation, the system's differential equations can subsequently be solved at an arbitrarily high spatial or temporal resolution, bypassing the constraints (e.g., small temporal and spatial steps) of traditional numerical solvers. In inverse mode, PINNs additionally perform parameter inference (systems identification) by having the NN optimise one or more of the equation parameters (which here represent tissue EP properties) during the training process.

2.4. Cardiovascular Applications of PINNs

PINNs have recently been used in several areas of cardiovascular medicine, especially for applications related to blood flow. Examples include the estimation of myocardial perfusion and related physiological parameters from dynamic contrast enhanced MRI (20) and the estimation of haemodynamic parameters from microscopic images of aneurysms-on-a-chip (21).

In the field of cardiac EP, Sahli-Costabal et al. (22) used PINNs in forward mode to estimate activation time maps (ATs, i.e., the arrival times of the action potential) and conduction velocity (CV) maps in the left atrium at high spatial resolution. Sahli-Costabal's method uses PINNs to solve the (isotropic diffusion) eikonal equation, a simple relationship between ATs and the spatial gradient of CV. This effectively interpolates AT and CV across the left atrial surface. Although their PINNs implementation was exclusively deployed on simulated data, the proposed application is very clinically relevant, as it aims to mitigate the low spatial resolution of clinical AT measurements. Grandits et al. (23) subsequently extended this PINNs-eikonal equation framework to anisotropic conduction, using it to estimate high-resolution AT maps and fibre directions from in silico and patient data. The PINNs method performance was nevertheless lower than that of a traditional (variational) inverse solver (24). As they rely on the eikonal equation, these tools are not well suited to the study of arrhythmic conditions or to the inference of EP parameters other than AT and CV.

2.5. Optical Mapping for Experimental EP-PINNs Testing

Maps of transmural electrical potential (V), similar to those simulated using the Aliev-Panfilov model, can be experimentally recorded using optical mapping. Optical mapping is a technique in which voltage-sensitive fluorescent dyes are added to cardiomyocyte preparations before imaging at high spatio-temporal resolution (25). It can be used to effectively obtain uncalibrated measurements of $V(\overrightarrow{x},t)$ in cardiac tissue across time. Although optical mapping can be challenging in vivo, in vitro experiments can provide very detailed insights into AP properties and cellular-level EP properties and gain insights into arrhythmic mechanisms (26). In particular, optical mapping can be used to study the effect of anti-arrhythmic drugs on cardiomyocyte preparations (27). These data were used to test EP-PINNs in an experimental setting.

3. Materials and Methods

In this section, we provide details about the finite differences (FD) model used, in a variety of settings, to generate training and test data for EP-PINNs. We then introduce the EP-PINNs architecture, before giving details about each of the in silico experiments in which EP-PINNs were deployed. We end this section by introducing the experimental data (in vitro optical mapping) used for testing EP-PINNs. The code used for EP-PINNs implementation and the generation of in silico GT data is freely available from github.com/martavarela/EP-PINNs.

3.1. Generation of in silico EP Data

We use an FD solver to generate in silico cardiac EP data, V_GT, which we use to train and evaluate the performance of EP-PINNs. Simulations were carried out in two different geometries: in 1D, in a 2-cm 1D domain cable; and in 2D, in a square with a 1-cm side. We use in-house code written in Matlab 2020b (Mathworks, Natick, MA, USA) that relies on central FD and an explicit 4-stage Runge-Kutta method to solve the isotropic monodomain Aliev-Panfilov model, as defined in Equations (1, 2) and model parameters listed in Supplementary Table 1. All simulations use Neumann boundary conditions and set dt = 5 × 10⁻³TU and dx = 100μm as temporal and spatial steps. All simulations were run for 300 TU, with the calculated $V_{GT}(\overrightarrow{x},t)$ field saved at every 1 TU and at every spatial step (every 100μm). We thus generate in total 1.4 × 10⁴ and 7.0 × 10⁵ data points for $V_{GT}(\overrightarrow{x},t)$, respectively in 1D and 2D.

Figure 1 shows example time frames of the generated V maps in 2D. All APs were initialised by adding an external stimulus current V_stim = 0.12 AU to the right-hand side of Equation (1) for 1 TU in a sub-domain of the studied geometry (Figure 1). In 2D simulations, this includes both planar (Figure 1a) waves (emanating from a rectangular stimulus) and centrifugal waves (from a point-like excitation) (Figure 1b). Spiral waves were also created using the cross-field protocol: as an initial planar excitation propagates, a second planar excitation wave, orthogonal to the first, is initiated. When timed appropriately (42 TU after the first stimulus, in our model), the second planar wave continuously curves as it moves to non-refractory tissue, giving rise to a sustained spiral wave (Figure 1c).

The data generated by the FD model is treated as GT data and is used to both train and test EP-PINNs. As in past studies (18, 28), we train EP-PINNs on a case-by-case basis, with a small subset of the data for which parameter inference is going to be performed. This is in contrast to most supervised NNs, which are usually trained with large amounts of data acquired in varied circumstances.

Training data for EP-PINNs, $V_{GT}(\overrightarrow{x},t)$, are provided to the network as an input and used to minimise the data discrepancy loss function, L_data in Equation (6). They amount to 10–20% of the total generated data, corresponding to a variable number of data points as detailed below and in Supplementary Table 2. The remaining data are withheld from EP-PINNs and used, in a post-processing step, to assess EP-PINNs' ability to reproduce cardiac APs, as detailed below.

3.2. Architecture and Training of EP-PINNs

EP-PINNs are designed, trained and deployed using the Python DeepXDE library (28). As in past successful implementations of PINNs (18, 28), we use a fully connected network architecture.

As detailed in Figure 2, EP-PINNs take as inputs the spatio-temporal points, $(\overrightarrow{x},t)$, where they will estimate the main outputs: V (and W). Experimental measurements of V, $V_{GT}(\overrightarrow{x},t)$, are also provided to EP-PINNs as inputs in a (training) subset of $(\overrightarrow{x},t)$. EP-PINNs minimise the hybrid physics-informed loss function described before (see Equations 5 and 6), by adjusting the network's weights and biases (collectively named θ in Figure 2). In inverse mode, EP-PINNs additionally adjust one or more parameters of the Aliev-Panfilov model (generically λ in Figure 2). The number of layers and neurons used by EP-PINNs is adjusted to the domain size and the complexity of the problem at hand, as described below and detailed in Supplementary Table 2.

FIGURE 2

Figure 2. Network architectures and training schemes for EP-PINNs used in this study. Neurons are represented by σ. In forward simulations only network parameters θ were estimated, whereas in the inverse setting one or more EP parameters (λ) were also calculated. (A) Architecture used in all forward simulations and all inverse simulations (1D and 2D) with homogeneous EP parameters. (B) Architecture used in all 2D forward simulations and all inverse simulations where D was a spatially-varying field. (C) Used training schemes. 1 was used in in silico 1D problems, whereas training scheme 2 was used in 2D problems and for optical mapping experimental data. The number of NN neurons and layers varied across different experiments, as did the learning rates (lr) and number of iterations (iter) in training scheme 2. Further details about the architecture and parameterisation of the NNs and training schemes can be found in Supplementary Table 2.

The optimisation approach used for EP-PINNs also depends on problem size. For 1D problems, we use Adam optimisation (29)—see Figure 2C. We empirically determined that, in 2D, EP-PINNs' performance improved when initially using Adam optimisation for only the data agreement term (L_data in Equation 6), followed by Adam optimisation for the full loss function and ending in a final phase of L-BFGS optimisation (30)—see Figure 2C. The initial Adam training phases are used to rapidly approach the desired minimum and the final L-BFGS optimisation phase helps the network converge faster towards it (28).

In the presence of spatially-varying EP parameters (D in the current study), we use network architecture B (see Figure 2B). Here, $D(\overrightarrow{x},t)$ is estimated by a parallel NN, NN_D, with the same number of layers and neurons as the main NN. In this setup, D is treated as a system variable (on par with V and W) and its estimates directly contribute to the loss term that ensures the agreement with the EP equations: L_fV in Equation (6).

The hyperbolic tangent function (tanh) is used throughout as the differentiable activation function and Glorot initialisation from a uniform distribution is used for all weights (31). Additionally, to minimise convergence problems caused by explosive gradients (28) and enhance the NN's stability, we implemented an automatic reset of the training process when the losses at the first epoch of the training exceeded a predefined threshold.

EP-PINNs were trained on a high performance machine with 1 RTX6000 GPU and 4 AMD EPYC 7742 CPUs. Typical training times varied between 15 min (for 1D problems) and 16h (for heterogeneous spiral wave problems), as detailed in Supplementary Table 2.

3.2.1. Assessment of EP-PINNs Performance

To assess the performance of EP-PINNs, we calculate, across all test points N_test, the root mean squared error (RMSE) for estimates of V:

$\begin{array}{ll}RMSE=\sqrt{\frac{1}{N_{test}}{\displaystyle \sum _{i=1}^{N_{test}}}{(V(i)-V_{GT}(i))}^2}& (7)\end{array}$

The RMSE is, by construction, adimensional and in the same scale range as V. All experiments were repeated at least 5 times to probe the variability in RMSE. In inverse mode, we additionally calculate the precision of the estimated model parameters using the standard deviation of the parameters estimated by EP-PINNs in these 5 different runs.

3.3. Forward Solution of EP Models

3.3.1. 1D Cable Geometry

We assessed the accuracy of EP-PINNs when solving the monodomain equation with the Aliev-Panfilov ionic model (forward problem) in the 20-cm cable during 70 TU (corresponding to 903 ms). We divided the GT data from the FD solver, V_GT, into test and training datasets using a 90–10% split. This corresponds to 140 randomly chosen points across the temporal and spatial domains for training EP-PINNs and 1,260 points for testing. As for inverse 1D problems, EP-PINNs were implemented in this instance using architecture A and training scheme 1 (see Supplementary Table 2).

We used two different training setups:

1. using only in silico experimental V_GT measurements as GT, as in Equation (6).

2. using in silico experimental points for V_GT and W_GT in the loss function, by adding an extra term: $L_{W_{GT}}=\frac{1}{N}\sum _{i=1}^N{(W(x_i,t_i)-W_{GTi})}^2$ to Equation (5).

Setup 1 more closely resembles an experimental setup, as the latent variable W_GT is not usually measurable.

To gauge whether EP-PINNs performance depended on model parameter choice, we used EP-PINNs on GT data synthesised with two different sets of model parameters, as detailed in Supplementary Table 1.

We additionally evaluated the performance of EP-PINNs in in silico data corrupted by noise. For this, we added to V_GT zero-mean Gaussian noise with standard deviations of 0.05, 0.10, 0.50 or 1.00 AU (with 1 AU being the approximate amplitude of an AP).

We also tested the performance of EP-PINNs in the presence of a reduced number of training data points. For this, we provided the network with V_GT at 1 × 10⁴, 5 × 10³, 1 × 10³ or 100 random training points within the 1D space-time domain (compared to 1.40 × 10⁷ in usual conditions). We used architecture A with training scheme 1 in all 1D problems (see Supplementary Table 2).

3.3.2. 2D Rectangular Geometry

We solved the forward problem in 2D, using the Aliev-Panfilov model in an isotropic and homogeneous 10-cm side square over 70 TU. We used 1.4 × 10⁵ randomly chosen data points (across the temporal and spatial domains) from the FD solver for the training of EP-PINNs (corresponding to 20% of the total generated points) and reserved the remaining 5.6 × 10⁵ points to assess EP-PINNs' performance. Simulations were carried out in three scenarios: planar wave (Figure 1a), centrifugal wave, emanating from point-like excitation in a corner (Figure 1b) and a spiral wave (Figure 1c). Here and in the equivalent inverse setup, EP-PINNs were implemented using architecture A and training scheme 2 (see Supplementary Table 2).

The point-like excitation and spiral wave scenarios were also simulated in heterogeneous conditions, where a 2-mm side square within the spatial domain was assigned a permanent diffusion coefficient ($D_{lesion}=0.02m{m}^2/TU$) lower than that of background tissue ($D_0=0.1m{m}^2/TU$). Architecture B with training scheme 2 was used for all (forward and inverse) heterogeneous problems (see Supplementary Table 2). A sigmoid function was used in the NN dedicated to estimating $D(\overrightarrow{x},t)$ (NN_D in Figure 2B) to account for the fact that $D(\overrightarrow{x},t)$ follows a binomial distribution: D₀ in healthy tissue and D_lesion otherwise.

3.4. Inverse Estimation of EP Parameters

We withheld the value of one or more of the model parameters from EP-PINNs, which were instead estimated by it. These parameters were chosen for their comparatively simple biophysical interpretation, known susceptibility to both disease remodelling and pharmacological action, and the limited degree of mathematical coupling between them. They were:

• a, which is related to the tissue excitation threshold (see Equation 1);

• D, the scalar diffusion coefficient (proportional to the electrical conductivity of the tissue, see Equation 1);

• b, which controls APD, see Equation 2.

• a and D simultaneously;

• b and D simultaneously.

3.4.1. Homogeneous 1D and 2D Geometries

We solved the inverse problem in the same (1D or 2D) setup, architecture and training scheme as for the forward problems and using a similar division of randomly chosen data points for training and testing. We assumed that none of the selected parameters varied across time or space (except for D, in the heterogeneous problem described below) and used the values listed in Supplementary Table 1 as GT values. In 1D, we additionally investigated how EP parameter estimation was affected by experimental noise. For this, Gaussian noise (σ = 0.05 or 0.10 AU) was added to the in silico data as described for the forward mode in section 3.3.

In addition to EP-PINNs robustness in the presence of experimental noise, we are also interested in its ability to cope with model uncertainty. Therefore, to assess EP-PINNs' ability to generalise beyond the model it is trained on, we additionally tested it on APs generated on a much more complex canine atrial EP model (8). These atrial APs are markedly different from those of the Aliev Panfilov model the EP-PINNs assumes, both in morphology and restitution properties. The canine atrial model data were synthesised using Matlab (https://models.cellml.org/workspace/47c) with central FD and explicit forward Euler schemes (dt = 5μs and dx = 100μm), with data saved at every spatial step and at every ms. Using this model, we created GT data for left atrial cells at 3 stages of AF-induced remodelling, which differed in APD. We tested, using the 1D model in inverse mode, whether EP-PINNs could identify the reduction in APD (detected as an increase in parameter b) in left atrial APs caused by increasing amounts of AF remodelling.

3.4.2. Estimation of EP Parameter Heterogeneities

We assessed EP-PINNs' ability to recognise spatial heterogeneities in model parameters in 2D, as a test for EP-PINNs potential for identifying spatially-varying lesions such as fibrosis. For this, we used the same setup as in section 3.3.2, with $D_0=0.1m{m}^2/TU$ reduced to $D_{lesion}=0.02m{m}^2/TU$ in a similar square region. As before, we estimated $D(\overrightarrow{x},t)$ on its own and simultaneously with either the a or b global model parameters.

3.5. Parameter Estimation Using Optical Mapping Data

We tested EP-PINNs performance on in vitro datasets using optical mapping data from neonatal ventricular rat myocyte preparations stained with a voltage-sensitive dye, as described in detail by Chowdhury et al. (27). Briefly, we used four time series (movies) of optical mapping images (field of view: 4.1 × 0.1 mm, spatial resolution: 1.2 μm, temporal resolution: 2 ms, duration: 300 ms). In two of these image series, ionic channel modulating drugs (E-4031 or nifedipine) had been administered at half maximal inhibitory concentration (IC₅₀): 772.2 nM for nifedipine and 243.4 nM for E-4031. The other two temporal image series consisted of matched control (baseline) preparations, to which no drug had been given.

We manually selected two square regions of interest (ROIs) with a side of 2.3 μm and with their centres 0.7 mm apart, in the same image location for each time series. We spatially averaged the optical signal over each ROI to obtain a signal trace across time. From this signal, we manually selected two consecutive APs and normalised the signal to the [0, 1] interval for consistency with the Aliev-Panfilov model. To improve the signal to noise ratio (SNR) of this trace, we applied a mean average filter twice, aligned and averaged the two APs over time to obtain a single higher-SNR AP. These pairs of post-processed APs were used as inputs to 1D EP-PINNs in inverse mode, with b as the variable to be estimated. All model parameters were unchanged from those in Supplementary Table 1. We used EP-PINNs' architecture A with training scheme 2 to estimate b 10 times for each setting. 116-200 points were used for training EP-PINNs and 29–50 to test it, as detailed in Supplementary Table 2.

We investigated in particular whether EP-PINNs could detect the effect on APD of E-4031 and nifedipine, which, respectively, block the hERG voltage-gated potassium channel (I_Kr current) and the L-type calcium channel (I_CaL current). Whereas, E-4031 will extend APD (and thus shorten b in the Aliev-Panfilov model—see Equations 1, 2), nifedipine will have the opposite effect, decreasing APD (increasing b). We compared EP-PINNs' b estimates for i) E-4031 vs. baseline and ii) nifedipine vs. baseline and assessed, using a t-test, whether the administered drugs significantly changed b values.

4. Results

EP-PINNs successfully solved the monodomain equation coupled with the Aliev-Panfilov model in 1D and 2D, including in the presence of heterogeneities in D. In inverse mode, the proposed setup could also successfully perform EP parameter estimation from in silico and in vitro data, as described below.

4.1. Forward Solution of EP Models

4.1.1. 1D Cable Geometry

EP-PINNs accurately reproduced the features, morphology and conduction properties of the APs generated by the Aliev-Panfilov model (Figure 3A). We found that EP-PINNs could accurately simulate APs even in the absence of GT values for the latent variable W. Although the error was slightly increased in the absence of W_GT (RMSE = 6.0 × 10⁻³ ± 2.0 × 10⁻³ vs. 9.0 × 10⁻³ ± 4.0 × 10⁻³), the RMSE was still minimal (see top left inset in Figure 3A) and the estimated $V(\overrightarrow{x},t)$ were visually indistinguishable from GT traces in both cases. As a consequence, in all other experiments described in this study, GT data for W was not used to train EP-PINNs, whose training relied solely on V_GT.

FIGURE 3

Figure 3. Effect of the number of sampled experimental points (A) and experimental noise (B) on EP-PINNs V estimates. The error in estimation is measured using the root mean square error (RMSE, see Equation 7) in forward mode in 1D. Representative V(t) plots sampled at a random spatial location are shown as insets for some of the probed conditions.

We found that EP-PINNs could solve the model accurately in the presence of even small numbers of V_GT points for training, with RMSE ≤ 2.5 × 10⁻² even when trained with only 100 points (Figure 3A). As expected, increasing noise in V_GT led to increasing levels of error in V estimates (Figure 3B), but EP-PINNs were able to converge in the presence of Gaussian noise with a standard deviation below 0.5 AU (approximately 0.5 times the amplitude of an AP). When using a different set of parameters for the Aliev-Panfilov model (see Supplementary Table 1), EP-PINNs' accuracy was comparable, suggesting that the obtained EP-PINNs performance is robust to different biophysical model settings.

4.1.2. 2D Rectangular Geometry

In homogeneous conditions, EP-PINNs were also able to reproduce AP propagation in 2D for planar, centrifugal and spiral waves (see Figures 4a–c), with excellent accuracy (RMSE < 3.0 × 10⁻² throughout). In the presence of heterogeneities in D, EP-PINNs were also able to accurately simulate APs with an RMSE of 5.6 × 10⁻³ ± 6.6 × 10⁻⁴ for centrifugal waves, which increased slightly to 2.7 × 10⁻² ± 4.3 × 10⁻³ for spiral ones (Figures 4a–c). In the spiral wave scenario, EP-PINNs found it most difficult to reproduce V in the high wavefront curvature regions close to the spiral wave tip.

FIGURE 4

Figure 4. EP-PINNs 2D forward solutions to the Aliev-Panfilov monodomain system for the same conditions as the GT data depicted in Figure 1. (a) Planar wave. (b) Centrifugal wave. (c) Spiral wave. (d) Centrifugal wave in the presence of a square heterogeneity in D. (e) Spiral wave in the presence of a square heterogeneity in D.

Movies showing the propagation of APs in the 2D rectangular domain across time (for both GT and EP-PINNs solvers) can be seen in Supplementary Videos 1–5.

4.2. Inverse Estimation of EP Parameters

4.2.1. Homogeneous 1D and 2D Geometries

EP-PINNs were able to estimate global model parameters in 1D the presence of varying degrees of noise, as detailed in Figure 5. As for the forward problem (Figure 3), AP morphology and main properties were well reproduced even in the presence of large amounts of noise, with RMSE < 9.0 × 10⁻³ overall. Example plots of V(t) in inverse mode in 1D are shown in Supplementary Figure 1.

FIGURE 5

Figure 5. Error in global EP parameter estimates by EP-PINNs in the inverse setting in 1D in the presence of different amounts of experimental noise. We show the relative error for a, b and D, when estimated separately and in pairs.

When estimating only one model parameter, relative errors (RE) did not exceed on average 27%, even in the presence of noise. Estimates of b, which determines AP duration, were the most accurate (|RE| < 3%), followed by D (|RE| < 35%), which EP-PINNs tended to underestimate. a estimation was considerably more difficult. Joint estimation of two parameters led in general to less accurate parameter estimates (Figure 5), with an error as high as 100% for a when estimated jointly with D in the presence of noise (see Figure 5). When performing simultaneous estimation of two parameters, no evidence of coupling between them was observed.

EP-PINNs were additionally able to perform robust parameter estimation on synthetic experimental data generated by a different EP model (8). When estimating b in APs generated by a different atrial EP model, it correctly inferred that APD is reduced (i.e., b is increased) for increasing degrees of AF remodelling, as shown in Figure 6. Moreover, the main AP features were well reproduced, with small discrepancies related to the differences between the two models (Figure 6). Interestingly, the solution proposed by EP-PINNs consistently shows a less steep depolarisation than expected from either Aliev-Panfilov model and the detailed canine model. This mismatch is likely to be a consequence of the EP-PINNs' adjustment to slightly different AP morphologies from those in the Aliev-Panfilov model in its loss function. This suggests that the cross-model estimation of parameters related to excitability (e.g., a) may not be as successful as parameters related to APD (such as b).

FIGURE 6

Figure 6. 1D inverse EP-PINNs solution for a detailed canine left atrial model in the conditions of: (A) no remodelling. (B) moderate AF remodelling. (C) severe AF remodelling. Representative V(t) plots are shown throughout, accompanied by the models 90% APD and EP-PINNS estimates for b, a parameter inversely proportional to APD.

Global parameter estimation in 2D was again successful, for both unidirectional propagation and spiral wave conditions, as demonstrated in Figure 7. As in 1D, V was also correctly reproduced in the different analysed conditions (Figures 7b–d), with RMSE < 2.3 × 10⁻² throughout. Across all experiments, parameter estimation and V reconstruction were most successful for planar wave conditions followed by centrifugal waves and less accurate for spiral waves, as shown in Figure 7a. The comparatively worse performance of EP-PINNs in spiral wave conditions may be caused by the spatial dependency of wave front curvature in this setting.

FIGURE 7

Figure 7. EP-PINNs inverse solution in homogeneous conditions in 2D. (a) Relative error for global estimates of: a and D and b and D, when estimated separately or simultaneously. (b–d) Corresponding representative V maps for: planar wave (b), centrifugal wave (c), and spiral wave (d). Compare (b–d) to the corresponding GT in Figures 1a-c and the forward solutions in Figures 4a–c.

Errors were largest when estimating two parameters simultaneously, with EP-PINNs again struggling to estimate a, especially when in conjunction with D (|RE|_a < 100%). Estimates of b were once again the most accurate (|RE|_b < 8%) and, as in 1D, EP-PINNs tended to underestimate D across all settings and to underestimate all parameters in the spiral wave scenario.

4.2.2. Estimation of EP Parameter Heterogeneities

EP-PINNs were able to estimate D on a pixel-by-pixel basis with remarkable accuracy, as demonstrated in Figures 8c,d, accurately identifying the low D region. RMSE_D was consistently below 3.5 × 10⁻² (Figure 8e) and was lower for the centrifugal wave case than the spiral wave one. As before, V was similarly well reproduced (Figures 8a,b), especially in the centrifugal wave scenario (see Figures 8b,d), with RMSE < 3.0 × 10⁻². In the spiral wave scenario, EP-PINNs found the estimation of D hardest near the spiral tip (see Figure 8d), where the high wavefront curvature may resemble the wavefront bending caused by low D regions.

FIGURE 8

Figure 8. EP-PINNs 2D inverse solution in the presence of heterogeneities in D. Maps showing representative V and D estimates for: (a,c) centrifugal wave and (b,d) spiral wave. (e) Error for global estimates of a and b and RMSE for estimates of D across the 2D domain, for all estimated parameter combinations. Compare (a,b) to the corresponding GT in Figures 1d,e and the forward solutions in Figures 4d,e.

Using architecture B, global estimation of a and b was also possible in the presence of the heterogeneous D field, keeping the same trends as in the 1D and 2D homogeneous cases (Figure 8e). As before, no evidence of coupling between the simultaneously estimated parameters was found. In detail:

• b was very accurately estimated (|RE| < 10%), whereas a estimates had a larger error (|RE| < 100%);

• estimating a and b simultaneously with the D field led to small decreases in accuracy overall;

• estimating V and D in the spiral wave setting was more difficult than in the presence of a centrifugal wave.

4.3. Parameter Estimation using Optical Mapping Data

Using optical mapping signals, EP-PINNs were able to accurately reproduce experimental APs and identify the actions of nifedipine and E-4031, correctly estimating that they, respectively, reduce and increase APD (Figure 9). The reduction in APD caused by nifedipine, an I_CaL blocker, was detected by EP-PINNs as a significant increase in b in the Aliev-Panfilov model [Δb = (2.7 ± 0.9) × 10⁻², p < 10⁻⁵]. The E-4031-induced increase in APD was more subtle [Δb = (−2.1 ± 3.0) × 10⁻²] and non-significant (p = 0.38). The reduced effect of E-4031 in these data is consistent with the modest role I_Kr, the current blocked by E-4031, is expected to play in rodent APs (32).

FIGURE 9

Figure 9. 1D inverse EP-PINNs solution for experimental optical mapping data, in the presence of (A) I_CaL channel blocker nifedipine and (B) I_Kr channel blocker E-4031. (Δb refers to the change in model parameter b in the presence of the drug when compared to the baseline value. Each AP was acquired in separate datasets and juxtaposed in the figure to allow visual comparisons.

5. Discussion

We present EP-PINNs, a successful framework to estimate EP parameters from measurements of trans-membrane potential V. We demonstrate EP-PINNs ability to accurately reproduce AP propagation in 1D and 2D in the presence of very sparse experimental measurements, experimental noise and model uncertainty. EP-PINNs can also estimate, for 1D and 2D in silico and in vitro data, global markers of APD, excitation threshold and/or conductivity (diffusion coefficient, D). We additionally show that EP-PINNs are further capable of identifying heterogeneities in EP parameters, such as D, even in arrhythmic conditions, showcasing their potential for clinically useful applications.

5.1. Forward Solution of EP Models

EP-PINNs offer a flexible and easy to implement framework for parameter estimation in EP. Sahli-Costabal et al. (22) and Grandits et al. (23) had already demonstrated PINNs' potential in cardiac EP by estimating high-resolution left atrial AT and CV maps in sinus rhythm conditions using a simple activation-only biophysical model. We extend PINNs' applications in EP by applying them to a more complex biophysical model, the monodomain Aliev-Panfilov model (7), which also captures restitution properties through the inclusion of a latent (non-measurable) variable, W. For the first time, we use a PINNs framework for the simulation of arrhythmic conditions, such as spiral waves, and for the estimation of parameters unrelated to the AT of V.

Importantly, we show EP-PINNs' are able to reproduce APs and perform parameter estimation in the absence of any data for W, which is not available experimentally. This bodes well for the deployment of PINNs for even more complex EP models, which use a higher number of latent variables to model individual ionic channels. This possibility is also supported by the work of Yazdani et al. (33), who successfully used PINNs for parameter estimation across several biological systems described by large sets of coupled ODEs.

We demonstrate PINNs' ability to describe AP dynamics in several circumstances. In 1D, EP-PINNs were able to reproduce APs even in the presence of very reduced amounts of experimental data (Figure 3A) and large amounts of noise (Figure 3B). PINNs' incorporation of explicit biophysical equations in the NN's loss function acts as an effective regulariser in EP problems, as demonstrated before in many other physical systems (17, 18). We note that these inherent regularisation properties allow PINNs to be trained with much lower amounts of training data than conventional NNs. The main drawback is the need for a comparatively time-intensive training on a case-by-case basis, compared to the global training usually employed with supervised NNs.

In 2D, we were able to accurately replicate AP dynamics for planar, centrifugal and planar waves, even in the presence of heterogeneities in the diffusion coefficient (see Figures 4, 8). We found that a more sophisticated training scheme and, for spiral waves, an increased NN capacity (5 layers of 64 neurons vs. 4 layers of 32 neurons, see Supplementary Table 2) were necessary for convergence in these large and complex problems. This more complex setup could, of course, have been used to solve the simpler 1D problems, through a trade-off between computational time and the convenience of a one-size-fits-all EP-PINNs approach.

5.2. Inverse Estimation of EP Parameters

It is in inverse mode, when estimating model properties from sparse measurements of V, that the EP-PINNs framework showcases its usefulness. Parameter estimation is an important topic in EP, as it is essential for both the personalisation of models and for the understanding of the effect of pathology and drugs on APs. Although parameter estimation in EP has been extensively discussed as a means of reducing the uncertainty associated with current biophysical models (34), NNs had not yet been used for dedicated model parameter estimation in EP. This contrasts with the more common use of NNs as efficient solvers of EP systems (in a similar fashion to the EP-PINNs forward mode in the current study) (15, 35, 36).

Using EP-PINNs, we were able to estimate, with different degrees of accuracy, three different biophysical parameters, each controlling, in an almost uncoupled manner, different observable properties of the system: APD (through b), excitability (through a) and conduction velocity (through D). In both 1D and 2D (homogeneous and heterogeneous) problems, the network was highly successful at estimating b, but struggled with a, especially when estimating it in tandem with D. This is likely to reflect a dependency between EP-PINNs' inverse mode effectiveness and the solution type that is probed experimentally. Indeed, when compared to b, the experimental inputs to the network ($V_{GT}(\overrightarrow{x},t)$), depend little on a providing the initial stimulus is supra-threshold. An exception could have been the spiral wave scenario, whose properties (e.g., the distance between spiral arms) depend strongly on model parameters such as a. Model parameters are more strongly coupled in the properties of spiral solutions, however, explaining the consistently lower accuracy of EP-PINNs' estimates in this scenario (Figures 7a,d). Parameter estimation in the spiral wave scenario may be improved when EP-PINNs are trained in longer time series, in which the spiral wave tip samples more of the spatial domain.

These issues underlie the difficulties of finding a single experimental design that allows for simultaneous accurate estimation of several EP parameters. A solution for this problem may be the training of PINNs using data from the same system acquired in different experimental conditions, perhaps by training separate NNs in parallel with a combined loss function, in an analogous manner to architecture B in this study (see Figure 2).

EP-PINNs were additionally able to identify heterogeneities in D across the 2D domain (see Figure 8), by estimating D in a separate parallel NN which shared terms of the loss function (L_data, see Equation 5) with the main NN. Characterising heterogeneities in EP parameters from electrical measurements is an interesting problem from a clinical point of view, as arrhythmias such as AF are often accompanied by heterogeneous changes in EP properties. Of these, regions of dense fibrosis [often modelled as areas with a reduced D (10)] are a promising candidate for personalised ablation sites (37), which may increase the overall efficacy of these procedures. EP-PINNs are thus well placed to help locate these putative ablation sites by identifying spatial heterogeneities in EP parameters such as D.

The inputs to the current EP-PINNs implementation, $V(\overrightarrow{x},t)$ are not, however, measurable clinically. Future work will modify EP-PINNs to instead perform parameter inference from extracellular electrical potentials, ϕ_e, which are regularly measured during clinical procedures using contact electrodes. ϕ_e can be interpreted as a weighted spatial convolution of the $\overrightarrow{\nabla }.(D\overrightarrow{\nabla }V)$ term in Equation (2) (38, 39), making the identification of localised EP changes more difficult. To take this into account, EP-PINNs designed for ϕ_e analysis may benefit from a move away from the current fully-connected architecture to incorporate, for example, convolutional layers.

5.3. Parameter Estimation Using Optical Mapping Data

An important point for future clinical applications of EP-PINNs is its ability to generalise beyond the details of the setup it is trained on. We showed that the proposed EP-PINNs implementation is model-agnostic, as it was able to perform robust parameter inference on in silico data generated by a much more complex atrial EP model (8) than the 2-variable ventricular one used in its loss function. In particular, EP-PINNs were able to correctly identify the decrease in APD (manifest as an increase in b) that is associated with increasing degrees of AF-induced remodelling in this model (see Figure 6).

In contrast to most previous PINN studies (17, 22, 28, 33, 40), we complemented the in silico studies with an assessment of the EP-PINNs performance on experimental biological data. Despite requiring a proportionally higher amount of training data than in silico experiments, EP-PINNs were able to cope well with the noise and artefacts unavoidably present in experimental data to identify the effect on APD of two different drugs: an I_Kr blocker and an I_CaL blocker. As for the data generated by a different mathematical model (Figure 6), EP-PINNs coped well with differences between the experimental data and the Aliev-Panfilov model, namely in resting membrane potential (see Figure 9). This ability to generalise well to data with different characteristics could be due to the use of PDE as a soft constraint (a term in the loss function) in the EP-PINNs framework, as well as the lack of assumptions about the distributions from which data come from.

As demonstrated for the in silico tests, the EP-PINNs framework can easily be extended to simultaneously infer the effect of drugs on more than one EP parameter, which may be useful for the characterisation and safety assessments of anti-arrhythmic drugs. These applications may further benefit from the training of EP-PINNs on more complex biophysical models, to obtain a more fine-grained characterisation of potential pharmacological (or pathological) effects.

5.4. Limitations and Future Plans

The current study aims to demonstrate the potential of PINNs within EP, as an initial necessary step towards clinical applications of this method. As such, we only trained EP-PINNs using one comparatively simple EP model, in 1D and 2D scenarios and for the estimation of a small number of EP parameters. We additionally did not test EP-PINNs in the chaotic or pseudo-chaotic scenarios of spiral wave break-up (9), which may be relevant for some arrhythmias.

Generalisations of the proposed framework to 2D/3D geometries representative of cardiac chambers, anisotropic conditions and more detailed EP biophysical models can be achieved by further increasing the capacity of the deployed EP-PINNs, with a concurrent increase in computational resources. However, promising and less resource-intensive applications for EP-PINNs may be the characterisation of pharmacological effects on AP or the identification of heterogeneities in EP properties from EGM signals.

Data Availability Statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Author Contributions

CH and AO wrote the EP-PINNs software, conducted the experiments, and analysed the results. RC acquired the experimental data and oversaw its analysis. AB, EU, RC, and NP contributed to the study design and the revision of the results. MV designed the study, wrote the GT data generation software, critically analysed the results, and drafted the manuscript. All authors critically reviewed the manuscript and gave final approval for publication.

Funding

This work was supported by the British Heart Foundation (RE/18/4/34215, RG/16/3/32175, PG/16/17/32069, and Centre of Research Excellence), the Imperial-TUM Seed Funding, the National Institute for Health Research (UK) Biomedical Research Centre and the Rosetrees Trust through the interdisciplinary award Atrial Fibrillation: A Major Clinical Challenge and the Generalitat Valenciana Conselleria d'Educació, Investigació, Cultura i Esport (ACIF/2021/205).

Conflict of Interest

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.

Publisher's Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We thank the Imperial College Research Computing Service (DOI: 10.14469/hpc/2232) for help with the computational setup of the project and technical support.

Supplementary Material

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

References

1. Hindricks G, Potpara T, Dagres N, Arbelo E, Bax JJ, Blomström-Lundqvist C, et al. 2020. ESC Guidelines for the diagnosis and management of atrial fibrillation developed in collaboration with the European Association for Cardio-Thoracic Surgery (EACTS). Eur Heart J. (2020). 42:373–498. doi: 10.1093/eurheartj/ehaa945

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Ganesan AN, Shipp NJ, Brooks AG, Kuklik P, Lau DH, Lim H, et al. Long-term outcomes of catheter ablation of atrial fibrillation: a systematic review and meta-analysis. J Am Heart Assoc. (2013) 2:e004549. doi: 10.1161/JAHA.112.004549

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Nattel S, Dobrev D. Controversies about atrial fibrillation mechanisms. Circ Res. (2017) 120:1396–9. doi: 10.1161/CIRCRESAHA.116.310489

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Calkins H, Kuck KH, Cappato R, Brugada J, John Camm A, Chen SA, et al. 2012 HRS/EHRA/ECAS expert consensus statement on catheter and surgical ablation of atrial fibrillation. J Intervent Cardiac Electrophysiol. (2012). 14:171–257.

PubMed Abstract | Google Scholar

5. Clayton RH, Bernus O, Cherry E, Dierckx H, Fenton F, Mirabella L, et al. Models of cardiac tissue electrophysiology: progress, challenges and open questions. Prog Biophys Mol Biol. (2011) 104:22–48. doi: 10.1016/j.pbiomolbio.2010.05.008

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Bueno-Orovio A, Kay D, Grau V, Rodriguez B, Burrage K. Fractional diffusion models of cardiac electrical propagation: role of structural heterogeneity in dispersion of repolarization. J R Soc Interface. (2014) 11:20140352. doi: 10.1098/rsif.2014.0352

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Aliev RR, Panfilov A. A simple two-variable model of cardiac excitation. Chaos Solitons Fractals. (1996) 7:293–301. doi: 10.1016/0960-0779(95)00089-5

CrossRef Full Text | Google Scholar

8. Varela M, Colman M, Hancox J, Aslanidi O. atrial heterogeneity generates re-entrant substrate during atrial fibrillation and anti-arrhythmic drug action: mechanistic insights from canine atrial models. PLoS Comput Biol. (2016) 12:e1005245. doi: 10.1371/journal.pcbi.1005245

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Fenton F, Karma A. Vortex dynamics in three-dimensional continuous myocardium with fiber rotation: filament instability and fibrillation. Chaos. (1998) 8:879. doi: 10.1063/1.166311

PubMed Abstract | CrossRef Full Text

10. Roy A, Varela M, Aslanidi O. Image-Based computational evaluation of the effects of atrial wall thickness and fibrosis on re-entrant drivers for atrial fibrillation. Front Physiol. (2018) 9:1352. doi: 10.3389/fphys.2018.01352

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Nelles O. Nonlinear System Identification From Classical Approaches to Neural Networks, Fuzzy Models, and Gaussian Processes. 2nd ed. Cham: Springer International Publishing (2020).

Google Scholar

12. Hoffman MJ, Cherry EM. Sensitivity of a data-assimilation system for reconstructing three-dimensional cardiac electrical dynamics. Philos Trans R SocA. (2020) 378:20190388. doi: 10.1098/rsta.2019.0388

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Drovandi CC, Cusimano N, Psaltis S, Lawson BAJ, Pettitt AN, Burrage P, et al. Sampling methods for exploring between-subject variability in cardiac electrophysiology experiments. J R Soc Interface. (2016) 13:20160214. doi: 10.1098/rsif.2016.0214

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Dhamala J, Arevalo HJ, Sapp J, Horácek BM, Wu KC, Trayanova NA, et al. Quantifying the uncertainty in model parameters using Gaussian process-based Markov chain Monte Carlo in cardiac electrophysiology. Med Image Anal. (2018) 48:43–57. doi: 10.1016/j.media.2018.05.007

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Fresca S, Manzoni A, Dedè L, Quarteroni A. Deep learning-based reduced order models in cardiac electrophysiology. PLoS ONE. (2020) 15:e0239416. doi: 10.1371/journal.pone.0239416

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Sahli Costabal F, Perdikaris P, Kuhl E, Hurtado DE. Multi-fidelity classification using Gaussian processes: Accelerating the prediction of large-scale computational models. Comput Methods Appl Mech Eng. (2019) 357:112602. doi: 10.1016/j.cma.2019.112602

CrossRef Full Text | Google Scholar

17. Raissi M, Perdikaris P, Karniadakis GE. Physics-informed neural networks: a deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations. J Comput Phys. (2019) 378:686–707. doi: 10.1016/j.jcp.2018.10.045

CrossRef Full Text | Google Scholar

18. Karniadakis GE, Kevrekidis IG, Lu L, Perdikaris P, Wang S, Yang L. Physics-informed machine learning. Nat Rev Phys. (2021) 3:422–40. doi: 10.1038/s42254-021-00314-5

CrossRef Full Text | Google Scholar

19. Baydin AG, Pearlmutter BA, Radul AA, Siskind JM. Automatic differentiation in machine learning: a survey. J Mach Learn Res. (2018) 18:1–43. Available online at: https://w.jmlr.org/papers/volume18/17-468/17-468.pdf

PubMed Abstract | Google Scholar

20. van Herten RLM, Chiribiri A, Breeuwer M, Veta M, Scannell CM. Physics-informed neural networks for myocardial perfusion MRI quantification. (2020).

Google Scholar

21. Cai S, Li H, Zheng F, Kong F, Dao M, Karniadakis GE, et al. Artificial intelligence velocimetry and microaneurysm-on-a-chip for three-dimensional analysis of blood flow in physiology and disease. Proc Natl Acad Sci USA. (2021) 118:e2100697118. doi: 10.1073/pnas.2100697118

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Sahli Costabal F, Yang Y, Perdikaris P, Hurtado DE, Kuhl E. Physics-informed neural networks for cardiac activation mapping. Front Phys. (2020) 8:42. doi: 10.3389/fphy.2020.00042

CrossRef Full Text | Google Scholar

23. Grandits T, Pezzuto S, Costabal FS, Perdikaris P, Pock T, Plank G, et al. Learning atrial fiber orientations and conductivity tensors from intracardiac maps using physics-informed neural networks. Lecture Notes Compu Sci. (2021) 12738:650–8. doi: 10.1007/978-3-030-78710-3_62

CrossRef Full Text | Google Scholar

24. Grandits T, Pezzuto S, Lubrecht JM, Pock T, Plank G, Krause R. PIEMAP: personalized inverse eikonal model from cardiac electro-anatomical maps. Lecture Notes Comput Sci. (2020) 12592:76–86. doi: 10.1007/978-3-030-68107-4_8

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Efimov IR, Nikolski VP, Salama G. Optical imaging of the heart. Circ Res. (2004) 95:21–33. doi: 10.1161/01.RES.0000130529.18016.35

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Hansen B, Zhao J, Li N, Zolotarev A, Zakharkin S, Wang Y, et al. Human atrial fibrillation drivers resolved with integrated functional and structural imaging to benefit clinical mapping. JACC: Clin Electrophysiol. (2018) 4:1501–15. doi: 10.1016/j.jacep.2018.08.024

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Chowdhury RA, Tzortzis KN, Dupont E, Selvadurai S, Perbellini F, Cantwell CD, et al. Concurrent micro-to macro-cardiac electrophysiology in myocyte cultures and human heart slices. Sci Rep. (2018) 8:1–13. doi: 10.1038/s41598-018-25170-9

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Lu L, Meng X, Mao Z, Karniadakis GE. DeepXDE: a deep learning library for solving differential equations. SIAM Rev. (2021) 63:208–28. doi: 10.1137/19M1274067

CrossRef Full Text | Google Scholar

29. Kingma DP, Ba J. Adam: a method for stochastic optimization. In: 3rd International Conference on Learning Representations. San Diego, CA: ICLR. (2015).

Google Scholar

30. Liu DC, Nocedal J. On the limited memory BFGS method for large scale optimization. Math Programm. (1989) 45:503–28. doi: 10.1007/BF01589116

CrossRef Full Text | Google Scholar

31. Glorot X, Bengio Y. Understanding the difficulty of training deep feedforward neural networks. In: Teh YW, Titterington M, editors. Proceedings of the Thirteenth International Conference on Artificial Intelligence and Statistics. Sardinia: PMLR (2010). p. 246–56. Available Online at: http://proceedings.mlr.press/v9/glorot10a/glorot10a.pdf

32. Odening KE, Gomez AM, Dobrev D, Fabritz L, Heinzel FR, Mangoni ME, et al. ESC working group on cardiac cellular electrophysiology position paper: relevance, opportunities, and limitations of experimental models for cardiac electrophysiology research. Europace. (2021) 23:1795–814. doi: 10.1093/europace/euab142

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Yazdani A, Lu L, Raissi M, Karniadakis GE. Systems biology informed deep learning for inferring parameters and hidden dynamics. PLoS Comput Biol. (2020) 16:e1007575. doi: 10.1371/journal.pcbi.1007575

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Clayton RH, Aboelkassem Y, Cantwell CD, Corrado C, Delhaas T, Huberts W, et al. An audit of uncertainty in multi-scale cardiac electrophysiology models. PhilosTrans R Soc A. (2020) 378: 20190335. doi: 10.1098/rsta.2019.0335

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Cantwell CD, Mohamied Y, Tzortzis KN, Garasto S, Houston C, Chowdhury RA, et al. Rethinking multiscale cardiac electrophysiology with machine learning and predictive modelling. Comput Biol Med. (2019) 104:339–51. doi: 10.1016/j.compbiomed.2018.10.015

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Kashtanova V, Ayed I, Cedilnik N, Gallinari P, Sermesant M. EP-Net 2.0: Out-of-domain generalisation for deep learning models of cardiac electrophysiology. In: Ennis DB, Perotti LE, Wang VY, editors. Functional Imaging and Modeling of the Heart. Cham: Springer International Publishing (2021). p. 482–92. doi: 10.1007/978-3-030-78710-3_46

CrossRef Full Text | Google Scholar

37. Roy A, Varela M, Chubb H, MacLeod R, Hancox JC, Schaeffter T, et al. Identifying locations of re-entrant drivers from patient-specific distribution of fibrosis in the left atrium. PLoS Comput Biol. (2020) 16:e1008086. doi: 10.1371/journal.pcbi.1008086

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Varela M, Aslanidi O. Role of atrial tissue substrate and electrical activation pattern in fractionation of atrial electrograms: a computational study. IEEE Eng Med Biol Soc Ann Conf. (2014) 2014:1587–90. doi: 10.1109/EMBC.2014.6943907

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Plonsey R, Barr R. Bioelectricity: a Quantitative Approach. 2nd ed. Berlin: Springer Science Business Media (2007).

Google Scholar

40. Arthurs CJ, King AP. Active training of physics-informed neural networks to aggregate and interpolate parametric solutions to the Navier-Stokes equations. J Comput Phys. (2021) 438:110364. doi: 10.1016/j.jcp.2021.110364

PubMed Abstract | CrossRef Full Text | Google Scholar