Direct Automatic Coronary Calcium Scoring in Cardiac and Chest CT

TL;DR

This paper introduces a fast, accurate, and computationally efficient deep learning method for automatic coronary calcium scoring in cardiac and chest CT scans, enabling rapid CVD risk assessment.

Contribution

A novel ConvNet-based approach that directly regresses calcium scores from CT scans, bypassing segmentation and significantly reducing computation time.

Findings

Achieved correlation coefficient of 0.98 with manual scores.

Predicted scores in less than 0.3 seconds per scan.

Almost perfect agreement with manual CVD risk categorization.

Abstract

Cardiovascular disease (CVD) is the global leading cause of death. A strong risk factor for CVD events is the amount of coronary artery calcium (CAC). To meet demands of the increasing interest in quantification of CAC, i.e. coronary calcium scoring, especially as an unrequested finding for screening and research, automatic methods have been proposed. Current automatic calcium scoring methods are relatively computationally expensive and only provide scores for one type of CT. To address this, we propose a computationally efficient method that employs two ConvNets: the first performs registration to align the fields of view of input CTs and the second performs direct regression of the calcium score, thereby circumventing time-consuming intermediate CAC segmentation. Optional decision feedback provides insight in the regions that contributed to the calcium score. Experiments were performed using 903 cardiac CT and 1,687 chest CT scans. The method predicted calcium scores in less than 0.3 s. Intra-class correlation coefficient between predicted and manual calcium scores was 0.98 for both cardiac and chest CT. The method showed almost perfect agreement between automatic and manual CVD risk categorization in both datasets, with a linearly weighted Cohen's kappa of 0.95 in cardiac CT and 0.93 in chest CT. Performance is similar to that of state-of-the-art methods, but the proposed method is hundreds of times faster. By providing visual feedback, insight is given in the decision process, making it readily implementable in clinical and research settings.

Tables

Table 1. TABLE I: Number of scans per CVD risk category for training, validation, and test sets. CVD risk categorization is based on the total Agatston score per scan: I: very low < 1 absent 1 <1 , II: low [ 1 , 10 ) 1 10 [1,10) , III: moderate [ 10 , 100 ) 10 100 [10,100) , IV: moderately high [ 100 , 400 ) 100 400 [100,400) , V: high ≥ 400 absent 400 \geq 400

		I	II	III	IV	V
Cardiac CT	Training	120	14	33	29	41
	Validation	68	14	28	15	11
	Test	260	49	89	70	62
Chest CT	Training	272	76	207	205	252
	Validation	39	14	46	30	40
	Test	128	42	99	112	125

Table 2. TABLE II: Efficient ConvNet architectures were used for atlas-registration as well as calcium scoring.

Atlas-Registration ConvNet	Calcium Scoring ConvNet
512$\times$512$\times$N 3-D input	256$\times$256 2-D input
6$\times$6$\times${1,2} Avg. Pooling	224$\times$224 cropping
32*3$\times$3$\times$3 Convolutions	32*3$\times$3 Convolutions
32*2$\times$2$\times$2 Avg. Pooling	32*2$\times$2 Max Pooling
32*3$\times$3$\times$3 Convolutions	32*3$\times$3 Convolutions
32*2$\times$2$\times$2 Avg. Pooling	32*2$\times$2 Max Pooling
32*3$\times$3$\times$3 Convolutions	32*3$\times$3 Convolutions
32*2$\times$2$\times$2 Avg. Pooling	32*2$\times$2 Max Pooling
32*3$\times$3$\times$3 Convolutions	32*3$\times$3 Convolutions
32*3$\times$3$\times$3 Convolutions	32*2$\times$2 Max Pooling
Global Avg. Pooling	32*3$\times$3 Convolutions
	32*2$\times$2 Max Pooling
	32*3$\times$3 Convolutions
	32*2$\times$2 Max Pooling
64 Fully Connected Nodes	64 Fully Connected Nodes
64 Fully Connected Nodes	64 Fully Connected Nodes
6 Output Nodes	1 Output Node

Table 3. TABLE III: Confusion matrices showing agreement in CVD risk categorization based on the total Agatston scores: I: very low < 1 absent 1 <1 , II: low [ 1 , 10 ) 1 10 [1,10) , III: moderate [ 10 , 100 ) 10 100 [10,100) , IV: moderately high [ 100 , 400 ) 100 400 [100,400) , V: high ≥ 400 absent 400 \geq 400 . The method is evaluated separately on the test sets of cardiac CTs (left) and chest CTs (right). The corresponding linearly weighted κ 𝜅 \kappa is shown below the confusion matrices.

		Predicted
		I	II	III	IV	V
Reference	I	259	0	1	0	0
	II	9	36	4	0	0
	III	2	3	82	2	0
	IV	0	1	2	65	2
	V	0	0	0	11	51
		Cardiac CT $\kappa =0.95$

Table 4. TABLE IV: Quantitative evaluation of visual feedback. Evaluation was performed segmenting CAC lesions with the visualization feedback. An optimal threshold was selected using precision recall analysis on the validation data shown in Figure  8 . Final results show that visualization by the heatmap is is as accurate on the validation as on the test set.

	Cardiac CT	Chest CT
Precision	0.77	0.78
Recall	0.85	0.86
Accuracy	0.99	0.99
F1 (Dice) score	0.81	0.82

Table 5. TABLE V: Comparison of direct calcium scoring experiments using various datasets with and without atlas-registration. The number of CTs in the training set as well as the number of slices are given for each experiment. Additionally we provide the fraction of slices having a calcium score > 0 absent 0 >0 . The results indicate that a calcium scoring ConvNet can be trained on one type of data and evaluated on another.

						Evaluated on:
						Cardiac CT			Chest CT
		Data	CTs	Slices	Fraction CAC	$\kappa$	Acc.	ICC	$\kappa$	Acc.	ICC
Trained on:	Non-Registered	Cardiac CT	237	10,468	10.4%	0.92	0.89	0.89	0.46	0.41	0.24
		Chest CT	1,012	211,353	06.6%	0.48	0.59	0.24	0.91	0.86	0.93
		Cardiac + Chest CT	1,239	221,821	06.7%	0.90	0.86	0.87	0.92	0.88	0.94
	Registered	Cardiac CT	237	10,016	10.9%	0.92	0.88	0.97	0.86	0.79	0.90
		Chest CT subset	237	11,716	14.8%	0.91	0.86	0.95	0.90	0.85	0.93
		Cardiac + Chest CT subset	574	21,732	13.0%	0.94	0.92	0.99	0.91	0.88	0.97
		Chest CT	1,012	100,379	13.8%	0.94	0.91	0.98	0.93	0.89	0.98
		Cardiac + Chest CT	1,239	110,395	13.5%	0.95	0.93	0.98	0.93	0.90	0.98

Table 6. TABLE VI: Results of state-of-the-art automatic calcium scoring methods in cardiac CT and chest CT, direct calcium scoring methods in chest CT, and our proposed method. For each of the studies the number of scans are given that were used for evaluation. To allow better comparison, similar statistics are reported as described in the other studies: ICC quantifies agreement between automatic scores and manual reference scores, and correlation is reported with Pearson’s ρ 𝜌 \rho . Additionally, linearly weighted κ 𝜅 \kappa and accuracy are reported for three different stratifications into risk categories: (A) five categories as used in [ 17 ] and [ 19 ] { < 1 absent 1 <1 , [ 1 , 10 ) 1 10 [1,10) , [ 10 , 100 ) 10 100 [10,100) , [ 100 , 400 ) 100 400 [100,400) , ≥ 400 absent 400 \geq 400 }; (B) four categories as used in  [ 24 ] { ≤ 10 absent 10 \leq 10 , ( 10 , 100 ] 10 100 (10,100] , ( 100 , 1000 ] 100 1000 (100,1000] , > 1000 absent 1000 >1000 }; and (C) five categories as used in  [ 15 ] and [ 46 ] { [ 0 , 10 ) 0 10 [0,10) , [ 10 , 100 ) 10 100 [10,100) , [ 100 , 400 ) 100 400 [100,400) , [ 400 , 1000 ) 400 1000 [400,1000) ≥ 1000 absent 1000 \geq 1000 }. Methods evaluated on identical datasets can be compared directly. The methods by Wolterink et al. [ 19 ] , Lessmann et al. [ 24 ] , and our proposed method were evaluated on systems with an Intel Xeon E5-1620 CPU, 32 GB of internal memory, and an NVIDIA Titan X GPU.

Data

Correlation

A

B

C

Execution time

Source

Number

ICC

$\rho$

$\kappa$

acc.

$\kappa$

acc.

$\kappa$

acc.

CPU

GPU

Cardiac CT

Wolterink et al.[19]

UMCU

530

0.96

–

0.95

0.91

–

–

–

–

20 min

–

Proposed method

UMCU

530

0.97

0.99

0.95

0.93

0.95

0.96

0.94

0.93

5 s

0.15 s

Chest CT

Cano-Espinosa et al. [46]

COPDGene

1,000

–

0.93

–

–

–

–

0.80

0.76

–

–

Lessmann et al. [24]

NLST

506

–

–

–

–

0.91

0.91

–

–

–

7 min

Proposed method

NLST

506

0.98

0.97

0.93

0.90

0.92

0.91

0.93

0.90

11 s

0.30 s

Table 7. TABLE VII: Results of the proposed method on orCaScore challenge data. Left: The confusion matrix shows agreement in CVD risk categorization based on the total Agatston scores: I: 0 0 , II: [ 1 , 100 ) 1 100 [1,100) , III: [ 100 , 300 ) 100 300 [100,300) , IV: > 300 absent 300 >300 . The corresponding linearly weighted κ 𝜅 \kappa is shown below the confusion matrix. Right: Comparison with other methods evaluated in the challenge  [ 25 ] .

		Predicted
		I	II	III	IV
Reference	I	8	0	0	0
	II	0	12	0	0
	III	0	0	8	0
	IV	0	0	1	11
	$\kappa =0.98$

Table 8. TABLE VIII: Intraclass correlation coefficient (ICC) for per-artery calcium scores. Since CVD risk categories are not defined for per-artery scores, CVD risk categorization was evaluated with linearly weighted κ 𝜅 \kappa and accuracy (Acc.) on the total calcium scores obtained by summation.

	Per-artery ICC			Total scores
	LAD	LCX	RCA	$\kappa$	Acc.	ICC
Cardiac CT	0.93	0.88	0.97	0.94	0.91	0.97
Chest CT	0.91	0.80	0.98	0.92	0.88	0.96

Equations

$$T_{\text{3D}}=\left[\begin{array}{cccc}1& 0& 0& t_x\\ 0& 1& 0& t_y\\ 0& 0& 1& t_z\\ 0& 0& 0& 1\end{array}\right]\left[\begin{array}{cccc}\cos {\theta }_z& -\sin {\theta }_z& 0& 0\\ \sin {\theta }_z& \cos {\theta }_z& 0& 0\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right]\left[\begin{array}{cccc}{s}_{xy}& 0& 0& 0\\ 0& s_{xy}& 0& 0\\ 0& 0& s_z& 0\\ 0& 0& 0& 1\end{array}\right]$$

$$T_{\text{2D}}=\left[\begin{array}{ccc}1& 0& t_x\\ 0& 1& t_y\\ 0& 0& 1\end{array}\right]\left[\begin{array}{ccc}\cos {\theta }_z& -\sin {\theta }_z& 0\\ \sin {\theta }_z& \cos {\theta }_z& 0\\ 0& 0& 1\end{array}\right]\left[\begin{array}{ccc}{s}_{xy}& 0& 0\\ 0& s_{xy}& 0\\ 0& 0& 1\end{array}\right]$$

$$\text{Agatston Score}=\sum _{S\in V}\sum _{l\in S}{A}_l{w}_l\frac{i_S}{t_S}.$$

$$\text{Pseudo-Agatston Score}=\sum _{l\in S}{n}_l\cdot s_{xy}\cdot w_l.$$

$$L=|\hat{y}-\ln (y+1)|,$$

Full text

I Introduction

Cardiovascular disease (CVD) is the global leading cause of death[1]. To reduce the burden of cardiovascular disease the World Health Organization underlines the need for early detection and treatment of individuals with CVD or those who are at high cardiovascular risk due to the presence of one or more risk factors [2]. A strong and independent risk factor for CVD events, e.g. myocardial infarction, is the quantity of coronary artery calcium (CAC) [3, 4, 5]. Quantification of CAC, i.e. calcium scoring, is typically performed in dedicated non-contrast-enhanced ECG-synchronized cardiac CT scans[4]. Alternatively, calcium scoring can be performed in other non-contrast-enhanced CTs that visualize the heart; e.g. in low-dose CT attenuation correction scans acquired in hybrid PET/CT and SPECT/CT [6, 7], or in radiation therapy planning CTs of breast cancer patients [8]. Furthermore, it has been shown that calcium scoring in lung screening low-dose chest CT scans is a predictor for all-cause mortality [9, 10]. In fact, in the National Lung Screening Trial (NLST) CVD was the leading cause of mortality [11]. Thus, CAC quantification, especially as an unrequested finding, has garnered much attention.

Clinically, calcium scoring is performed by experts who manually identify CAC in CT image slices. This is a tedious process of finding and selecting high density voxels in the coronary arteries—commonly defined as two or more connected voxels above 130 Hounsfield Units (HU). In scans not dedicated to calcium scoring this can be particularly cumbersome because of high noise, low resolution, and motion artifacts. Subsequently, when lesions are identified, region growing is used to fully segment the calcified lesions. Finally, after all CAC lesions have been segmented, CAC is quantified using the Agatston score [12]. The Agatston score takes into account the lesion area and the weighted maximum density of the lesion. This score can be used to stratify patients into risk categories [13].

The additional cost involved with manual calcium scoring makes the process prohibitive in settings where it is not the primary request. To simplify the task, qualitative stratification into CVD risk groups was proposed [14, 10]. Qualitative calcium scoring is faster and it demonstrates good inter-rater agreement. However, such an analysis still demands experts who closely inspect the scans. With the ever-increasing amount of scans and the increasing interest in calcium scoring, especially as an unrequested finding, the use of fully-automatic methods might be the preferred direction.

Several automatic methods have been introduced for calcium scoring in non-contrast-enhanced CT, ranging from rule-based approaches [15, 16], to the better performing conventional machine learning approaches[17, 18, 19, 20] and recent deep learning approaches [21, 22, 23, 24]. The main difficulty in automatic calcium scoring is to differentiate CAC from other dense structures. Obviously, CAC exclusively resides in the walls of the coronary arteries, thus most of the automatic methods exploit this prior knowledge.

Išgum et al. [17] introduced the first method for automatic calcium scoring in chest CT. CAC lesions were described with features and subsequently classified using a two-stage classification approach of k-nearest neighbor and support vector classification. Among texture, size, and shape features, highly important for CAC identification, were the location features. Location features were determined by registering an input image to an atlas image and by extracting the location features from a map of a priori spatial probabilities of CAC. The probability map was created from known CAC locations in 237 chest CTs that were registered to a single priorly chosen atlas image. Shahzad et al. [18] used a similar machine learning approach for calcium scoring in cardiac CT, but they employed pair-wise deformable image registration to ten atlases that encoded the coronary arteries. The atlases were made from 85 contrast enhanced CT angiography scans with annotated coronary arteries. The methods of Išgum et al. [17] and Shahzad et al. [18] relied on feature selection methods to reduce dimensionality. Wolterink et al. [19] circumvented feature selection by using an extremely randomized trees classifier. Their method also depended on location features that were obtained by deformable image registration of ten atlases with encoded coronary arteries, but these were obtanied from non-contrast-enhanced CTs. Durlak et al. [20] combined the principles of the aforedescribed methods: they employed a random forest and made an a priori probability map of coronary arteries locations, made from automatically extracted coronary arteries from cardiac CT angiography images. Instead of using time-consuming deformable image registration to align input images and atlas images, they achieved a speed-up by using affine registration. Similarly, other methods employed information from CTA to aid calcium scoring in cardiac CT. These methods were specifically designed for the coronary calcium score (orCaScore) challenge, and employed rule-based image analysis or conventional machine learning [25].

Most recently proposed methods employ deep learning methods for automatic calcium scoring, in particular convolutional neural networks (ConvNets). ConvNets are known for their automatic feature extracting capabilities and alleviate the need for handcrafting features. Wolterink et al. [22] used ConvNets to classify CAC in cardiac CT angiography scans. All voxels were classified using a pair of ConvNets. One ConvNet identified voxels likely to be CAC and discarded the majority of non-CAC-like voxels such as lung and fatty tissue. The other ConvNet more precisely discriminated between CAC and CAC-like negatives. In the method of Lessmann et al.  [23] a single ConvNet was used that classified candidate CAC lesions in lung screening chest CTs. To simplify the classification tasks, both these deep learning methods used an additional ConvNet that localized the heart with a bounding box [26]. More recently, the method of Lessmann et al.[24] fully exploited the feature extraction capabilities of ConvNets without dedicated localization methods. They employed two sequential ConvNets to classify CAC as well as aortic valve, mitral valve, and aorta calcifications in chest CT. The first ConvNet identified candidate calcifications based on their location, and the second ConvNet refined the classification results by reducing false positive errors.

While all aforementioned methods use different strategies, they all follow a workflow similar to current clinical calcium scoring: CAC is first identified and thereafter quantified. The automatic methods show high accuracy, but often at considerable computational cost. Employing these methods on large datasets would require dedicated servers. To alleviate computational cost, we propose a workflow that circumvents intermediate identification and that performs direct quantification (see Figure 1). Direct quantification has proven to be useful for atrial and ventricle volume quantification [27, 28, 29]. Furthermore, attempts are being made to use it for calcium scoring. In our preliminary study we presented a direct calcium scoring method that uses 2-D ConvNet regression [30, 31]. The method performs direct calcium scoring in extracted image slices from bounding boxes cropped around the heart. In a recently proposed method, Cano-Espinosa et al. used a 3-D regression ConvNet for direct calcium scoring in downsampled CT volumes also cropped around the heart. However, their method could not be used in 14% of the scans, because heart localization failed. Furthermore, previously proposed automatic calcium scoring methods are dedicated to either cardiac CT or chest CT. These methods required retraining for application in other types of CT [8, 32].

We present an automatic method that performs real-time direct calcium scoring in different types of non-contrast-enhanced CT. Unlike previous methods that focused on a single type of CT, the proposed method is able to perform calcium scoring directly in multiple types of CT by using an unsupervised deep learning atlas-registration method to align their fields of view (FOVs). For this we employ two ConvNets: one for atlas-registration and one for calcium scoring, as shown in Figure 2. The atlas-registration ConvNet makes the FOV of input CT images alike using Deep Learning Image Registration (DLIR) [33, 34] further developed to facilitate atlas-registration. Subsequently, a calcium scoring ConvNet predicts the calcium score in image slices mimicking clinical calcium scoring with the Agatston score. When desired, decision feedback can be queried for every slice with a predicted calcium score. For this purpose, a visual attention heatmap accurately reveals the regions that contributed to the calcium score. The method provides robust and accurate predictions of calcium scores and it is computationally efficient, obtaining an Agatston score in less than 0.3 s in cardiac and chest CT.

II Data

This study included two datasets used in previous studies that presented automatic coronary calcium scoring in cardiac CT [19] and in chest CT [24]. To allow a direct comparison of methods, the original training, validation, and test set distributions were used.

II-A Cardiac CT

The set of 903 cardiac CT scans (age range: 18 to 88 years, 31% women) originates from a set of routinely acquired scans for clinical calcium scoring of the University Medical Center Utrecht, Utrecht, The Netherlands. The need for informed consent was waived by the local Medical Research Ethics Committee. Scans were acquired with a 256-detector row Philips Brilliance iCT scanner (tube voltage 120 kVp, tube current 55 mAs) during a single breath-hold, with ECG-triggering and without contrast enhancement. The images were reconstructed to 3 mm slice thickness and slice increment with in-plane resolution ranging from 0.29 mm to 0.49 mm, depending on patient size. The dataset was divided into 237 scans for training, 136 scans for validation, and 530 scans were in the hold-out test set only used for final evaluation.

II-B Chest CT

The set of 1,687 chest CT scans (age range: 43 to 74 years, 39% women) originates from a set of 6,000 available baseline scans from the National Lung Screening Trial (NLST) [11]. All scans were acquired during inspiratory breath-hold without contrast enhancement. Scans were acquired in 31 different hospitals with 120 or 140 kVp tube voltage and 30-160 mAs tube current. Axial images slices were reconstructed with varying kernels, varying slice thickness (1.00-3.00 mm), varying slice increments (0.63-3.00 mm), and with varying in-plane resolutions (0.49-0.98 mm per voxel). In our study, scans with less than 100 slices or slices thicker than 3.00 mm were not considered, because they were not adequate for calcium scoring. Furthermore, the scans were resampled to 3.00 mm slice thickness and 1.50 mm slice increment to make the scans suitable for calcium scoring [35]. The dataset was divided into 1,012 scans for training, 169 scans for validation, and 506 scans were in the hold-out test set only used for final evaluation.

II-C Reference standard

The reference standard was defined by experts who manually identified CAC lesions in the scans. CAC lesions were segmented following a standard procedure: region growing was used to select 26-connected voxels $\geq$130 HU. In the chest CTs with low radiation dose this procedure could lead to faulty segmentations (i.e. leakage) because of excessive noise. In such cases annotations were manually corrected by voxel painting [24]. Agatston scores were calculated in each axial slice for training. Total Agatston scores for each scan were calculated for final evaluation. Additionally, each subject was assigned to one of five CVD risk categories [13] based on the Agatston score: very low: $<$1; low: [1, 10), moderate: [10, 100), high: [100, 400), very high: $\geq$400. Table I provides an overview of the number of scans per risk category per dataset.

III Methods

The method employs two ConvNets in sequence (Figure 2). The first ConvNet registers input CTs to an cardiac CT atlas-image. The second ConvNet performs calcium scoring. When desired, visual feedback can be queried for image slices with a score. For this purpose an attention heatmap reveals the regions that contributed to the calcium score.

III-A Atlas-registration strategy

An atlas-registration ConvNet ensures that all input images have a similar FOV and resemble a cardiac CT. The ConvNet is trained with a modified version of our framework for Deep Learning Image Registration (DLIR) [33]. The DLIR framework uses an end-to-end unsupervised approach that trains a ConvNet for image registration. Similar to a conventional intensity-based image registration framework it exploits optimization of an image similarity metric. Figure 3 shows the schematics of training an atlas-registration ConvNet using the atlas image as a static fixed image. The task of the ConvNet is to analyze moving images and predict the transformation parameters that warp the moving images to the atlas-image. Image similarity between the atlas and the warped image, is used for backpropagation during training. By optimizing image similarity (e.g. minimizing negative cross correlation) with gradient descent, the atlas-registration ConvNet learns the registration task in an unsupervised manner. After training, the ConvNet can register unseen moving images in one shot.

A cardiac atlas-image is created using an iterative inter-subject registration strategy [36]. With this strategy an initial atlas image is made by averaging multiple images. The atlas image is iteratively refined by registering the individual images to the atlas. Subsequently, the final atlas image is used to train the atlas-registration ConvNets for cardiac and chest CT alignment used for subsequent calcium scoring.

III-B Atlas-registration ConvNet training

For registration we propose a global 3-D rigid registration model with six degrees of freedom (shown in Figure 4). The model allows translations in any direction, but rotations are restricted to the axial ($z$) axis. Furthermore, scaling in the axial plane is isotropic and independent from scaling along the axial axis. These restrictions preserve the relation of reference Agatston scores that are defined on the original (unregistered) axial image slices. This facilitates training of the subsequent calcium scoring ConvNet.

We use a computationally efficient ConvNet architecture that is listed in Table II. For fast analysis, images are downsampled close to 3 mm isotropic voxel dimensions; i.e. 6$\times$6$\times$1 downsampling for cardiac CT, and 6$\times$6$\times$2 downsampling for chest CT using average pooling. The ConvNet has three alternating layers of 3$\times$3$\times$3 convolutions and 2$\times$2$\times$2 average pooling and those are followed by two layers of 3$\times$3$\times$3 convolution. To facilitate a fixed output, global average pooling is applied before connection with two fully connected layers. The final output layer has six nodes, one for each transformation parameter. Throughout the network exponential linear units are used for activation, except in the output nodes. Three output nodes are unconstrained translation parameters ($t_{x}$, $t_{y}$, $t_{z}$), the rotation parameter ($\theta_{z}$) is constrained with a hyperbolic tangent between $-\pi$ and $\pi$, and the two scaling parameters ($s_{xy}$, $s_{z}$) are constrained with a hyperbolic tangent between $0.25$ and $4$ scaling factors. These output parameters are used to constitute the following 3-D transformation matrix:

[TABLE]

III-C Atlas-registration ConvNet inference

We train an atlas-registration ConvNet for 3-D registration, but we use it for slice selection and 2-D warping. As a consequence, correspondence is guaranteed between warped axial slices and the per-slice calcium scores. Axial image slices are extracted from the original image from $t_{z}$ to $t_{z}+d_{z}/s_{z}$, where $d_{z}$ is the depth of the atlas image along the axial axis. These slices are resampled using bi-linear interpolation to a 256$\times$256 grid with the following 2-D transformation matrix:

[TABLE]

III-D Calcium scoring ConvNet

The calcium scoring ConvNet employs direct regression to predict an Agatston score from input axial image slices. The choice of 2-D ConvNets, in favor of 3-D ConvNets, is based on the number of samples that are available for training. There are more image slices available than image volumes. Furthermore, 2-D image analysis mimics clinical calculation of the Agatston calcium score that is performed in 2-D axial slices:

[TABLE]

where $l$ is a 2-D CAC lesion in a slice $S$ of a CT volume $V$. $A_{l}$ is the area of the lesion. The weighted intensity $w_{l}$ is based on the maximum radio-density in HU of a 2-D lesion in the following manner: 1 = [130, 200), 2 = [200, 300), 3 = [300, 400), and 4 = $\geq$400. The Agatston score is corrected when image slices are overlapping, thus when slice increment $i_{S}$ is not equal to slice thickness $t_{S}$ [37].

Agatston scores are dependent on the CAC lesion area. Given that input images have different voxels sizes, we chose to simplify the prediction task by determining a pseudo-Agatston score. This score is obtained by cancelling out the axial pixel dimensions, the slice increment, and the slice thickness of the original Agatston score. The resulting target is the product of the number of voxels in a lesion $n_{l}$, the predicted slice scaling factor $s_{xy}$, and the weighted intensity $w_{l}$:

[TABLE]

The calcium scoring ConvNet uses an efficient architecture that is listed in Table II. It analyzes random image croppings of 224$\times$224 pixels during training and center croppings during application. It has alternating layers of 3$\times$3 convolutions and $2\times 2$ max pooling, followed by two fully connected layers, and an output layer of one node. Throughout the network batch normalization [38] is used and exponential linear units are used for activation[39]. The final output node has a linear output to facilitate continuous prediction. However, given that clinically used CVD risk categories are exponentially increasing, the task of the calcium scoring ConvNet was modified to learn a log-transform of the pseudo-Agatston score:

[TABLE]

where $\hat{y}$ is the predicted score, and $y$ is the reference pseudo-Agatston score. The log-transform induces relatively high penalties for erroneous low calcium score predictions, and relatively low penalties for erroneous high calcium score predictions. Consequently, higher precision is forced for lower calcium burden, which is favorable for CVD risk stratification. During application of the calcium scoring ConvNet, the predicted outputs are converted to the original Agatston scores.

III-E Decision feedback

By employing regression of calcium scores, we circumvent time-consuming intermediate segmentation. On the other hand, it may be desirable to visualize regions in image slice that contributed to the calcium score. Inspired by the study of Zeiler and Fergus [40], we provide such visualization by using a deConvNet. The deConvNet uses the same operations of filtering and pooling as a ConvNet, but in reverse order from output to input. The reverse operations map the activities back to the input pixel space, and it shows which input patterns originally contributed to the activations in the feature maps. To obtain a smooth visual attention heatmap, the deConvNet is applied until the third convolutional layer, by taking the absolute value per feature of this layer, and by summing these features along the feature map dimension to get 2-D matrix. Using third order interpolation we obtain a smooth map that can be superpositioned on the image slice as a heatmap. This resulting heatmap visualizes attention by highlighting the regions that contributed to the Agatston score.

IV Evaluation

Automatically predicted per-subject Agatston scores were compared with manually determined reference scores. Evaluations were performed on the hold-out test sets which were not used during method development. Two-way mixed intra-class correlation coefficient (ICC) for absolute agreement was computed and Bland-Altman analysis was performed to evaluate bias between predicted and reference Agatston scores. In addition, for each subject, CVD risk category was determined based on the Agatston score as defined in section II-C. Agreement between predicted and reference CVD risk categories was determined using accuracy and Cohen’s linearly weighted kappa ($\kappa$).

V Experiments and results

In this section we evaluate the atlas-registration ConvNet, the calcium scoring ConvNet, and the quality of decision feedback. In addition, we will evaluate whether the calcium scoring ConvNet requires to be trained on all data, or whether it can be trained on one dataset and applied to the other. Finally, we will compare state-of-the-art automatic calcium scoring methods with the proposed method. All experiments were performed with Theano [41], Lasagne [42], and OpenCV [43] on an Intel Xeon E5-1620 3.60 GHz CPU with an NVIDIA Titan X GPU.

V-A Atlas-registration ConvNet

Figure 5a shows the initial atlas image that was created by aligning all cardiac training images using their geometric centroids. We chose the median dimensions and voxels sizes of all the cardiac training images define the atlas image space. The atlas can be iteratively refined, but given the constraints of the global registration model used here, only one update was sufficient. The final atlas image, shown in Figure 5d, was used to train the atlas-registration ConvNets for cardiac and chest CT alignment. Thus, in total three ConvNet instances were trained: one to create an atlas image, one for cardiac CT alignment, and one for chest CT alignment. All ConvNets were trained in 15,000 iterations with mini-batches containing 32 randomly selected images. Training took about 40 hours per ConvNet. Adam [44] was used with a learning rate of 0.001 for mini-batch gradient descent. To illustrate performance of the atlas-registration ConvNets, Figure 5 shows images before and after registration. Figure 5b shows the average image of the 530 cardiac CT images from the test set before registration and Figure 5e shows these images after registration. Similarly, Figure 5c shows an average image of the 506 chest CTs before registration and Figure 5f shows these after registration. Note the similarity of the registered image with the refined atlas image shown in Figure 5d.

Quantitative evaluation of registration results revealed that registration erroneously cropped CAC out of the selected slices. Between one and four image slices containing CAC were not selected in three cardiac CTs and three chest CTs. Upon closer inspection, two of the chest CTs had calcifications in the aortic arch and descending aorta incorrectly labeled as CAC in the reference, thereby affecting CVD risk categorization. Nevertheless, these annotations were left uncorrected in further analysis to facilitate a fair comparison with previously developed methods. The registration errors did not have an adverse effect on CVD risk categorization in the other cases.

V-B Calcium scoring ConvNet

The calcium scoring ConvNet was trained in 150,000 iterations using Adam [44]. Training took 21 hours with 100 image slices per mini-batch randomly selected from the registered image slices taken from the cardiac and chest CT training sets. High imbalance between the minority of slices with a calcium score and the majority of slices with zero calcium score prevented convergence during ConvNet training. To ensure convergence, the amount of image slices with CAC (Agatston score $>0$) and without CAC (Agatston score $=0$) were balanced during training. To prevent bias, training continued on the full imbalanced training set after 10,000 iterations. Additionally, we ensured stable convergence by decreasing the learning rate to 10% of its previous value every 50,000 iterations.

After training, the test sets were used to evaluate the calcium scoring ConvNet. Per-subject scores show high intraclass correlation coefficients (ICC); the ICC for cardiac CT and chest CT were both 0.98 with 95% confidence intervals of 0.98 to 0.99. Slight positive bias in cardiac and chest CT is visualized with the Bland-Altman plots shown in Figure 6. This was mainly caused by overestimations of the higher Agatston scores. However, this was not noticeable in CVD risk stratification. Table III shows confusion matrices of predicted risk categories vs. the manual reference standard. In cardiac CT calcium scoring only four scans were two categories off, and in chest CT calcium scoring eight scans were two categories off. The scan that was three categories off was a scan with incorrectly annotated aorta calcium, as discussed in the previous section. Nonetheless, overall agreement was almost perfect [45] with Cohen’s linearly weighted $\kappa$s of 0.95 in cardiac CT and 0.93 in chest CT. Accuracy in CVD risk categorization was 0.93 for cardiac CT and 0.90 for chest CT. Because efficient network architectures are used, the method is able to achieve high speed when used on a single CPU core: within 5 s a score for cardiac CT is obtained and within 11 s a score for chest CT is obtained. When using a GPU, calcium scoring can be performed in real-time. Including image registration and image resampling, a calcium score for cardiac CT is obtained in less than 0.15 s and for chest CT in less than 0.30 s.

V-C Decision feedback

Decision feedback visualizes attention of the calcium scoring ConvNet. This feedback informs and end-user about the regions that contributed to the calcium score. Figure 7 shows examples of such feedback. The feedback helps an expert to quickly navigate and evaluate the image slices containing CAC.

We propose visual feedback as an optional qualitative tool, but we have performed a quantitative analysis to provide insight in its accuracy. To obtain quantitative results we analyzed heatmaps for slices with predicted calcium scores. The heatmaps were warped to the original image spaces by using the inverse transformation matrices. The values of the heatmaps were scaled between 0 and 1 to mimic probability maps for CAC candidate voxels. CAC candidates were defined as high density 26-connected voxels with a volume between 1.5 and 1,500 mm3[19]. For evaluation of these maps we performed precision-recall analysis (Figure 8). We have defined an optimal threshold by selecting the maximum F1 (i.e. Dice) score on the validation set. Table IV shows the obtained scores using the selected threshold on the test sets. The results show that detection performance is very accurate on the validation set as well as the test set.

Additionally, decision feedback aided our analysis by clarifying incorrect calcium scores. Decision feedback revealed that the largest CVD miscategorizations were not caused by incorrect quantification but by incorrect recognition of CAC. Figure 9 shows six examples of the largest miscategorizations made by the calcium scoring ConvNet. The majority of errors were made in identification of calcifications near the coronary artery ostia. Calcifications near the ostia can be partly in the aorta and partly in the coronary artery. These calcifications are difficult to distinguish, especially when no information of neighboring slices is available.

V-D Influence of training data and registration

For clinical application it would be useful to investigate whether the method needs training data from both datasets or if data from one set would suffice, and we investigated the influence of atlas-registration is required. Thus, we performed experiments using different combinations of training data with and without atlas-registration, as listed in Table V. The calcium scoring ConvNets were trained with either cardiac CT images, chest CT images, or a combination thereof. To balance cardiac and chest CT data, a subset of chest CT images was created by taking images from 237 randomly selected subjects and by removing every other slice in the chest CT images. Additionally, the histograms shown in Figure 10 provide insight in the distribution of calcium amount in the training data. Note that the chest CT subset has a very similar distribution compared to the cardiac CT training set.

The best performance was achieved using atlas-registration with a calcium scoring ConvNet trained on all cardiac and chest CT images. Lower scores are found when a calcium scoring ConvNet is only trained with cardiac CT or the subset of chest CTs. However, combining the two datasets increased the scores notably, giving a performance close to the ConvNet trained with all images. Furthermore, the results show that atlas-registration facilitated training on one type of data and high performance on the other: the ConvNet trained with the full set of chest CTs achieved a high performance on the cardiac CT test images that was very close to the best results.

V-E Comparison with other methods

Table VI shows a comparison with other state-of-the-art calcium scoring methods by Wolterink et al. [19] and Lessmann et al. [24] using the same datasets. The proposed method achieves similar performance compared to these methods, but it is hundreds of times faster. Even when ran on a single core of a CPU, the method achieves high speed. Additionally, we listed results from other direct calcium scoring methods by González et al. [15] and Cano-Espinosa et al. [46] using chest CT data from the COPDGene study [47]. We provide similar performance metrics to give an indication, but please note that a direct comparison between these methods and ours was not possible.

V-F Performance on orCaScore data

We evaluated our method on data from the orCaScore challenge [25]. This challenge provides data to evaluate a method for coronary calcium scoring. The data consists of non-contrast enhanced ECG-triggered cardiac CT acquired on CT scanners from four different vendors from four different hospitals. Training data is provided, but we evaluated our method on the test set of 40 patients without retraining. Table VII shows the obtained confusion matrix and lists the results of dedicated cardiac CT calcium scoring methods that competed in the challenge. Given that our method does not differentiate between location of CAC, we only provide total calcium scoring results.

V-G Per-artery calcium scores

Routine coronary artery calcium scoring is typically performed per artery. Currently, only total coronary calcium scores are reported and used for CVD risk prediction. For research purposes, per-artery calcium scores might provide interesting additional information. Hence, we evaluated performance of the proposed method for per-artery calcium scoring, i.e. scoring in the the LAD, LCX, and RCA. We chose to combine CAC scores in the LM and LAD, since it is difficult, if not impossible, to differentiate them in chest CT scans. The direct scoring ConvNet was adapted by changing the number of output nodes from one to to three. Similar to the experiment described in Section V-B, training started with a balanced set of image slices with and without calcium scores for the first 10,000 iterations and continued with the full set of image slices thereafter. Additionally, each mini-batch had at least three image slices containing each type of arterial calcification. Risk categories are clinically not defined for per-artery calcium scores, but they are obtained for total calcium scores by summation of per-artery scores. The results are listed in Table VIII.

VI Discussion

We have presented a method for automatic coronary calcium scoring in cardiac CT and chest CT. The method uses an atlas-registration ConvNet to align FOVs making input images alike. The atlas-registration ConvNet is trained for 3-D registration, but its rigid model is constrained to enable 2-D slice selection and 2-D image warping. Selected and warped input image slices are presented to a calcium scoring ConvNet that directly predicts the Agatston score in these slices. The method circumvents time-consuming CAC segmentation. To provide decision feedback, a visual attention heatmap can be generated that shows the regions in an image contributing to the calcium score. The method achieves excellent agreement for calcium score prediction for CVD risk categorization compared to manual calcium scoring. The method achieves similar performance compared to state-of-the-art methods, but achieves it hundreds of times faster.

In preliminary experiments we found that only a small ConvNet architecture was able to learn direct calcium scoring. Large ConvNet architectures architectures were unstable and failed to converge during training. By limiting the degrees of freedom of a ConvNet, i.e. by using a small architecture, we were able to train a ConvNet that learned to differentiate coronary calcification from other types of calcification e.g. aorta calcification, pericardium calcification, and heart valve calcification.

To simplify the problem we extracted bounding boxes around the heart in our preliminary work [30, 31]. However, this was a supervised method that classified presence of the heart in image slices. In case of noisy images, consecutive image slices could have discontinuous predictions. Discontinuous predictions resulted in an incorrect bounding box extracting a partial heart. For atlas-registration used in our current work this is not an issue.

The atlas-registration ConvNets were highly successful in pre-alignment of input CTs, i.e. in slice selection and image warping. Only 4 out of 1,036 test images had slices containing CAC that were missed by erroneous slice selection. Erroneous slice selection was likely caused by incorrect focus of the atlas-registration ConvNet on high contrast areas like the diaphragm. A mask drawn around the heart might steer focus of the ConvNet and might increase registration performance. Alternatively, a simple adjustment could be made by padding slice selection with some slices. Nonetheless, the errors caused by registration had negligible impact on calcium scoring and did not affect CVD risk categorization. Calcium scoring is better with atlas-registration than without it. Moreover, registration allows training and application of direct calcium scoring on datasets with different FOVs.

In general accuracy of predicted Agatston scores was high. Although Bland-Altman analysis showed that the method underestimated subjects with high Agatston scores. In fact, this was by design, because the method estimates a log transformed Agatston score, which induces relatively low precision for higher scores, and high precision for lower scores. Because the clinically used CVD risk categories are based on exponentially increasing Agatston scores, it is obviously more important to differentiate between subjects at low to moderate risk, than to differentiate between subjects at high risk. Thus, we imposed this higher precision on lower Agatston scores. Still, the largest CVD miscategorizations were found in the lower risk categories. Miscategorization was predominantly caused by incorrect identification of CAC and aortic calcifications near the coronary artery ostia. Even manual classification of these calcifications can be very difficult when they spread from the aorta through an ostium into the coronary artery. It often involves inspecting multiple adjacent slices in 3-D. Thus, performance of the method might be improve by exploiting additional 3-D information in future work. Additionaly, performance might improve by increasing input image resolution. The current resolution was chosen based on the majority of chest CT images, being roughly half the resolution of cardiac CTs. Nevertheless, even though all cardiac CTs were downsampled a high performance was obtained in these CTs.

The proposed method shows near perfect agreement in CVD risk categorization compared to manual calcium scoring, even when trained with a relatively low number of scans from a single dataset. Interestingly, training on one type of data allowed the model to be applied to the other type of data without any modifications or transfer learning. However, we found that a model trained on only chest CT led to better results than a model trained only on cardiac CT. One potential reason for this may be the distribution of CAC in the datasets: the population of ex-heavy smokers typically have more CAC [9] than the population undergoing calcium scoring cardiac CT. However, Figure 10 shows that the distribution of CAC in equally sized datasets of cardiac CT and chest CT is similar. An alternative reason could be the presence of motion artifacts, which are nearly absent in ECG-synchronized cardiac CT, but abundant in non-ECG-synchronized chest CT. Therefore, a model trained on chest CT may be more robust to such artifacts. While our experiments indicated that a ConvNet trained on the cardiac and chest CT datasets supplement each other, a calcium scoring ConvNet trained with only chest CTs almost matched performance of the best performing ConvNet. Additionally, we have shown that the method obtained near perfect CVD risk categorization results on cardiac CTs from the orCaScore challenge. The method did not require retraining on representative data from the different hospitals and vendors. Having a single system that can handle potentially any CT scan that visualizes the heart would be very practical in a routine radiology setting. In future work we will investigate whether the method could be readily applied on other types CTs, without requiring retraining or fine-tuning.

Additionally, we have shown that the method can provide per-artery calcium scores. While this is not required for CVD risk categorization, it might be interesting for clinical research. In terms of ICC [48], per-artery calcium scoring achieved good reliability ($>0.75$) in the LCX, and excellent reliability ($>0.90$) in the LAD and the RCA. In addition, determination of CVD risk using combined per-artery scores led to almost perfect agreement (${\kappa>0.90}$) [45]. Nevertheless, performance was slightly better when total calcium was directly determined. This difference in performance may be a consequence of increased complexity of the per-artery scoring task while using the same number of samples for training.

The proposed method can achieve a calcium score hundreds of times faster than previously proposed methods. This is mainly due to one-shot (i.e. non-iterative) registration, and direct quantification using regression. The direct calcium scoring method circumvents time-consuming intermediate segmentation. The method might also be suitable for e.g. determination of volume, (pseudo-)mass, or number of CAC; and for quantification of other lesions or diverse anatomical structures. However, the benefit of using a segmentation approach over direct scoring is that it provides immediate insight to the end-user. We mitigate this shortcoming of direct scoring, by providing decision feedback with a visual attention heatmap. In this way valuable feedback is still provided whenever an end-user requires it.

VII Conclusion

We have presented an automatic method for direct calcium scoring in cardiac CT and chest CT. The method employs two ConvNets, one for atlas-registration to align the FOV of input images to an atlas image made from cardiac CTs and one for direct calcium scoring of input image slices using regression. The method achieves robust and accurate predictions of calcium scores in real-time. By providing visual feedback, insight is given in the decision process, making it readily implementable in a clinical and research settings.