Deep learning analysis of coronary arteries in cardiac CT angiography for detection of patients requiring invasive coronary angiography

TL;DR

This study introduces a deep learning-based method for non-invasively identifying patients with coronary artery stenosis requiring invasive angiography, using cardiac CT scans and autoencoders to analyze coronary arteries.

Contribution

It presents a novel deep unsupervised analysis approach combining autoencoders and SVM classification for detecting significant coronary stenosis from CCTA images.

Findings

Achieved an AUC of 0.81 at artery level and 0.87 at patient level.

Demonstrated feasibility of automatic non-invasive detection of patients needing ICA.

Potential to reduce unnecessary invasive procedures.

Abstract

In patients with obstructive coronary artery disease, the functional significance of a coronary artery stenosis needs to be determined to guide treatment. This is typically established through fractional flow reserve (FFR) measurement, performed during invasive coronary angiography (ICA). We present a method for automatic and non-invasive detection of patients requiring ICA, employing deep unsupervised analysis of complete coronary arteries in cardiac CT angiography (CCTA) images. We retrospectively collected CCTA scans of 187 patients, 137 of them underwent invasive FFR measurement in 192 different coronary arteries. These FFR measurements served as a reference standard for the functional significance of the coronary stenosis. The centerlines of the coronary arteries were extracted and used to reconstruct straightened multi-planar reformatted (MPR) volumes. To automatically identify arteries with functionally significant stenosis that require ICA, each MPR volume was encoded into a fixed number of encodings using two disjoint 3D and 1D convolutional autoencoders performing spatial and sequential encodings, respectively. Thereafter, these encodings were employed to classify arteries using a support vector machine classifier. The detection of coronary arteries requiring invasive evaluation, evaluated using repeated cross-validation experiments, resulted in an area under the receiver operating characteristic curve of $0.81 \pm 0.02$ on the artery-level, and $0.87 \pm 0.02$ on the patient-level. The results demonstrate the feasibility of automatic non-invasive detection of patients that require ICA and possibly subsequent coronary artery intervention. This could potentially reduce the number of patients that unnecessarily undergo ICA.

Tables

Table 1. Table I: Average reconstruction MAPE, (within lumen HU range), total size of the final encoding used for classification and the achieved AUC for artery-level classification across a range of different input sizes of first CAE (CAE 1), different encoding sizes of first and second CAE (CAE 2) or different encoding strategies: PCA: using two consecutive principle component analyses; 1D: using a 1D-CAE (as CAE 2); 2D: using a 2D-CAE (as CAE 2); Global: using a single 3D-VCAE applied to the complete artery. Please note that the proposed configuration (input size of 40x40x5 with 16/64 encoding sizes) is listed multiple times for easy comparison. ∗ Number of principle components.

Input size of CAE 1	Encoding size of CAE 1	Encoding size of CAE 2	Encoding strategy	Reconstruction MAPE (within lumen)	Total encoding size	Average AUC
40x40x5	8	64	1D	$32.7\pm 12.6$	512	$0.74\pm 0.03$
40x40x5	16	64	1D	$29.1\pm 9.8$	1024	$0.81\pm 0.02$
40x40x5	32	64	1D	$17.7\pm 5.2$	2048	$0.63\pm 0.02$
40x40x5	${16}^{\ast }$	${64}^{\ast }$	PCA	$51.4\pm 13.0$	1024	$0.51\pm 0.03$
40x40x5	16	64	1D	$29.1\pm 9.8$	1024	$0.81\pm 0.02$
40x40x5	16	1024	2D	$40.7\pm 13.6$	1024	$0.68\pm 0.03$
40x40x800	1024	-	Global	$79.0\pm 18.1$	1024	$0.52\pm 0.01$
40x40x5	16	32	1D	$29.8\pm 9.6$	512	$0.78\pm 0.04$
40x40x5	16	64	1D	$29.1\pm 9.8$	1024	$0.81\pm 0.02$
40x40x5	16	128	1D	$23.6\pm 8.4$	2048	$0.61\pm 0.03$
40x40x5	16	64	1D	$29.1\pm 9.8$	1024	$0.81\pm 0.02$
40x40x13	16	64	1D	$28.5\pm 8.8$	1024	$0.80\pm 0.03$
40x40x23	16	64	1D	$32.1\pm 10.5$	1024	$0.72\pm 0.04$

Table 2. Table II: Average diagnostic accuracy for the detection of arteries and patients requiring ICA on the artery- and patient-levels shown in four different subgroups corresponding to four ranges of FFR measurements. N indicates the number of arteries and patients in each subgroup.

FFR Range	N Arteries	Artery-level Accuracy	N Patients	Patient-level Accuracy
$FFR\le 0.7$	32	0.66	26	0.70
	52	0.75	41	0.76
	73	0.79	52	0.78
$FFR\ge 0.9$	35	0.73	18	0.80

Table 3. Table III: Average diagnostic accuracy for the detection of arteries requiring ICA shown in three subgroups, corresponding to the three main coronary arteries: Left circumflex artery (LCX), right coronary artery (RCA) and left anterior descending (LAD). N indicates the number of arteries in the data set.

Artery	N	Accuracy
LCX	52	0.65
RCA	36	0.71
LAD	104	0.80

Table 4. Table IV: Performance comparison with previous work. Table lists number of evaluated patients and arteries, achieved accuracy and the area under the ROC curve (AUC) per-patient and per-artery for classification according to FFR as reported in the original studies. Please note that these methods use different FFR thresholds and perform different analyses: either analyzing of the blood flow in the coronary arteries ( Flow ), detecting ischemic changes directly in LV myocardium ( Myo. ), or, as proposed in this work; classifying coronary arteries with features extracted by CAEs.

					Per-patient		Per-artery
	Study	Patients	Arteries	$FFR\le$	Accuracy	AUC	Accuracy	AUC	Requirements
Artery	Nørgaard et al.[18]	254	484	0.8	0.81	0.90	0.86	0.93	Artery lumen segmentation
	Coenen at al.[44]	106	189	0.8	-	-	0.74	-	Artery lumen segmentation
	Coenen at al.[20]	303	525	0.8	0.71	-	0.78	0.84	Artery lumen segmentation
Myo.	Zreik et al.[12]	126	-	0.8	0.64	0.66	-	-	LV myocardium segmentation
	Han et al.[25]	252	407	0.8	0.63	-	0.57	-	LV myocardium segmentation
	Proposed	137	192	0.9	0.80	0.87	0.78	0.81	Artery centerline tracking

Full text

††publicationid: pubid: Accepted in IEEE Transactions on Medical Imaging, 2019††publicationid: pubid:

Accepted in IEEE Transactions on Medical Imaging, 2019

I Introduction

Obstructive coronary artery disease (CAD) is the most common type of cardiovascular disease [1]. Obstructive CAD develops when atherosclerotic plaque builds up in the wall of the coronary arteries, narrowing the coronary artery lumen [2]. This is defined as coronary stenosis, which can potentially limit blood supply to the myocardium, and could lead to ischemia and irreversible damage [3]. Only functionally significant stenoses, i.e. those stenoses which significantly limit blood flow, need to be invasively treated in order to reduce CAD morbidity [3, 4, 5, 6]. Contrarily, invasively treating a functionally non-significant stenosis may lead to harmful output [5, 7]. Therefore, it is crucial to assess the functional significance of a coronary stenosis to guide treatment.

Cardiac CT angiography (CCTA) is typically used to noninvasively identify patients with suspected CAD and visually detect coronary artery stenosis [8]. Although CCTA has high sensitivity in determining the functional significance of the stenosis, its specificity for this task is low [9, 10, 11]. Therefore, to determine whether a coronary artery stenosis is functionally significant, patients with obstructive CAD typically undergo invasive coronary angiography (ICA) to measure the fractional flow reserve (FFR) in the coronary arteries. FFR is currently the reference standard for establishing the functional significance of a coronary stenosis and it is used to guide treatment [3, 4]. However, because of the low specificity of CCTA, up to 50% of patients undergo invasive FFR measurement unnecessarily [11]. To reduce the number of unnecessary invasive procedures, noninvasive determination of the functional significance of stenoses based on CCTA images has been intensively investigated. Several automatic methods for determination the functional significance of coronary artery stenosis in CCTA have been proposed [12]. These methods can be divided into those that simulate and analyze blood flow in the coronary arteries [13, 14, 15, 16], and those that analyze and characterize the left ventricle (LV) myocardium [12, 17].

Methods that simulate and analyze the blood flow in the coronary arteries in CCTA images estimate FFR values along the coronary artery, which can be used to determine the functional significance of coronary artery stenosis. Taylor et al. [13] were the first to propose noninvasive flow-based FFR estimation from CCTA images, which was later validated in multiple clinical studies[18, 19]. To determine FFR values along the coronary artery, computational fluid dynamics, coupled with assumptions of physiological boundary conditions, were used. Itu et al. [14] also presented a method to estimate FFR in the coronary artery tree in CCTA images by simulating blood flow. This method uses a parametric lumped heart model, while modeling the patient-specific hemodynamics in both healthy and diseased coronary arteries. Nickisch et al. [15] determined FFR values along the coronary artery by simulating blood flow and pressure along the coronary artery arteries using an electrical patient-specific parametric lumped model. Moreover, Itu et al. [16] presented a machine-learning-based model for estimating FFR along the coronary artery. The model is trained on a large number of synthetically generated coronary anatomies, where the target values are computed using a blood flow-based model [14]. This method was further evaluated in [20]. While these techniques [13, 14, 15, 16] achieved high accuracy, they are remarkably dependent on the accuracy of coronary artery lumen segmentation [21]. Manual annotation of the coronary artery lumen is a time consuming and a complex task, where commercially available automatic software tools typically require substantial manual interaction and correction, especially in CCTA scans with excessive atherosclerotic calcifications or imaging artefacts due to stents and cardiac motion [22].

Recently, methods that do not model the blood flow in the coronary arteries but employ characteristics extracted from the myocardium in CCTA scans, have shown to be feasible. Our recent work [12, 23] presented a deep learning approach to automatically identify patients with a functionally significant coronary artery stenosis using analysis of the LV myocardium in CCTA. The method first characterizes the LV myocardium using a convolutional autoencoder (CAE). Thereafter, using the extracted characteristics, patients are classified according to the presence of functionally significant stenosis using an SVM classifier. Previously, Xiong et al. [17] presented a machine learning based approach to detect patients with anatomically significant stenosis using characteristics of the LV myocardium derived from a CCTA scan. In this method, the LV myocardium is aligned with the standard 17-segments model [24] to relate each myocardial segment to its perfusing coronary artery. Then, hand-crafted features, describing each myocardial segment, are extracted and used for supervised classification of patients according to the presence of anatomical significant stenosis. Thereafter, Han et al. [25] employed the technique described in [17] to detect patients with functionally significant stenosis, as defined by the invasively measured FFR. Although these new methods [12, 17] have presented promising results without the need for accurate coronary artery lumen segmentation, they still need to be validated in large and diverse patients cohorts.

Moreover, in our recent work [26], we have analyzed the coronary arteries employing a recurrent convolutional neural network (RCNN) for detecting and classifying the anatomical significance of the coronary artery stenosis. The RCNN employs a 3D convolutional neural network to extract local features along the coronary artery. Subsequently, a recurrent neural network aggregates the features to perform the classification tasks. However, such an approach cannot be directly employed for the detection of the functional significance of a coronary stenosis for two reasons. First, such RCNN only performs a local analysis of the artery, were the complete artery is not taken into account. Second, to train such an RCNN, local reference labels are required. Such a requirement is not practical in the case of the functional significance of a coronary stenosis, where FFR is used as the reference and is usually provided on the artery level only.

Here, we present a method to automatically and non-invasively identify coronary arteries and patients requiring further invasive evaluation, i.e. ICA, as determined by the invasively measured FFR. Blood flow in the coronary artery may be affected by multiple coronary artery stenoses and arterial plaques [3, 27]. Therefore, to classify an artery according to the functional significance of the coronary artery stenosis, local analysis of a single stenosis may be insufficient. Hence analysis of the complete artery should be performed. Moreover, in clinical practice, usually a single, i.e. lowest, FFR value per coronary artery is reported. Consequently, employing supervised machine learning methods to directly analyze a whole volume of an artery (e.g. with 3D-CNN or RCNN [26]) to detect the functional significance of each stenosis or estimating the invasively measured FFR values at every point along the coronary artery is unfeasible. Therefore, in the proposed work, a complete artery is analyzed in an unsupervised manner to extract lower-dimensional encoding, and thereafter to determine the presence of abnormal FFR. First, using the extracted coronary artery centerline [28], the straightened 3D multi-planar reformatted (MPR) volume is reconstructed. Then, an MPR volume of a complete artery is characterized with a fixed number of encodings using convolutional autoencoders (CAEs) [29, 30, 31], which serve as unsupervised feature extractors. As MPR volumes of complete coronary arteries have large volumetric sizes and variable lengths and shapes, a single traditional CAE cannot be successfully and directly applied to efficiently encode a complete artery. Therefore, in the here proposed work, two disjoint CAEs are employed. The first CAE performs spatial encoding of local sub-volumes along the artery. Then, a second CAE encodes the output of the first CAE - which depends on the artery length - into a fixed-length encoding. Finally, a support vector machine (SVM) [32] classifies arteries based on these encodings according to presence of functionally significant stenosis, as defined by the invasively measured FFR. The proposed approach is illustrated in Fig. 1. Our contributions are twofold. Firstly, we propose to jointly employ two disjoint CAEs that perform spatial and sequential encoding of large volumes with varying lengths. Secondly, in contrast to previous methods that detect the presence of functionally significant stenosis or determine FFR values non-invasively, our method does not require accurate and difficult to obtain segmentation of the coronary artery lumen or LV myocardium. Instead, it only requires the coronary artery centerline, which can be obtained automatically or semi-automatically [33].

The remainder of the manuscript is organized as follows. Section II describes the data and reference standard. Section III describes the method. Section IV reports our experimental results, which are then discussed in Section V.

II Data

II-A Patient and Image Data

This study includes retrospectively collected CCTA scans of 187 patients (age: $58.6\pm 8.7$ years, 145 males) acquired between 2012 and 2016. The Institutional Ethical Review Board waived the need for informed consent.

All CCTA scans were acquired using an ECG-triggered step-and-shoot protocol on a 256-detector row scanner (Philips Brilliance iCT, Philips Medical, Best, The Netherlands). A tube voltage of 120 kVp and tube current between 210 and 300 mAs were used. For patients $\leq 80$ kg contrast medium was injected using a flow rate of 6 mL/s for a total of 70 mL iopromide (Ultravist 300 mg I/mL, Bayer Healthcare, Berlin, Germany), followed by a 50 mL mixed contrast medium and saline (50:50) flush, and next a 30 mL saline flush. For patients $>80$ kg the flow rate was 6.7 mL/s and the volumes of the boluses were 80, 67 and 40 mL, respectively. Images were reconstructed to an in-plane resolution ranging from 0.38 to 0.56 mm, and 0.9 mm thick slices with 0.45 mm spacing.

In each CCTA scan, coronary arteries were tracked and their centerlines were extracted using the method previously described by Wolterink et al. [28]. The method tracks the visible coronary arteries, where the arterial centerlines are extracted between the ostia and the most distal visible locations. Using the extracted centerlines, a 3D straightened MPR volumes with 0.3 $mm^{3}$ isotropic resolution were reconstructed for all coronary arteries and used for further analysis. Note that we define an artery as the vessel starting from the ostium until the most distal location visible in the CCTA.

II-B FFR Measurements

Out of the 187 patients, 137 patients suspected of obstructive CAD underwent invasive FFR measurements ($0.81\pm 0.10$, interquartile range: 0.74-0.89), up to one year after the acquisition of the CCTA scan. In these patients, FFR was measured in 192 different arteries. FFR was recorded with a coronary pressure guidewire (Certus Pressure Wire, St. Jude Medical, St. Paul, Minnesota) at maximal hyperemia conditions. Maximal hyperemia was induced by administration of intravenous adenosine (at a rate of 140 $\mu$g/kg per minute) through a central vein. The FFR wire was placed at the most distal part possible in the target artery. Using manual pullback, a single minimal FFR value was assessed and recorded for each artery.

III Methods

Blood flow in the coronary artery may be affected by a single or multiple coronary artery stenoses [3, 27], located anywhere along the coronary artery; starting from the ostium until the most distal location visible in the CCTA. Therefore, to classify an artery according to the functional significance of a stenosis, local analysis of a single stenosis may be insufficient, but the analysis of the complete artery is needed. Moreover, in clinical practice, determining invasive FFR values for each voxel within the artery lumen, or recording an FFR value for each point on the coronary artery centerline, is impractical and typically not performed. Instead, the minimal single FFR value per coronary artery is recorded, resulting in a single reference label per artery. Hence, given the sparsity of reference labels along the artery, the large input dimensions and the limited dataset size, employing a supervised end-to-end machine learning methods (e.g. with 3D-CNN or RCNN [26]) to directly detect the functional significance of each stenosis or estimating FFR values at every point along the coronary artery would be prone to overfitting. Therefore, in the proposed work, an MPR of a complete artery is analyzed to determine the presence of abnormal FFR. First, to extract robust features of complete arteries, MPR volumes are characterized by a fixed number of encodings using convolutional autoencoders (CAEs) [29, 30, 31], regardless of the artery length. Then, the extracted encodings are used as input to an SVM classifier that determines whether the artery needs further invasive evaluation, in the form of ICA, to establish the need of intervention.

III-A Encoding the artery

The main purpose of a convolutional autoencoder (CAE) is to extract robust compact features from unlabeled data, while removing input redundancies and preserving essential aspects of the data [29, 30, 34]. A CAE consists of two main parts, an encoder and a decoder [29, 30]. The encoder compresses the data to a lower dimensional latent space by convolutions and down-sampling. The decoder expands the compressed form to reconstruct the input data by deconvolutions and upsampling. A CAE is trained to minimize a difference loss between the encoder input and decoder output. This ensures that the encodings contain sufficient information to reconstruct inputs with low error [30]. Once the CAE is trained, the decoder is removed and the encoder is used to generate encodings for unseen data.

Coronary arteries are complex anatomical 3D structures, with varying lengths and anomalies across patients [35]. The resolution of modern CT scanners is high and a large number of voxels (millions) is contained in an MPR volume of a single artery. Therefore, following the straightforward approach of training a single CAE, applied directly to the complete artery volume without a large reconstruction error, is infeasible. Therefore, in this work, we propose a two-stage encoding approach to encode a complete MPR volume of the coronary artery, regardless of its length. Fig. 2 illustrates the proposed encoding flow. First, a 3D variational convolutional autoencoder (3D-VCAE) is applied to local sub-volumes extracted from the MPR along the artery centerline. As the 3D-VCAE is only applied to small input volumes, the number of its trainable parameters is relatively low. The 3D-VCAE encodes each sub-volume into a set of small number of encodings. When applied to all sequential sub-volumes along the artery, the result is a feature map of the same height as the number of encodings and the same length as the artery length. This feature map is then represented as a set of individual 1D sequences of encodings. Each sequence contains an individual encoding out of the set of encodings, running along the artery (colored signals in Fig. 2). This allows the analysis of complete arteries with varying length by a 1D convolutional autoencoder (1D-CAE), with low number of trainable parameters which decreases the chance of overfitting. Hence, the 1D-CAE encodes the varying length sequences of encodings further into a fixed number of encodings, that represent the complete artery, regardless of its length.

III-A1 Spatial encoding with 3D variational convolutional autoencoder

VAEs are generative models, which approximate data generating distributions [31]. Through approximation and compression, the resulting models have been shown to capture the underlying data manifold; a constrained, smooth, continuous, lower dimensional latent (feature) space where data is distributed [36, 37]. Having in mind possible reconstruction errors of the encoding sequences by the second CAE, a variational CAE is chosen as the first CAE for the ability of its decoder to handle small variations in the encodings [36]. Inspired by these advantageous properties of the latent space, a VCAE is employed to compress and encode local volumes along the artery. To capture local volumetric characteristics of the artery, the input to the 3D-VCAE is set to a volume of 40x40x5 voxels, centered around a coronary artery centerline point. The size of the input is chosen so that it contains the whole arterial lumen and the vicinity of the artery [2]. The output of the encoder in the 3D-VCAE is set to 16; i.e. an encoding in a $R^{16}$ latent space. The dimension of the input volume and encoding size are determined in preliminary experiments to balance between the compactness (i.e. size of the encoding) and the expressiveness of the encodings (i.e. reconstruction error). Table I lists these findings. To encode the complete artery, overlapping volumes with stride of 1 are extracted and encoded with 3D-VCAE (Fig. 2). This results in 16xL encodings, where L is the length of the artery. To reconstruct a complete artery, the middle slices of each overlapping reconstructed volume are used. The 3D-VCAE architecture used in this work is shown in Fig. 3(a). In the 3D-VCAE , batch normalization [38] layers and rectified linear units (ReLUs) are used after all convolution layers except the encoder and decoder output layers.

III-A2 Sequential encoding with 1D convolutional autoencoder

When representing the coronary artery to determine the functionally significant stenosis according to FFR, characteristics along the artery, starting from the ostium to the most distal part of the artery, need to be taken into account [13, 14, 15]. Therefore, to analyze the complete artery at once, the local encodings extracted previously by the 3D-VCAE along the length of the artery need to be merged. To accomplish this, the feature map, consisting of L sets of 16 values of encoding generated by 3D-VAE at each coronary artery center point, is represented as L 1D sequences. As in the 3D-VCAE design, the size of the 1D-CAE encoding is determined in preliminary experiments to balance between the compactness (i.e. size of the encoding), the expressiveness of the encodings (i.e. reconstruction error) and the classification performance (AUC). Table I lists these findings. Each sequence consists of 1xL values, where L represents the length of the artery, i.e. number of coronary artery centerline points. To encode arteries with different lengths, sequences of encodings of short arteries were padded into a maximum length of 800, which corresponds to the number of centerline points in the longest artery in the dataset. This representation leads each sequence to represent a specific member of the encoding in the $R^{16}$ latent space along the artery (colored 1D signals in Fig. 2). This, consequently, allows us to apply a 1D-CAE to each of the 16 sequences separately. The weights of the 16 1D-CAEs are shared, where each 1D-CAE encodes one of the 16 sequences into an encoding of a second latent space of 64 dimensions ($R^{64}$). This results in 1024 (16x64) features that represent the complete artery. The 1D-CAE architecture used in this work is shown in Fig. 3(b). In the 1D-CAE, the exponential linear units (ELUs) are used after all convolutions layers except the encoder and decoder output layers.

III-B Classification of arteries and patients

Based on the extracted encodings from the encoding stage, arteries are classified according to the need of further invasive evaluation, in the form of ICA. This was defined by the invasively measured FFR. As the standard deviation of the differences in repeated FFR measurements can reach up to $5\%$ [39, 40, 41, 42], especially in the so called ”gray-zone” [41], where the measured FFR is between 0.75 and 0.85, in our experiments, the threshold on FFR value was set to 0.9. This results in a positive class with $FFR\leq 0.9$ representing arteries requiring ICA to establish the need of intervention, and a negative class with $FFR>0.9$ representing absence of functionally significant stenosis, where ICA is not necessary. The classification is performed using an SVM classifier with a linear kernel and an $L_{1}$ regularization. For each classified artery, the continuous output of the trained SVM is used to assign a predicted class.

As patients with suspected obstructive CAD undergo ICA to measure the FFR in all diseased coronary arteries, in this work, classification of patients is also performed. To classify patients, the highest output value of all classified arteries in a patient is used to assign a predicted class to the patient. The minimal FFR across the arteries of a patient is taken as a reference.

Classification performance is evaluated using a receiver operating characteristic (ROC) curve and the corresponding area under the ROC curve (AUC).

IV Experiments and Results

IV-A Encoding the artery

To train the 3D-VCAE and the 1D-CAE, a set of CCTA images of 50 patients, who did not undergo ICA and hence had no FFR measurements, were used. From these, 38 CCTA images were randomly selected for training, and the remaining 12 images were used for validation. In both sets, MPR volumes of the extracted arteries were used to train and validate the autoencoders. Both autoencoders’ hyperparameters were determined in preliminary experiments using the validation set. Please note that no data augmentation was performed for training the autoencoders. As a baseline reconstruction, principal component analysis (PCA) was employed in a similar way as the two disjoint autoencoders: A first PCA was applied to all 40x40x5 voxels volumes along the artery to reduce each into 16 principal components. Then, a second PCA was applied to all outputs of the first PCA to reduce the dimensionality of the outputs into 1024 components. Table I lists the findings.

To train and validate the 3D-VCAE, 40x40x5 voxels volumes were randomly extracted along centerlines of arteries in the training and validation sets, respectively. Mini-batches of 32 volumes were used to minimize the loss function with Adam optimizer [43] with a learning rate $0.001$. The mean squared error between the input and the reconstructed volumes, and the Kullback-Leibler (KL) divergence with the reparameterization trick [31] were employed as a loss function for the variational autoencoder. L2 regularization was used with $\gamma=0.001$ for all layers. Training was performed until convergence. Fig. 4(a)-(b) show an example of a complete artery which was encoded and reconstructed with the trained 3D-VCAE. This was performed by extracting, encoding and reconstructing input volumes around each point along the MPR centerline. Fig. 4(c) shows the absolute reconstruction error, i.e. the absolute difference between the input and the reconstructed artery. Fig. 4(d) shows all 16 sets of encodings, presented as continuous sequences running along the artery. These sequences are to be encoded further in a later stage using the 1D-CAE into a fixed number of encodings.

To train and validate the 1D-CAE, arteries with the corresponding sets of 16 encodings sequences, obtained by the 3D-VCAE, were randomly chosen from the training and validation sets, respectively. Sequences of encodings of short arteries were padded into a maximum length 800, which corresponds to the longest artery in the dataset. Mini-batches of 32 sets of encodings sequences were used to minimize the loss function with Adam optimizer with a learning rate $0.001$. The masked mean squared error was employed as a loss function, where padded values in the input sequences did not contribute to the loss value or its gradients and were therefore ignored. L2 regularization was used with $\gamma=0.001$ for all layers. Training was performed until convergence. Fig. 4(e) shows an example of 3 randomly chosen encoding sequences of a complete artery which were encoded and reconstructed with the trained 1D-CAE.

To demonstrate the effectiveness of the proposed combined two-stage encoding approach in preserving the original shape and appearance of the artery, both disjoint trained autoencoders were combined and tested on complete arteries. To accomplish this, an inference with four steps was performed. First, the encoder of the 3D-VCAE was applied to local volumes along the MPR volume of a complete artery, resulting in 16 sequences of encodings with L values each. Second, the encoder of the 1D-CAE encoded the sequences into a single encoding vector of 1024 values. Third, the decoder of the 1D-CAE decoded the encoding vector back to 16 encodings sequences. Last, the decoder of the 3D-VCAE reconstructed those reconstructed encodings sequences to the original MPR volume size. Fig. 4(a),(f),(g) show an example of a complete artery that was encoded with the combined strategy, reconstructed back to the original volume dimensions, and the corresponding reconstruction error. Fig. 5 compares the average mean absolute reconstruction percentage errors (MAPE) between the local reconstructions made by only the 3D-VCAE, the baseline PCA reconstruction approach and the combined approach reconstruction, across a range of CT Hounsfield units (HU). As the MAPE might be misleading around small image values or image values equal to zero, the range of intensity values characteristic for coronary artery lumen (250-450 HU) is highlighted. Fig. 4 and Fig. 5 demonstrate the high resemblance and the low reconstruction error between the results of the local and the combined approaches compared to the original volume.

IV-B Evaluation of alternative encoding strategies

To demonstrate that the proposed combined encoding strategy is advantageous, two additional encoding strategies were evaluated and compared to the proposed sequential disjoint autoencoders.

First, the most straightforward approach was evaluated, where a single autoencoder analyzes the complete MPR volume, encodes it to a fixed number of encodings (1024), and reconstructs it back to the input size. To handle arteries with different lengths, shorter arteries were padded to the maximal artery length in our dataset (800). Hence, the input of the autoencoder was defined as 40x40x800. The architecture of the evaluated VCAE is shown in Fig. 6.

Second, the 1D-CAE used in the combined encoding strategy (Fig. 2) was replaced by a 2D-CAE to jointly process the encoding sequences. While the proposed 1D-CAE encoded each sequence of encodings separately, the here evaluated 2D-CAE mutually encoded all sequences of encodings. This was performed by representing the encodings map as a 2D image and applying two-dimensional convolutional kernels. The architecture of the evaluated 2D-CAE is identical to the 1D-CAE (Fig. 3(b)), but the convolutions were performed by applying two-dimensional kernels of the same size (3x3). Additionally, dropout of 0.1 was applied between fully connected layers to avoid overfitting.

The two additional autoencoders were trained and validated in a similar manner as in the combined encoding strategy. In the case of the 2D-CAE, the trained decoder of the 3D-VCAE was used to reconstruct the MPR volume. Fig. 7 shows an example of an artery that was encoded and reconstructed using the two evaluated auto-encoding strategies, and compared with the reconstruction of the proposed combined encoding strategy. Fig. 5 also shows the average MAPE across a range of different HUs. Both figures demonstrate a clear advantage for the reconstruction of the combined strategy over the two additionally evaluated approaches.

IV-C Classification of arteries and patients

Classification of arteries was performed using the arteries’ encodings, extracted by the two disjoint autoencoders (Section IV-A), and an SVM classifier. In preliminary experiments, different classifiers (logistic regression, random forest) and various SVM configurations, including different regularizations and kernels types, were tested. The best performance was achieved using a linear $L_{1}$-regularized SVM. All 50 CCTA images of patients used in training and validation of the autoencoders were excluded. Thus, CCTA images of 137 patients and the corresponding reference FFR measurements in 192 different arteries were used for this analysis. To assess the performance and the robustness of the classification, 1000 stratified Monte-Carlo cross-validation experiments were performed. In each experiment, 10 arteries were used as a test set, and the remaining arteries were used as a training set. The assignment of arteries to test or training sets was random, however it insured that arteries from the same patient were included either in the training or the test set. Optimal SVM parameters were selected in every experiment using a grid search on the training set only.

The obtained results are shown in Fig. 8. On the artery-level, an average AUC of $0.81\pm 0.02$ was achieved, while on the patient-level, an average AUC of $0.87\pm 0.02$ was achieved. Table II lists the average diagnostic accuracy on the artery- and patient-levels across four different ranges of FFR measurements, and Table III lists the achieved performance in the three main coronary arteries. Moreover, implemented in Keras with TensorFlow, the runtime of encoding and classifying a single MPR volume of a complete coronary artery was on average 11 seconds, while using a single NVIDIA TITAN X (Pascal) GPU with an Intel Xeon machine with 256 GB RAM.

IV-D Evaluation of alternative classification strategies

To demonstrate the effectiveness of the proposed classification scheme, we have performed several additional classification experiments that can be divided into two categories.

First, the influence of the FFR threshold was investigated. As in Section III-B, arteries were classified into positive or negative class with a binary SVM classifier using the extracted encodings. Different thresholds on reference FFR values were applied, resulting in different class interpretations: When a threshold of 0.7 was applied, positive classes represented arteries in need for an invasive intervention without the need of establishing the FFR value first. When a 0.8 threshold was applied, positive class represented arteries with a functionally significant stenosis.

Second, the influence of regression vs. classification was investigated. In contrast to the former classification scenarios where retraining the SVM classifier was needed for each different FFR threshold, here, a single SVM regressor was trained to estimate continuous values of FFR (i.e. regression) and then an FFR threshold (0,7, 0.8 or 0.9) was applied to output the binary classes.

The ROC curves showing the results are given in Fig. 9. The results show that binary classifications outperform the corresponding regression experiments, regardless of FFR threshold. Moreover, when an FFR threshold of 0.7 was used, performance of the both classification and regression approaches were moderate, while the experiments using threshold of 0.8 on FFR values showed lowest areas under the ROC curves.

IV-E Comparison with other FFR classification methods

We compare our classification performance with the reported results of previous methods. These methods either analyzed the blood flow in the coronary arteries [18, 44, 20], requiring a highly accurate segmentation of the arterial lumen, or analyzed the LV myocardium [25, 12], requiring a segmentation of the LV myocardium. Table IV lists the results as originally reported. However, all of the compared methods were evaluated on different datasets that included different patients cohorts. Moreover, the compared methods employed an FFR cut-off value of 0.8 to define the functional significance of a stenosis, while in this study, a 0.9 cut-off point was used. Therefore, these results only indicate the differences in approaches and should not be directly compared.

V Discussion

A method for automatic and non-invasive identification of patients requiring evaluation with invasive coronary angiography has been presented. The method analyzes complete coronary arteries with two convolutional autoencoders that characterize the MPR volume of each artery with general robust features, and encode the complete artery into a fixed number of encodings to reduce the dimensions of the input. Then, these encodings are used with an SVM classifier to identify arteries with functionally significant stenosis in an supervised manner. Unlike previous methods that detect functionally significant stenosis by relying either on the coronary artery lumen segmentation [13, 14, 15, 16] or the left ventricle myocardium segmentation [12, 17], the proposed method requires only the coronary artery centerline as an input along with the CCTA scan. Artery centerline extraction is a simplified task compared to myocardium segmentation and to the arterial lumen segmentation, where the latter occasionally requires substantial manual interaction, especially in diseased population with heavily calcified arteries. In this work, to extract the coronary artery centerlines, we have employed our previously designed method for artery centerline extraction [28]. However, any other manual, semi-automatic or automatic method could be employed instead.

As the dimensions of MPR volumes of complete arteries are large, and the reference labels are only provided on the artery level, employing a straight-forward supervised 3D-CNN to detect the functional significance of a stenosis is far from feasible. Therefore, here, we have used unsupervised learning to characterize and encode each artery before employing a supervised classifier to detected abnormal FFR. To do so, two disjoint CAEs, that were applied sequentially, were employed. This is contrary to the more common approach of using a single CAE that encodes the complete artery volume at once. The results show that the learned encodings were able to represent the artery shape and appearance accurately, as was demonstrated qualitatively (Fig. 4) and quantitatively by the relatively small mean absolute error between the input and the reconstructed volumes (Fig. 5).

The output of the combined encoding strategy was examined and compared to the output of the 3D-VCAE on local volumes. Fig. 4 and Fig. 5 show that, in the range of CT values of the artery lumen (250-450 Hounsfield units), both the local and the combined encoders (with 1D-CAE) achieved satisfactory results, where the local approach was slightly advantageous. However, the lower error achieved with the local approach could be explained by the larger number of encodings used per artery compared to the combined approach. It can be noticed (Fig. 4(b) and (f)) that the reconstructions of both approaches preserved the shape and the morphology of the artery, while not being able to preserve the texture of neither the lumen nor the background. Although the lumen of the coronary artery was accurately reconstructed by the proposed encoding method, some small calcifications within the artery were entirely or partially lost (Fig. 4 and Fig. 7, respectively). This might be due to the low number of encodings used in the 3D-VCAE (16). Beside increasing the size of the encoding, future work could address this by modifying the loss function of the 3D-VCAE to penalize such errors in reconstruction, or by over-sampling such calcifications in the training process.

In the proposed combined encoding strategy, both autoencoders were disjoint during training and were combined only during inference, which could lead to error propagation. To overcome this, one alternative would be training both autoencoders simultaneously and end-to-end. However, in preliminary experiments, this was proven difficult, mainly due to hardware limitations. Another alternative would be training the 3D-VCAE separately, and then, using its trained decoder during training the 1D-CAE. This could be done by directly minimizing the mean squared error between the original and the reconstructed MPR volumes of the complete artery instead of between its original and reconstructed encodings sequences. Such training might potentially compensate for errors or prevent error propagation between the two disjoint training processes. Future work might address this.

Additional encoding strategies were performed and compared with the proposed one. Although other studies showed that a single 3D-CAE might be successfully used on large volumetric input [45], in the proposed work training a single 3D-VCAE (Fig. 6), applied directly to the complete artery volume without a large reconstruction error, was proven infeasible (Fig. 5 and Fig. 7). This could be due to a number of reasons: The very large number of trainable parameters ($\sim 65\times 10^{6}$) of the CAE, the high variability among the shapes and lengths of the arteries, the atypical aspect ratio of the input (40x40x800), or the lack of large set of training data. Moreover, treating all encoding sequences as a single 2D image and encoding this image with a 2D-CAE proved inferior when compared with the proposed 1D-CAE (Fig. 5 and Fig. 7). This might be explained by the lack of local spatial relations between the different encodings ($\mu_{0}-\mu_{15}$) in Fig. 2 at a given location along the coronary artery. These local spatial relations among the encodings motivate the use of 2D kernels in a typical 2D-CAE, when analyzing natural or medical images.

The artery was represented by multiple sequences of encodings, obtained after applying the 3D-VCAE to local sub-volumes along the artery. To further encode the artery to a fixed number of encodings, a 1D-CAE was applied to each sequence of encodings separately. As the fully connected layer in 1D-CAE (layer $e$ in Fig. 3(b)) expects a fixed number of inputs, the input sequences were padded to the maximal length of an artery in the dataset. A masked loss function was employed during training the 1D-CAE to minimize the effect of such padding on the reconstruction error. Despite this masking, the proposed padding could affect the classification performance as a function of the artery length. To enable the autoencoder to handle variable length sequences, without the need of padding it, a recurrent autoencoder could be employed [46]. In such a recurrent autoencoder, known as sequence-to-sequence autoencoders, a recurrent layer, with Gated Recurrent Units (GRUs) [47] or long short-term memory (LSTM) units [48], replaces the fully connected layer in the proposed 1D-CAE, to recursively process and encode a sequential varying length input. Future work might investigate such a recurrent autoencoder and its affect on the classification performance.

Our experiments show that moderate classification performance on both the artery- and patient-levels was achieved, while using only features derived in an unsupervised manner from a CCTA of coronary arteries (Fig. 8). These results show that the proposed approach could potentially lead to a reduction in the number of patients that unnecessarily undergo invasive coronary angiography. For example, as seen in Fig. 8(b), at the sensitivity of 80% or 90% in detecting patients requiring ICA, i.e. those having $FFR\leq 0.9$, unnecessary ICA could have been prevented in 76% or 53% of the negative patients, i.e. those having $FFR>0.9$, respectively. Moreover, we have compared the classification results across different ranges of FFR measurements (Table II). The comparison shows slightly higher accuracies for $FFR>0.8$, on both the artery- and patient-levels. The reason might be that arteries with $FFR>0.8$ typically contain less plaque and therefore were better characterized by the autoencoders. A comparison of the diagnostic accuracy in the three main coronary arteries (Table III) shows that the highest accuracy was achieved for the LAD. This might be due the small number of available training examples for the LCX and RCA.

Unlike this study, most previous methods that analyze blood flow [13, 14, 15, 16] for detection of functionally significant stenosis as determined by invasive FFR estimate continuous values of FFR along the coronary artery and localize the functionally significant stenoses using a threshold of 0.8 on the determined FFR [41]. However, in preliminary experiments for estimation of continuous FFR values (i.e. with regression) or applying such a threshold, the proposed method has not achieved satisfactory results (refer to Section IV-D and Fig. 9). The reason may be threefold. First, unlike other methods [13, 14, 15, 16], the here proposed method analyzes only a single coronary artery at once, while not taking into account the other arteries in the complete coronary artery tree. Analyzing the entire coronary tree might be crucial to differentiate between arteries with functionally significant or non-significant stenoses. Second, as the complete coronary artery is characterized as a whole using the CAEs, spatial information about a specific stenosis can not be retained. Third, the small number of encodings used in this work preserves the coarse shape and morphology of the analyzed artery (see Fig. 4(f)). However, fine morphology might be crucial for differentiating the functionally significant stenoses with FFR measurements in the most difficult range around the FFR of 0.8 [41]. Although these methods that analyze blood flow reported better results [13, 14, 15, 16], they are heavily dependent on the accuracy of coronary artery lumen segmentation [21]. Highly accurate lumen segmentation is extremely challenging task especially in patients with excessive atherosclerotic calcifications or imaging artefacts or stents [22]. As a result, these patients are typically not eligible for such analysis and are excluded[19, 49, 50]. In contrast, our method does not require lumen segmentation and therefore heavily diseased patients were not excluded. With a larger data set, future work may further investigate estimation of FFR on the continuous scale, possible performance enhancement when complete coronary artery tree is taken into analysis, and investigate different encoding approaches that may preserve fine morphology of the arteries.

To conclude, this study presented an automatic and non-invasive analysis of the coronary arteries in CCTA for detection of patients requiring invasive coronary angiography to establish the need of coronary intervention. The method is based on two disjoint convolutional autoencoders that characterize and encode volumes of complete coronary arteries into a set of encodings. Thereafter, a support vector machine classifier classifies arteries, employing these encodings, according to the presence of abnormal invasively measured FFR. The achieved moderate classification performance shows the feasibility of reducing the number of patients that unnecessarily undergo invasive FFR measurements.