1 Introduction
Noninvasive electrophysiological (EP) imaging involves the reconstruction of cardiac electrical activity from highdensity bodysurface electrocardiograms (ECGs) [6]. It solves an illposed inverse problem that deteriorates as the imaging depth increases from the epicardium to the endocardium [9]. One type of increasingly utilized regularization considers knowledge about the welldefined physiological process of cardiac electrical propagation. This is often realized in a modelconstrained approach, where the optimization or statistical inference of cardiac electrical activity is constrained by a predefined model describing local activation/repolarization and its spatial propagation [4, 11, 12]. Earlier models include step jump functions [10], logistic functions [11], and 3D curve models [4] empirically parameterized to mimic the physiological process. Recently, more expressive cardiac EP simulation models have also been used [12, 7].
These modelconstrained approaches are afflicted with a common challenge: they are controlled by highdimensional parameters often associated with local tissue properties and the origin of electrical activation that are unknown a priori. The more expressive the model is, the more parameters it has. To fix these model parameters in optimization/inference, as is common in existing approaches [12]
, model errors may be introduced decreasing the accuracy of the estimated electrical activity
[12]. To adapt these model parameters to the observed data, as is desired for accurate inference, is however difficult due to their highdimensionality and nonlinear relationship with the observed ECG data [3].In this paper, we introduce a novel modelconstrained inference framework that replaces the conventional physiological models with a deep generative model that is trained to generate the spatiotemporal dynamics of transmembrane potential (TMP) from a lowdimensional set of generative factors. These generative factors can be viewed as a lowdimensional abstraction of the highdimensional physical parameters, which allows us to efficiently adapt the prior physiological knowledge to the observed ECG data (through inference of the generative factors) for an improved reconstruction of TMP dynamics.
In specific, the presented method consists of two novel contributions. First, to obtain a generative model that is sufficiently expressive to reproduce the temporal sequence of 3D spatial TMP distributions, we adopt a novel sequencetosequence variational autoencoder (VAE) [2]
with cascaded long shortterm memory (LSTM) networks. This VAE is trained on a large database of simulated TMP dynamics originating from various myocardial locations and with a wide range of local tissue properties. Second, once trained, the VAE decoder describes the likelihood of the TMP conditioned on a lowdimensional set of generative factors, while the encoder learns the posterior distributions of the generative factors conditioned on the training data. We utilize these two components within the Bayesian inference, and present a variation of the expectationmaximization (EM) algorithm to jointly estimate the generative factors and transmural TMP signals from observed ECG data. In a set of synthetic and realdata experiments, we demonstrate that the presented method is able to improve the accuracy of transmural EP imaging in comparison to statistical inference either constrained by a conventional physiological model
[12] or without physiological constraints.2 Generative Modeling of TMP via Sequential VAE
To learn to generate the spatiotemporal TMP sequences, we use a sequential variation of VAE [8] based on the use of LSTM networks [2].
VAE Architecture: The architecture of the sequential VAE is summarized in the red block in Fig. 1
. Both the encoder and the decoder consists of two layers of LSTM, where the second layer includes separate mean and variance networks. The spatial dimension decreases from the original TMP signal
U to the latent representation Z, while the temporal relationship is modeled by the LSTMs. Note that while the random variables in a standard VAE are vectors, a sequential VAE deals with matrices. By defining the conditional distribution of a matrix as the product of distributions over its columns, we obtained the likelihood distribution
and the variational posterior distribution as:(1) 
(2) 
where and are output from the mean and variance networks of the encoder parameterized by , and and are output from the mean and variance networks of the decoder parameterized by .
VAE Training: Training of the VAE is performed by maximizing the variational lower bound on the likelihood of the training data given as:
(3) 
where is an isotropic Gaussian prior. The calculation of the KL divergence and cross entropy loss for the presented sequential architecture is carried out in a manner similar to that described in [8].
The training data is generated by the AlievPanfilov (AP) model [1], simulating spatiotemporal TMP sequences originated from different ventricular locations with different tissue properties.
3 Transmural EP Imaging
The biophysical relationship between cardiac TMP, and bodysurface ECG, can be described by a a linear measurement model: , where is specific to the hearttorso model of an individual. To estimate U from Y is severely illposed and requires the regularization from additional knowledge about U.
3.0.1 Probabilistic Modeling of the Inverse Problem:
We formulate the inverse problem in the form of statistical inference. We define the likelihood distribution of Y given U by assuming zeromean measurement errors with variance :
(4) 
To incorporate physiological knowledge about U, we model its prior distribution conditioned on Z using the VAE decoder with trained parameter :
(5) 
To further utilize the knowledge about the generative factor Z learned by the VAE from a large training dataset, we also utilize the VAEencoded marginal posterior distribution of Z as its prior distribution in Bayesian inference. In specific, we approximate samples from this marginalized distribution to be Gaussian:
(6) 
With this, we complete the statistical formulation of our problem. Our goal is to estimate the joint posterior distributions
3.0.2 Inference:
Due to the presence of a deep neural network, the posterior
is analytically intractable. To address this issue, we note that conditioned on Z, the distribution of U is Gaussian in each column; thus, is analytically available. We leverage this fact and employ a variant of the expectation maximization (EM) algorithm to obtain the maximum a posteriori (MAP) estimate of Z along with the posterior distribution of U given the MAP estimate of Z .Estep: Conditioned on an estimated value of Z (say ), we calculate , with the covariance and mean of the column of U as:
(7) 
where , and and are the column output of the VAE decoder network when is input to it.
Mstep: Given , we update Z by maximizing
(8) 
Realizing that a complete optimization of with respect to Z
would be expensive, we instead take a few gradient descent steps towards the optimum. The gradient of the second term is analytically available. The gradient of the first term is calculated by backpropagation through the decoder network.
The EM steps iterate until convergence, at which we obtain both the MAP value of Z and the posterior distribution of U conditioned on Z and Y.
4 Results
4.0.1 Synthetic Experiments:
Synthetic experiments are carried out on two imagederived human hearttorso models. On each heart, the VAE is trained using around 850 simulated TMP signals considering approximately 50 different origins of ventricular activation in combination with 17 different tissue property configurations. As an initial study, here we focus on tissue properties representing local regions of myocardial scars with varying sizes and locations.
The presented method incorporating the trained VAE model is then tested on simulated 120lead ECG data from three different settings, each with 20 experiments. The three settings include 1) presence of myocardial scar not included in training data, 2) origin of ventricular activation different from those used in training, and 3) both myocardial scar and activation origin not seen in training. In all experiments, the performance of the presented method is compared to 0order Tikhonov regularization with temporal constraint (Greensite method) [5] and conventional EP model constrained inference with fixed parameters [12].
The reconstruction accuracy is measured with three metrics: 1) normalized RMSE given by the ratio of Frobenius norm of the error matrix to that of the truth TMP matrix, 2) Euclidean distance between the reconstructed and true origins of ventricular activation, and 3) Dice coefficient of the reconstructed and true regions of scar as =2/(). In the two physiologically constrained methods, region of scar is defined based on absence or delay of activation and shortening of action potential duration; in Greensite method, since the reconstructed signal no longer preserves the temporal shape of TMP, the region of scar is defined based on the peak amplitude of the signal.
Computational cost: Training of the VAE takes approximately 40 hours on a 4 GB Nvidia Quadro P1000 GPU. Generation of training data for each heart takes about 7 hours and inference around 30 minutes on Quadcore CPU.
TMP generation: Fig. 2 shows examples of local TMP signals generated by the trained VAE decoder against TMP signals simulated by the AP model [1]. Note that, when generating from a isotropic Gaussian (Fig. 2 right), noisy rather than meaningful TMP signals may also be generated. In comparison, when sampling from the approximated posterior distribution of Z as described in equation (6), the generated signals closely resemble the simulated TMP signals.



Imaging TMP from various origins: Fig. 3 shows a snapshot from the early stage of ventricular activation reconstructed by the three methods in comparison to the ground truth. Since the EP model constrained approach assumes general sinusrhythm activation, it introduces model error that incorrectly dominates the results. The simple Greensite method, free from erroneous model assumption, actually does a better job in comparison. By adapting model generative factors to the data, the presented method demonstrates a significantly improved ability to reconstruct TMP sequence resulting from unknown origins.
Imaging TMP at the presence of myocardial scar: Fig. 4 shows the spatial distribution of scar tissue obtained by the three different methods, along with temporal TMP signals reconstructed in healthy and scar regions, in comparison to the ground truth. Without prior physiological knowledge, the Greensite method is not able to preserve the temporal TMP shape, resulting in high RMSE error as shown in Table 1. By thresholding the maximum amplitude of the reconstructed signals, the identified region of scar has high false positives and resembles poorly with the ground truth. The EP model constrained approach does a better job in retaining the temporal TMP shape. However, without prior knowledge about the scar, the model error again affects the accuracy of TMP reconstruction, especially in the early stage of activation when a smaller amount of ECG data is available for correcting the model error. The presented method, in comparison, is able to recognize the presence of scar tissue, adapting the physiological constraint for improved TMP reconstructions and scar identifications.
Summary:
Table 1 summarizes the quantitative comparison of the three methods tested in the three settings as described earlier. Although the test cases were not seen by the VAE during training, the proposed method shows a significant improvement in inverse reconstruction (paired ttest, p
0.001) when compared with the other two methods in all settings and metrics except with Euclidean distance using Greensite method, where improvement is only marginal. It shows the importance of physiological knowledge and its adaptation to observed data during modelconstrained inference.4.0.2 Real data Experiments:
Two case studies are performed on realdata from patients who underwent catheter ablation due to scarrelated ventricular arrhythmia. Spatiotemporal TMP is reconstructed from 120lead ECG data using the presented method and the EP model constrained method. In Fig. 5, scar regions (red regions with low voltage) identified from the reconstructed TMP are compared with scar regions (red regions) in the invivo bipolar voltage data. In both cases, while the scar tissue identified by two methods are generally in similar locations, the presented method shows less false positives and higher qualitative consistency with bipolar voltage maps.
5 Discussion and Conclusions:
To our knowledge, this is the first work that integrates a generative network learned from numerous examples into a statistical inference framework to allow the adaptation of prior physiological knowledge via a small number of generative factors. The results show the ability of this concept to improve modelconstrained inference. Since the present formulation is in a personalized setting, we intend to extend this architecture to learn a geometryinvariant generative model that can be trained on multiple heart models and applied on a new subject.
Acknowledgement
This work is supported by the National Science Foundation under CAREER Award ACI1350374.
References
 [1] Aliev, R.R., Panfilov, A.V.: A simple twovariable model of cardiac excitation. Chaos, Solitons & Fractals 7(3), 293–301 (1996)
 [2] Bowman, S.R., Vilnis, L., Vinyals, O., Dai, A.M., Jozefowicz, R., Bengio, S.: Generating sentences from a continuous space. arXiv preprint arXiv:1511.06349 (2015)
 [3] Ghimire, S., Sapp, J.L., Horacek, M., Wang, L.: A variational approach to sparse model error estimation in cardiac electrophysiological imaging. In: International Conference on MICCAI. pp. 745–753. Springer (2017)
 [4] Ghodrati, A., Brooks, D.H., Tadmor, G., MacLeod, R.S.: Wavefrontbased models for inverse electrocardiography. IEEE TBME 53(9), 1821–1831 (2006)
 [5] Greensite, F., Huiskamp, G.: An improved method for estimating epicardial potentials from the body surface. IEEE TBME 45(1), 98–104 (1998)
 [6] Gulrajani, R.M.: The forward and inverse problems of electrocardiography. IEEE Engineering in Medicine and Biology Magazine 17(5), 84–101 (1998)
 [7] He, B., Li, G., Zhang, X.: Noninvasive imaging of cardiac transmembrane potentials within threedimensional myocardium by means of a realistic geometry anisotropic heart model. IEEE TBME 50(10), 1190–1202 (2003)
 [8] Kingma, D.P., Welling, M.: Autoencoding variational bayes. arXiv preprint arXiv:1312.6114 (2013)
 [9] Plonsey, R., Barr, R.C.: Bioelectricity: a quantitative approach. Springer Science & Business Media (2007)
 [10] Pullan, A., Cheng, L., Nash, M., Bradley, C., Paterson, D.: Noninvasive electrical imaging of the heart: theory and model development. Annals of biomedical engineering 29(10), 817–836 (2001)
 [11] Van Dam, P.M., Oostendorp, T.F., Linnenbank, A.C., Van Oosterom, A.: Noninvasive imaging of cardiac activation and recovery. Annals of biomedical engineering 37(9), 1739–1756 (2009)
 [12] Wang, L., Zhang, H., Wong, K.C., Liu, H., Shi, P.: Physiologicalmodelconstrained noninvasive reconstruction of volumetric myocardial transmembrane potentials. IEEE Transactions on Biomedical Engineering 57(2), 296–315 (2010)
Comments
There are no comments yet.