Tensor-based Multimodal Learning for Prediction of Pulmonary Arterial Wedge Pressure from Cardiac MRI

Abstract

Heart failure is a severe and life-threatening condition that can lead to elevated pressure in the left ventricle. Pulmonary Arterial Wedge Pressure (PAWP) is an important surrogate marker indicating high pressure in the left ventricle. PAWP is determined by Right Heart Catheterization (RHC) but it is an invasive procedure. A non-invasive method is useful in quickly identifying high-risk patients from a large population. In this work, we develop a tensor learning-based pipeline for identifying PAWP from multimodal cardiac Magnetic Resonance Imaging (MRI). This pipeline extracts spatial and temporal features from high-dimensional scans. For quality control, we incorporate an uncertainty-based binning strategy to identify poor-quality training samples. We leverage complementary information by integrating features from multimodal data: cardiac MRI with short-axis and four-chamber views, and cardiac measurements. The experimental analysis on a large cohort of $1346$ subjects who underwent the RHC procedure for PAWP estimation indicates that the proposed pipeline has a diagnostic value and can produce promising performance with significant improvement over the baseline in clinical practice (i.e., $\Delta$AUC $=0.10$, $\Delta$Accuracy $=0.06$, and $\Delta$MCC $=0.39$). The decision curve analysis further confirms the clinical utility of our method. The source code can be found at: https://github.com/prasunc/PAWP.

1 Introduction

Heart failure is usually characterized by the inability of the heart to supply enough oxygen and blood to other organs of the body [4]. It is a major cause of mortality and hospitalization  [14]. Elevated Pulmonary Arterial Wedge Pressure (PAWP) is indicative of raised left ventricular filling pressure and reduced contractility of the heart. In the absence of mitral valve or pulmonary vasculature disease, PAWP correlates with the severity of heart failure and risk of hospitalization [1]. While PAWP can be measured by invasive and expensive Right Heart Catheterization (RHC), simpler and non-invasive techniques could aid in better monitoring of heart failure patients. Cardiac Magnetic Resonance Imaging (MRI) is an effective tool for identifying various heart conditions and its ability to detect disease and predict outcome has been further improved by machine learning techniques [3]. For instance, Swift et al. [17] introduced a machine-learning pipeline for identifying Pulmonary Arterial Hypertension (PAH). Recently, Uthoff et al. [18] developed geodesically smoothed tensor features for predicting mortality in PAH.

Cardiac MRI scans contain high-dimensional spatial and temporal features generated throughout the cardiac cycle. The small number of samples compared to the high-dimensional features poses a challenge for machine learning classifiers. To address this issue, Multilinear Principal Component Analysis (MPCA) [11] utilizes a tensor-based approach to reduce feature dimensions while preserving the information for each mode, i.e. spatial and temporal information in cardiac MRI. Hence, the MPCA method is well-suited for analyzing cardiac MRI scans. The application of the MPCA method to predict PAWP might further increase the diagnostic yield of cardiac MRI in heart failure patients and help to establish cardiac MRI as a non-invasive alternative to RHC. Existing MPCA-based pipelines for cardiac MRI [17, 18, 2] rely on manually labeled landmarks that are used for aligning heart regions in cardiac MRI. The manual labeling of landmarks is a cumbersome task for physicians and impractical for analyzing large cohorts. Moreover, even small deviations in the landmark placement may significantly impact the classification performance of automatic pipelines [16]. To tackle this challenge, we leverage automated landmarks with uncertainty quantification [15] in our pipeline. We also extract complementary information from multimodal data from short-axis, four-chamber, and Cardiac Measurements (CM). We use CM features (i.e., left atrial volume and left ventricular mass) identified in the baseline work by Garg et al. [5] for PAWP prediction.

Our main contributions are summarized as follows: 1) Methodology: We developed a fully automatic pipeline for PAWP prediction using cardiac MRI data, which includes automatic landmark detection with uncertainty quantification, an uncertainty-based binning strategy for training sample selection, tensor feature learning, and multimodal feature integration. 2) Effectiveness: Extensive experiments on the cardiac MRI scans of $1346$ patients with various heart diseases validated our pipeline with a significant improvement ($\Delta$AUC $=0.1027$, $\Delta$Accuracy $=0.0628$, and $\Delta$MCC $=0.3917$) over the current clinical baseline. 3) Clinical utility: Decision curve analysis indicates the diagnostic value of our pipeline, which can be used in screening high-risk patients from a large population.

2 Methods

Figure 1: The schematic overview of the PAWP prediction pipeline including preprocessing, tensor feature learning, and performance analysis. The blocks in gray color are explained in more detail in Section 2.

As shown in Fig. 1, the proposed pipeline for PAWP prediction comprises three components: preprocessing, tensor feature learning, and performance analysis.

Cardiac MRI Preprocessing: The preprocessing of cardiac MRI contains ($1$) normalization of scans, ($2$) automatic landmark detection, ($3$) inter-subject registration, and ($4$) in-plane downsampling. We standardize cardiac MRI intensity levels using Z-score normalization [7] to eliminate inter-subject variations. Furthermore, we detect automatic landmarks which is explained in the next paragraph. We perform affine registration to align the heart regions of different subjects to a target image space. We then carry out in-plane scaling of scans by max-pooling at $2$, $4$, $8$, and $16$ times and obtain down-sampled resolutions of $128\times 128$, $64\times 64$, $32\times 32$, and $16\times 16$, respectively.

Landmark Detection and Uncertainty-based Sample Binning: We utilize supervised learning to automate landmark detection using an ensemble of Convolutional Neural Networks (CNNs) for each modality (short-axis and four-chamber). We use the U-Net-like architecture and utilize the same training regime implemented in [15]. We employ Ensemble Maximum Heatmap Activation (E-MHA) strategy [15] which incorporates an ensemble of five models for each modality. We utilize three landmarks for each modality, with the short-axis modality using the inferior hinge point, superior hinge point, and inferolateral inflection point of the right ventricular apex, and the four-chamber modality using the left ventricular apex and mitral and tricuspid annulus. E-MHA produces an associated uncertainty estimate for each landmark prediction, representing the model’s epistemic uncertainty as a continuous scalar value.

A minor error in landmark prediction can result in incorrect image registration [16]. To address this issue, we hypothesize that incorrectly preprocessed samples resulting from inaccurate landmarks can introduce ambiguity during model training. For quality control, it is crucial to identify and effectively handle such samples. In this study, we leverage predicted landmarks and epistemic uncertainties to tackle this problem using uncertainty-based binning. To this end, we partition the training scans based on the uncertainty values of the landmarks. The predicted landmarks are divided into $K$ quantiles, i.e., $Q=\{q_{1},q_{2},...,q_{K}\}$, based on the epistemic uncertainty values. We then iteratively filter out training samples starting from the highest uncertain quantile. A sample is discarded if the uncertainty of any of its landmarks lies in quantile $q_{k}$ where $k=\{1,2,...,K\}$. The samples are discarded iteratively until there is no improvement in the validation performance, as measured by the area under the curve (AUC), for two subsequent iterations.

Tensor Feature Learning: To extract features from processed cardiac scans, we employ tensor feature learning, i.e. Multilinear Principal Component Analysis (MPCA) [11], which learns multilinear bases from cardiac MRI stacks to obtain low-dimensional features for prediction. Suppose we have $M$ scans as third-order tensors in the form of $\{\mathcal{X}_{1},\mathcal{X}_{2},..,\mathcal{X}_{M}\in\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}\}$. The low-dimensional tensor features $\{\mathcal{Y}_{1},\mathcal{Y}_{2},..,\mathcal{Y}_{M}\in\mathbb{R}^{P_{1}\times P_{2}\times P_{3}}\}$ are extracted by learning three ($N=3$) projection matrices $\{U^{(n)}\in\mathbb{R}^{I_{n}\times P_{n}},n=1,2,3\}$ as follows:

$$\mathcal{Y}_{m}=\mathcal{X}_{m}\times_{1}U^{(1)^{T}}\times_{2}U^{(2)^{T}}\times_{3}U^{(3)^{T}},m=1,2,...,M,$$ (1)

where $P_{n}<I_{n}$, and $\times_{n}$ denotes a mode-wise product. Therefore, the feature dimensions are reduced from $I_{1}\times I_{2}\times I_{3}$ to $P_{1}\times P_{2}\times P_{3}$. We optimize the projection matrices $\{U^{(n)}\}$ by maximizing total scatter $\psi_{\mathcal{Y}}=\sum_{m=1}^{M}||\mathcal{Y}_{m}-\bar{\mathcal{Y}}||_{F}^{2}$, where $\bar{\mathcal{Y}}=\frac{1}{M}\sum_{m=1}^{M}\mathcal{Y}_{m}$ is the mean tensor feature and $||.||_{F}$ is the Frobenius norm [10]. We solve this problem using an iterative projection method. In MPCA, $\{P_{1},P_{2},P_{3}\}$ can be determined by the explained variance ratio, which is a hyperparameter. Furthermore, we apply Fisher discriminant analysis to select the most significant features based on their Fisher score [8]. We select the top $k$-ranked features and employ Support Vector Machine (SVM) for classification.

Multimodal Feature Integration: To enhance performance, we perform multimodal feature integration using features extracted from the short-axis, four-chamber, and Cardiac Measurements (CM). We adopt two strategies for feature integration, namely the early and late fusion of features [6]. In early fusion, the features are fused at the input level without doing any transformation. We concatenate features from the short-axis and four-chamber to perform this fusion. We then apply MPCA [11] on the concatenated tensor, enabling the selection of multimodal features. In late fusion, the integration of features is performed at the common latent space that allows the fusion of features that have different dimensionalities. In this way, we can perform a late fusion of CM features with short-axis and four-chamber features. However, we can not perform an early fusion of CM features with short-axis and four-chamber features.

Performance Evaluation: In this paper, we use three primary metrics: Area Under Curve (AUC), accuracy, and Matthew’s Correlation Coefficient (MCC), to evaluate the performance of the proposed pipeline. Decision Curve Analysis (DCA) is also conducted to demonstrate the clinical utility of our methodology.

Table 1: Baseline characteristics of included patients. $p$ values were obtained using $t$-test [20].

	Low PAWP($\leq 15$)	High PAWP($>15$)	$p$-value
Number of patients	$940$	$406$	-
Age (in years)	$64.8\pm 14.2$	$70.5\pm 10.6$	$<0.01$
Body Surface Area (BSA)	$1.88\pm 0.28$	$1.93\pm 0.24$	$<0.01$
Heart Rate (bpm)	$73.9\pm 15.5$	$67.6\pm 15.9$	$<0.01$
Left Ventricle Mass (LVM)	$92.3\pm 25$	$106\pm 33.1$	$<0.01$
Left Atrial Volume ($ml^{2}$)	$72.2\pm 33.7$	$132.2\pm 56.7$	$<0.01$
PAWP (mmHg)	$10.3\pm 3.1$	$21.7\pm 4.96$	$<0.01$

3 Experimental Results and Analysis

Study Population: Patients with suspected pulmonary hypertension were identified after institutional review board approval and ethics committee review. A total of $1346$ patients who underwent Right Heart Catheterization (RHC) and cardiac MRI scans within $24$ hours were included. Of these patients, $940$ had normal PAWP ($\leq 15$ mmHg), while $406$ had elevated PAWP ($>15$ mmHg). Table 1 summarizes baseline patient characteristics. RHC was performed using a balloon-tipped $7.5$ French thermodilution catheter.

Cardiac MRI and measurement: MRI scans were obtained using a $1.5$ Tesla whole-body GE HDx MRI scanner (GE Healthcare, Milwaukee, USA) equipped with $8$-channel cardiac coils and retrospective electrocardiogram gating. Two cardiac MRI protocols, short-axis and four-chamber, were employed, following standard clinical protocols to acquire cardiac-gated multi-slice steady-state sequences with a slice thickness of $8$ mm, a field of view of $48\times 43.2$, a matrix size of $512\times 512$, a bandwidth of $125$ kHz, and TR/TE of $3.7/1.6$ ms. Following  [5], left ventricle mass and left atrial volume were selected as cardiac measurements.

Experimental Design: We conducted experiments on short-axis and four-chamber scans across four scales. To determine the optimal parameters, we performed $10$-fold cross-validation on the training set. From MPCA, we selected the top $210$ features. We employed early and late fusion on short-axis and four-chamber scans, respectively, while CM features were only fused using the late fusion strategy. We divided the data into a training set of $1081$ cases and a testing set of $265$ cases. To simulate a real testing scenario, we designed the experiments such that patients diagnosed in the early years were part of the training set, while patients diagnosed in recent years were part of the testing set. We also partitioned the test into $5$ parts based on the diagnosis time to perform different runs of methods and report standard deviations of methods in comparison results. For SVM, we selected the optimal hyper-parameters from $\{0.001,0.01,0.1,1\}$ using the grid search technique. The code for the experiments has been implemented in Python (version $3.9$). We leveraged the cardiac MRI preprocessing pipeline and MPCA from the Python library PyKale [9] and SVM implementation is taken from scikit-learn [12].

Figure 2: Performance comparison of removing a different number of bins of training data on $10$-fold cross-validation.

(a) Short-axis versus short-axis and CM

(b) Four-chamber versus four-chamber and CM

Figure 3: The effect of combining CM features on short-axis and four-chamber. SA: Short-axis; FC: Four-chamber.

Table 2: Performance comparison using three metrics (with best in bold and second best underlined). FC: Four-Chamber features; SA: Short-Axis features; CM: Cardiac Measurement features. The standard deviations of methods were obtained by dividing the test set into $5$ parts based on the diagnosis time.

Modality	Resolution	AUC	Accuracy	MCC
Unimodal (CM) [5]	-	$0.7300\pm 0.04$	$0.7400\pm 0.03$	$0.1182\pm 0.03$
Unimodal (SA) [17]	$64\times 64$	$0.7391\pm 0.05$	$0.7312\pm 0.07$	$0.3604\pm 0.02$
	$128\times 128$	$0.7495\pm 0.05$	$0.7321\pm 0.04$	$0.3277\pm 0.01$
Unimodal (FC) [17]	$64\times 64$	$0.8034\pm 0.02$	$0.7509\pm 0.04$	$0.4240\pm 0.02$
	$128\times 128$	$0.8100\pm 0.04$	$0.7925\pm 0.05$	$0.4666\pm 0.02$
Bi-modal (SA and FC):	$64\times 64$	$0.7998\pm 0.01$	$0.7698\pm 0.03$	$0.4185\pm 0.03$
Early fusion	$128\times 128$	$0.7470\pm 0.02$	$0.7283\pm 0.02$	$0.3512\pm 0.02$
Bi-modal (SA and FC):	$64\times 64$	$0.8028\pm 0.04$	$0.7509\pm 0.03$	$0.3644\pm 0.01$
Late fusion	$128\times 128$	$0.8122\pm 0.03$	$0.7547\pm 0.03$	$0.3594\pm 0.02$
Bi-modal (SA and CM):	$64\times 64$	$0.7564\pm 0.04$	$0.7585\pm 0.02$	$0.3825\pm 0.02$
Late fusion	$128\times 128$	$0.7629\pm 0.03$	$0.7434\pm 0.03$	$0.3666\pm 0.03$
Bi-modal (FC and CM):	$64\times 64$	$0.8061\pm 0.03$	$0.7709\pm 0.02$	$0.4435\pm 0.02$
Late fusion	$128\times 128$	$0.8135\pm 0.02$	$0.7925\pm 0.02$	$0.4999\pm 0.03$
Tri-modal (FC, SA, and CM)	$64\times 64$	$0.8146\pm 0.04$	$0.7774\pm 0.03$	$0.4460\pm 0.02$
Hybrid fusion	$128\times 128$	$\mathbf{0.8327\pm 0.06}$	$\mathbf{0.8038\pm 0.05}$	$\mathbf{0.5099\pm 0.04}$
Tri-modal Hybrid fusion	$64\times 64$	$0.7892\pm 0.04$	$0.7513\pm 0.05$	$0.4278\pm 0.02$
without uncertainty binning	$128\times 128$	$0.8036\pm 0.03$	$0.7820\pm 0.04$	$0.4779\pm 0.01$

Uncertainty-Based Sample Binning: To improve the quality of training data, we used quantile binning to remove training samples with uncertain landmarks. The landmarks were divided into $50$ bins, and then removed one bin at a time in the descending order of their uncertainties. Figure 2 depicts the results of binning using $10$-fold cross-validation on the training set, where the performance improves consistently over the four scales when removed bins $\leq 5$. Based on the results, we removed $5$ bins ($129$ out of $1081$ samples) from the training set, and used the remaining $952$ training samples for the following experiments.

Unimodal Study: The performance of three models on single-modality is reported in Table 2, including short-axis (SA), four-chamber (FC), and cardiac measurements (CM), where the CM based unimodal is considered as the baseline. The results demonstrate an improvement of $\Delta$AUC $=0.0800$ $\Delta$Accuracy $=0.0527$, and $\Delta$MCC $=0.3484$ over the baseline obtained by FC based unimodal, which indicates that tensor-based features have a diagnostic value.

Bi-modal Study: In this experiment, we compared the performance of bi-modal models. As shown in Table 2, bimodal (four-chamber and CM) produces superior performance (i.e., AUC $=0.8135$, Accuracy=$0.7925$ and MCC $=0.4999$) among bi-modal models. Next, we investigated the effect of fusing CM features with short-axis and four-chamber modalities in Fig. 3. It can be observed from these figures that the fusion of CM features enhances the diagnostic power of cardiac MRI modalities at all scales. The bi-modal (four-chamber and CM) model achieved the improvement in the performance ($\Delta$AUC $=0.0035$ and $\Delta$MCC $=0.0333$) over the unimodal (four-chamber) model.

Effectiveness of Tri-modal: In this experiment, we performed a fusion of CM features with the bi-modal models to create two tri-modal models. The first tri-modal is tri-modal late (CM with a late fusion of short-axis and four-chamber) and the second tri-modal is a tri-modal hybrid (CM with an early fusion of short-axis and four-chamber). As shown in Fig. 4, CM features enhance the performance of bi-modal models and tri-modal hybrid outperforms all. The tri-modal hybrid obtained the best performance (Table 2, where AUC $=0.8327$, Accuracy $=0.8038$, and MCC $=0.5099$) and a significant improvement of $\Delta$AUC $=0.1027$, $\Delta$Accuracy $=0.0628$, and $\Delta$MCC $=0.3917$ over the baseline method.

Figure 4: The effect of combining CM features on the bi-modals including early and late fusion of four-chamber and short-axis. Early fusion: early fusion of short-axis and four-chamber; late fusion: late fusion of short-axis and four-chamber.

Figure 5: Evaluating clinical utility of our method using Decision Curve Analysis (DCA) [19].“Treat All” means treating all patients, regardless of their actual disease status, while “Treat None” means treating no patients at all. Our predictive model’s net benefit is compared with the net benefit of treating everyone or no one to determine its overall utility.

Decision Curve Analysis (DCA) [19, 13] on the performance suggests the potential clinical utility of the proposed method. As shown in Fig. 5, the Tri-modal model outperformed the baseline method for most possible benefit/harm preferences, where benefit indicates a positive net benefit (i.e. correct diagnosis) and harm indicates a negative net benefit (i.e. incorrect diagnosis). The tri-modal model (the best model) obtained a higher net benefit between decision threshold probabilities of $0.30$ and $0.70$ which implies that our method has a diagnostic value and can be used in screening high-risk patients from a large population.

Feature contributions: Our model is interpretable. The highly-weighted features were detected in the left ventricle and interventricular septum in cardiac MRI. For cardiac measurements, left atrial volume (0.778/1) contributed more than left ventricular mass (0.222/1) to the prediction.

4 Conclusions

This paper proposed a tensor learning-based pipeline for PAWP classification. We demonstrated that: $1$) tensor-based features have a diagnostic value for PAWP, $2$) the integration of CM features improved the performance of unimodal and bi-modal methods, $3$) the pipeline can be used to screen a large population, as shown using decision curve analysis. However, the current study is limited to single institutional data. In the future, we would like to explore the applicability of the method for multi-institutional data using domain adaptation techniques.