A retrospective cohort study design was adopted. A total of 346 patients who underwent HVR at our cardiovascular disease center from January 2019 to June 2024 were continuously included. Inclusion criteria were as follows: first-time HVR surgery (mechanical or biological valve) [8], age ≥18 years, types of surgery included isolated aortic valve replacement, isolated mitral valve replacement, combined aortic and mitral valve replacement, and valve replacement combined with coronary artery bypass grafting, elective or sub-emergency surgery under cardiopulmonary bypass at our center, complete and retrievable clinical medical records. Exclusion criteria were as follows: patients who underwent transcatheter aortic valve implantation, patients complicated with acute infective endocarditis, patients complicated with acute myocardial infarction requiring concurrent surgery, patients with severe missing clinical data (>20%). The flow of patient screening and cohort formation was summarized in Supplementary Fig. 1.

Patient data were collected through the electronic medical record system, including: preoperative data: age, gender, body mass index, hypertension, diabetes, history of previous cardiac surgery, etc., preoperative left ventricular ejection fraction (LVEF), New York Heart Association (NYHA) functional classification (class III/IV was defined as poor cardiac function), EuroSCORE II, hemoglobin, serum albumin, creatinine, C-reactive protein, total bilirubin, platelet count. Intraoperative data: type of surgery, whether combined with coronary artery bypass grafting, cardiopulmonary bypass time, aortic cross-clamp time, valve type (mechanical/biological valve). The primary outcome was the occurrence of major adverse cardiovascular events (MACE, a composite endpoint) in the hospital after surgery, defined as the incidence of any of the following events within 30 days postoperatively: all-cause death, low cardiac output syndrome, new-onset stroke, reoperation due to bleeding, acute kidney injury requiring renal replacement therapy.

According to the international guidelines of cardiothoracic surgery, cases were grouped according to whether postoperative major complications occurred [9]. The complication group (MACE group) was defined as patients who experienced major adverse cardiovascular events (including all-cause death, low cardiac output syndrome, new-onset stroke, reoperation due to bleeding, and acute kidney injury requiring renal replacement therapy) within 30 days after valve replacement. The non-complication group (Non-MACE group) was the control group of patients who did not experience the above events during the same period. The 30-day follow-up period was chosen to capture short-term surgical outcomes consistently.

The sample size for this study was estimated based on the Events Per Variable (EPV) principle to ensure sufficient outcome events and mitigate model overfitting. With an estimated in-hospital MACE incidence of approximately 22% and planning for 5–6 final predictor variables, a minimum of 25–30 events were required. Our final cohort of 346 patients, yielding 76 MACE events, satisfied this requirement.

Data analysis was performed using SPSS 26.0 (IBM Corp., Armonk, NY, USA), R 4.2.3 (R Foundation for Statistical Computing, Vienna, Austria), and Python 3.8.5 (Python Software Foundation, Beaverton, OR, USA). Continuous variables are presented as mean ± standard deviation or median (interquartile range), and categorical variables as counts (percentages). Group comparisons used Student’s t-test, Mann-Whitney U test, or Chi-square test as appropriate. Before model construction, comprehensive data preprocessing was performed: (1) Missing data: Variables with missing rates >20% were excluded per protocol. For remaining variables, missing values (<20%) were imputed using multiple imputation (5 iterations). (2) Outliers: Continuous variable outliers were identified using the 3σ rule and reviewed for clinical plausibility; non-physiological extreme values were winsorized at the 99th percentile. (3) Standardization: All continuous variables were standardized via z-score normalization to neutralize scale effects. The dataset was initially split into a training set (70%, n = 242) and a validation set (30%, n = 104) using stratified random sampling to preserve the MACE incidence ratio. In the training set, univariate analysis was first conducted to screen for variables. Least absolute shrinkage and selection operator (LASSO) regression was applied to the significant variables from the univariate analysis for further feature selection. The optimal penalty coefficient (λ) was determined via 10-fold cross-validation to select the core variables with the highest predictive value. Statistically significant variables were included in a multivariate logistic regression analysis to identify independent influencing factors. Prediction models were constructed using random forest (RF) Model, K-Nearest Neighbors (K Model), and gradient boosting (GB Model) algorithms implemented in Python (v3.8.5) with the scikit-learn library. Hyperparameters for each model were optimized via a 5-fold cross-validated grid search on the training set. Key tuned parameters included: for the RF Model (n_estimators, max_depth, min_samples_split, min_samples_leaf); for the K Model (n_neighbors); and for the GB Model (n_estimators, max_depth, learning_rate). Model performance was primarily evaluated by the area under the receiver operating characteristic (ROC) curve (AUC). The Delong test was used for statistical comparison of AUCs between models and against the EuroSCORE II. Calibration was assessed with calibration curves and the Brier score. To rigorously evaluate model stability and optimism, 5-fold cross-validation was applied during training, and bootstrap resampling (1000 iterations) was performed on the entire dataset to obtain optimism-corrected performance estimates. Decision curve analysis (DCA) was used to evaluate the clinical net benefit. To enhance interpretability of the optimal model, SHapley Additive exPlanations (SHAP) analysis was conducted. A two-sided p-value < 0.05 was considered statistically significant.

Among the 346 enrolled patients, 76 (22.0%) experienced at least one MACE component within 30 days postoperatively. According to the ratio of 7:3, 346 patients who underwent HVR were divided into a training dataset (including 242 cases) and a validation dataset (including 104 cases). There were no statistically significant differences in the comparison of general data between patients in the training set and the validation set (all p > 0.05) (Table 1).

In the training set, univariate analysis showed that there were statistically significant differences in age, LVEF, EuroSCORE II, serum albumin, cardiopulmonary bypass time, and aortic cross-clamp time between the non-complication group and the complication group (all p < 0.05) (Table 2).

The outcome after heart valve replacement was used as the dependent variable (non-complication group = 0, complication group = 1) (Supplementary Table 1). The indicators with statistical significance in the univariate analysis were included in the LASSO regression for variable screening. Variables were selected using the screening criterion of lambda.1se (Supplementary Figs. 2,3). The appropriate predictive variables were age, LVEF, EuroSCORE II, serum albumin, cardiopulmonary bypass time, and aortic cross-clamp time. The result of multivariate logistic regression analysis identified age, EuroSCORE II, cardiopulmonary bypass time, aortic cross-clamp time, LVEF, and serum albumin as independent influencing factors for MACE (all p < 0.05) (Fig. 1).

Fig. 1.

Forest plot of the multivariable logistic regression analysis of influencing factors for MACEs after heart valve replacement. X1, age; X2, left ventricular ejection fraction; X3, EuroSCORE II; X4, serum albumin; X5, cardiopulmonary bypass time; X6, aortic cross-clamp time; MACEs, major adverse cardiovascular events.

The RF Model, K Model, and GB Model were used for prediction in the training set and validation set. The AUC values of the three models in the training set were 0.847 (95% CI: 0.785–0.908), 0.837 (95% CI: 0.776–0.898), and 0.835 (95% CI: 0.771–0.898) respectively, and those in the validation set were 0.823 (95% CI: 0.715–0.930), 0.790 (95% CI: 0.672–0.908), and 0.771 (95% CI: 0.646–0.896) respectively. The model with the largest AUC value was selected as the best model for this study, which was the RF Model (Fig. 2). The predictive performance of the ML models and the EuroSCORE II was summarized in Table 3. The results of the DeLong test indicated that in the validation set, the AUC of the RF Model was significantly higher than that of the EuroSCORE II score and superior to the other compared machine learning models. Specifically, the AUC of the RF Model was 0.823, which was significantly higher than that of EuroSCORE II (0.723), K Model (0.790), and GB Model (0.771), with corresponding p-values of 0.028, 0.036, and 0.021, respectively. 5-fold cross-validation results showed that the average AUC of the RF Model was 0.839 (95% CI: 0.792–0.886), which was consistent with the 70/30 split validation result (0.823), indicating good model stability. Bootstrap correction (1000 iterations) yielded an optimism-corrected AUC of 0.811 (95% CI: 0.755–0.867) for the RF Model, which remains substantially higher than the EuroSCORE II, indicating robust performance despite the sample size constraints.

Fig. 2.

Area under the receiver operating characteristic curve of machine learning models. (A) The training set. (B) The validation set. K Model, K-Nearest Neighbors Model; GB Model, gradient boosting Model; AUC, area under the receiver operating characteristic curve; TPR, true positive rate.

Furthermore, calibration analysis was performed for the RF Model. The Brier score for the validation set was 0.124, where a lower score indicates better model calibration. The calibration curve demonstrated good agreement between the predicted probability of MACE and the actual observed frequency, without significant systematic overestimation or underestimation (Fig. 3).

Fig. 3.

Calibration curves of machine learning models. (A) Training set. (B) Validation set.

DCA showed that across a wide range of clinical decision thresholds (high-risk probability threshold approximately 0.1–0.5), the net clinical benefit curve of the RF Model was consistently higher than the “treat all” and “treat none” strategies, indicating that using it to guide clinical decisions (e.g., intensifying interventions for high-risk patients) could yield greater net clinical benefit (Fig. 4). The net benefits of the K Model and the GB Model were also superior to the simple strategies but lower than those of the RF Model. Compared with the currently used EuroSCORE II score, the RF Model demonstrated significant reclassification improvement, with a net reclassification improvement index of 0.250 (95% CI: 0.070–0.430, p = 0.006) and an integrated discrimination improvement index of 0.080 (95% CI: 0.030–0.130, p = 0.001). After stratifying patients into low-, medium-, and high-risk groups based on the model’s predicted probabilities, the actual incidence of MACE in the high-risk group (45.2%) was significantly higher than that in the medium-risk (17.3%) and low-risk (3.8%) groups (p for trend <0.001), confirming the model’s good risk stratification capability.

Fig. 4.

Decision curves of machine learning models. (A) Training set. (B) Validation set.

As the number of decision trees increased, the overall error gradually stabilized. This change trend reflected the dynamic change characteristics of the predictive performance of the model during the iterative construction of decision trees. This trend can be used to assist in judging the convergence of the model. When the error curve tended to be stable, it indicated that the model complexity reached a certain level, and the optimization effect of new decision trees on the error was limited. This provided a basis for parameter selection to determine the optimal number of decision trees and improve the predictive efficiency of the model, helping to screen out the configuration that could balance model complexity and prediction accuracy, so as to enhance the model’s performance in predicting the outcomes after heart valve replacement (Supplementary Fig. 4). Based on the RF Model, the importance scores of independent influencing factors for the outcomes after heart valve replacement were calculated. The importance ranking was as follows: EuroSCORE II, LVEF, cardiopulmonary bypass time, serum albumin, age, and aortic cross-clamp time (Fig. 5).

Fig. 5.

Importance ranking of the random forest model. X1, age; X2, left ventricular ejection fraction; X3, EuroSCORE II; X4, serum albumin; X5, cardiopulmonary bypass time; X6, aortic cross-clamp time.

SHAP analysis was conducted to interpret the predictions of the optimal RF Model. The global feature importance ranking derived from mean absolute SHAP values was largely consistent with the Gini importance ranking reported earlier, with EuroSCORE II being the most influential predictor, followed by cardiopulmonary bypass time. Analysis of SHAP values showed the directional impact of each predictor: higher EuroSCORE II scores, longer cardiopulmonary bypass time, longer aortic cross-clamp time, and advanced age were associated with increased risk of MACE (positive SHAP values). Conversely, higher LVEF and serum albumin levels demonstrated a protective contribution, associated with a decrease in the model’s predicted risk. Furthermore, the SHAP dependence plots illustrated the non-linear relationships between key continuous predictors (e.g., cardiopulmonary bypass time) and the model’s predicted risk, demonstrating how risk increased with prolonged cardiopulmonary bypass time (Fig. 6).

Fig. 6.

SHapley Additive exPlanations feature importance graph. X1, age; X2, left ventricular ejection fraction; X3, EuroSCORE II; X4, serum albumin; X5, cardiopulmonary bypass time; X6, aortic cross-clamp time.

This study has several limitations. First, as a single-center, retrospective study, the model is derived from a specific patient population, which carries inherent risks of selection bias and may limit generalizability to other settings with different demographics or clinical practices. Although internal validation using cross-validation and bootstrap resampling showed good performance, external validation with multi-center, prospective cohorts is essential before clinical implementation to confirm robustness and universality. Second, the model relies on routinely available clinical variables and did not incorporatepotentially prognostic emerging biomarkers (e.g., N-terminal pro-brain natriuretic peptide (NT-proBNP), troponin) or genetic data. Future studies integrating multi-omics information could enhance predictive power. Third, while the model demonstrated superior discrimination and calibration compared to EuroSCORE II, a formal assessment of its clinical utility. Evaluating the model’s impact on real-world clinical decisions and patient outcomes is a crucial next step for prospective implementation research. Fourth, although we employed SHAP analysis to improve interpretability and reported global feature importance, the intrinsic complexity of the RF Model may still pose a “black-box” challenge for individual case interpretations. Future work could leverage more advanced explainable AI techniques to provide even clearer, patient-specific decision pathways. Finally, while hyperparameter tuning was conducted via grid search with cross-validation, and the performance of the RF Model was validated against other classical algorithms (K Model, GB Model), exploring more advanced ensemble methods (e.g., XGBoost, LightGBM) and automated machine learning (AutoML) frameworks in future research may yield further performance optimizations.