We collected the data of elderly HF patients with malnutrition (female or male $\ge$65 years) who were admitted to Bengbu Central Hospital (BCH) between October 2022 to October 2024 and the first affiliated hospital of Bengbu medical university (BMU 1st AH) between January 2021 to October 2024. The inclusion criteria were as follows: (1) patients with a confirmed diagnosis of HF; (2) patients had presence of nutritional impairment evidenced by the GNRI (score $\le$98). The exclusion criteria were as follows: (1) history of acute myocardial infarction (AMI) (within 3 months); (2) active malignancy; (3) end-stage renal diseases on maintenance dialysis; (4) decompensated cirrhosis (Child-Pugh score $\ge$7). Following screening, 889 patients from BMU 1st AH and 191 patients from BCH were analyzed (Fig. 1). Laboratory data were collected within 24 hours of admission for non-transferred patients, or within 24 hours post-transfer for those moved between departments. The primary endpoint was all-cause in-hospital mortality, assessed as a binary outcome (survival or death) at discharge. This study received ethical approval from the ethics committee of the first affiliated hospital of Bengbu medical university, China (No.2025-551), and was granted an exemption by Bengbu central hospital. Due to the retrospective observational nature of this study, the requirement for informed consent was waived.

Fig. 1.

A flowchart illustrating the data screening process of this study. HF, heart failure; AMI, acute myocardial infarction; BMU 1st AH, first affiliated hospital of bengbu medical university; BCH, bengbu central hospital.

Guided by the underlying causes of HF in geriatric patients, established clinical practice, and contemporary literature, we selected 62 variables for this study. The primary aim of this study was to develop a mortality prediction model. Accordingly, we collected binary outcome data (death or survival) for all patients. In addition, we collected demographic characteristics including sex, age, weight, height and body mass index (BMI) was calculated for analytical purposes. Comorbidities such as diabetes mellitus (DM), coronary artery disease (CAD), arrhythmia, valvular heart disease (VHD), hypertension (HTN), pulmonary diseases, urinary tract infection (UTI), cerebral infarction and so on. Our study evaluated patients’ nutritional status using the GNRI, calculated as:

GNRI = [1.489 $\times$ albumin (g/dL)] + [41.7 $\times$ (current weight / ideal weight)].

Male ideal weight (kg) = height (cm) – 100 – [height (cm) – 150] / 4.

Female ideal weight (kg) = height (cm) – 100 – [height (cm) – 150] / 2.5.

We assessed nutritional risk using the GNRI, including all patients with a score $\le$98 (indicative of malnutrition risk) and recording their specific values.

Therefore, we recorded serum albumin (ALB) levels at both admission and discharge, along with the total dosage of human albumin administered during hospitalization. NT-proBNP and left ventricular ejection fraction (LVEF) were systematically included in our study which as gold-standard diagnostic markers for HF. Additionally, we collected the following admission laboratory parameters: white blood cell (WBC), neutrophil percentage (NEUTperc), lymphocyte percentage (LYperc), hemoglobin (Hb), platelet (PLT), creatinine (Cr), blood urea nitrogen (BUN), aspartate aminotransferase (AST), alanine aminotransferase (ALT), cholinesterase (AchE), electrolyte levels, C-reactive protein (CRP), total cholesterol (TC), triglyceride (TG), low-density lipoprotein cholesterol (LDL-C), high-density lipoprotein cholesterol (HDL-C), cardiac troponin I (cTnI), creatine kinase-MB (CK-MB) and so on. And we calculated six novel inflammatory markers for model construction: NLR, PLR, systemic immune inflammation index (SII), neutrophil percentage to albumin ratio (NPAR), lymphocyte to C-reactive protein ratio (LCR), Platelet to high density lipoprotein cholesterol ratio (PHR).

(1) NLR = Neutrophil count ($\times$10⁹/L) / Lymphocyte count ($\times$10⁹/L).

(2) PLR = Platelet count ($\times$10⁹/L) / Lymphocyte count ($\times$10⁹/L).

(3) SII = (Platelet count $\times$ Neutrophil count) / Lymphocyte count (all ($\times$10⁹/L)).

(4) NPAR = Neutrophil percentage (%) / Serum albumin (g/dL).

(5) LCR = Lymphocyte count ($\times$10⁹/L) / C-reactive protein (mg/L).

(6) PHR = Platelet count ($\times$10⁹/L) / High-density lipoprotein cholesterol (mmol/L).

Due to the inability to perform echocardiography in some acutely ill HF patients during hospitalization, the LVEF data had a significant missing rate of 11.12%. To ensure the robustness of statistical inferences, multiple imputation was performed using the mice package in R to generate 5 complete datasets. For the purpose of describing baseline characteristics in Table 1, the first imputed database was used, as is the conventional practice for descriptive statistics. And due to the strong correlations among some variables, the study employed the least absolute shrinkage and selection operator (LASSO) technique to pick out key clinical factors, while simultaneously discarding extraneous data.

In this research, 889 patients from the BMU 1st AH were divided into training and validation cohorts at a 7:3 ratio through stratified sampling, while 191 patients from BCH as the testing cohort. The training cohort, which constitutes the largest proportion of the dataset, is used for model learning and parameter optimization. The validation cohort monitors the training process to prevent overfitting and helps fine-tune model performance. Finally, a completely independent test cohort evaluates the model’s final performance after training and parameter adjustment are completed. The XGBoost model achieved higher predictive accuracy than other machine learning models in forecasting cardiovascular events based on laboratory test results [16, 17]. Given that the primary data we collected consisted of laboratory test results and that HF was either the primary or secondary cause of death in deceased patients, we choose to construct an XGBoost model to predict patient mortality. Hyperparameter tuning was performed using a systematic grid search approach via the train function from the caret package in R. The objective was to optimize the performance of the XGBoost model, with the primary evaluation metric set to the logarithmic loss (logloss) calculated through cross-validation. After an exhaustive search within the defined parameter space, the final optimal configuration was determined as follows: a learning rate of 0.002, a maximum tree depth of 2, 200 boosting rounds, a gamma value of 0, a column sampling rate of 0.6, a minimum child weight of 1, and a subsample ratio of 0.8. This specific combination of parameters was selected as it yielded the lowest cross-validated logloss, thereby providing the model with the best generalizability on the given dataset within the explored search space. Model performance was evaluated using the receiver operating characteristic (ROC) curve, calibration curve, and decision curve analysis (DCA). Additionally, metrics including accuracy, sensitivity, and specificity were calculated and reported to comprehensively assess the predictive performance of the model. To address the inherent opacity ML algorithms and improve model interpretability, Shapley Additive Explanations (SHAP) values were applied for feature importance analysis. A bar plot was generated to visualize the contribution of each feature to the model’s predictions.

All statistical analyses and data visualizations were performed using R4.3.0 (R Foundation for Statistical Computing, Vienna, Austria). Continuous variables were presented as mean $\pm$ standard deviation (normally distributed data) or median (IQR) (non-normally distributed data), while categorical variables were expressed as frequencies (percentage). For between-group comparisons, Chi-square tests were used for categorical variables, and independent-sample t-tests (normal distribution) or nonparametric tests (non-normal distribution) were applied for continuous variables, as appropriate.

Complete blood count analysis was performed at BMU 1st AH (Bengbu, Anhui, China) using a Sysmex XN-9000 hematology analyzer (Sysmex Corporation, Kobe, Japan), and at BCH (Bengbu, Anhui, China) using a Mindray BC-7500 hematology analyzer (Shenzhen Mindray Bio-Medical Electronics Co., Ltd., Shenzhen, China). Biochemical testing at BMU 1st AH was conducted on a Cobas c 701 analyzer (Roche Diagnostics, Basel, Switzerland). Similarly, biochemical testing at Bengbu Central Hospital was independently performed using a Cobas c 701 analyzer of the same model (Roche Diagnostics, Basel, Switzerland).

The study flowchart is shown in Fig. 1. After applying the inclusion and exclusion criteria, clinical data from 889 patients were ultimately utilized to generate both the training set and validation set. Table 1 presents the baseline characteristics of these patients, primarily including demographic information, laboratory test results and comorbidities. Among the 889 elderly patients (median age was 78 years) with HF and malnutrition in the training and validation sets, 564 patients (63.44%) had CAD, 520 patients (58.49%) had arrhythmia, 446 patients (50.16%) presented with pulmonary infection upon admission and 437 patients (49.15%) had HTN. We employed stratified sampling for data partitioning, ensuring balanced baseline characteristics between the training and validation sets. With the addition of 191 cases from BCH for external validation, the total study population comprised 836 survivors and 244 deaths. The endpoint event rates in the training set, validation set and testing set were 22.953%, 22.932% and 20.942%, respectively.

Pairwise correlations among all 62 candidate variables were calculated using spearman methods with significance levels (p-values) adjusted for multiple comparisons via the Benjamini-Hochberg procedure. Results were visualized in a correlation heatmap (Fig. 2), revealing significant collinearity among multiple biomarker pairs. To mitigate multicollinearity effects on model stability, we employed iterative LASSO with 10-fold cross-validation [18]. The optimized hyperparameter $\lambda$ was identified based on the bivariate deviation. Based on the optimal $\lambda$, 13 predictors are identified from the initial pool of 62 candidate variables, as seen in Fig. 3, in order: GNRI, NEUTperc, LYperc, NT-proBNP, BUN, NLR, ANC, serum potassium (K), PLR, serum phosphorus (P), PHR, Hyponatremia, and AST.

Fig. 2.

Variable importance and correlation analysis. (A) AUC values of each variable. (B) Heatmap of correlations among variables, with more intense colors indicating higher correlation. *, **, ***, and **** indicate statistical significance at different levels, with more asterisks representing stronger significance.

Fig. 3.

LASSO regression for variable selection. (A) LASSO coefficient path plot. (B) LASSO cross-validation curve, the vertical dashed lines indicate the optimal $\lambda$ values where the model achieves improved predictive performance while retaining the corresponding variables. LASSO, least absolute shrinkage and selection operator.

The XGBoost model demonstrated robust performance in predicting mortality among elderly patients with HF and malnutrition. Key metrics across training, validation and test datasets are summarized in Table 2 and Fig. 4. The model exhibited outstanding discriminative ability (training AUC = 0.979, testing AUC = 0.936). While specificity remained high ($>$94% across all sets), sensitivity dropped significantly in validation (52.4%) and test sets (58.5%). The F1-score reflected this imbalance (training: 0.877 vs testing: 0.675), suggesting potential under-detection of true mortality cases in external data.

Fig. 4.

ROC curves of the predictive model. (A) Training set. (B) Validation set. (C) Testing set. ROC, receiver operating characteristic.

The SHAP summary plot demonstrates the relative contributions of each feature to the predictive model, with variables ranked by their mean absolute SHAP values in descending order: GNRI, NEUTperc, LYperc, NT-proBNP, BUN, NLR, K, PLR, P, PHR, hyponatremia, and AST. Higher SHAP values for a given feature indicate greater association with increased mortality risk. To enhance interpretability, we generated a beeswarm plot (Fig. 5) that visually distinguishes features with positive versus negative impacts on mortality prediction: yellow data points (higher feature values) generally represent risk-increasing effects, while purple (lower values) shows protective effects, with the intensity of color reflecting the magnitude of SHAP value contribution.

Fig. 5.

Feature importance and impact analysis. (A) Feature importance plot. Features are ranked by mean absolute SHAP values in descending order. (B) Beeswarm plot of significant variables, demonstrating the impact of each feature on model output. Each dot represents a patient, with color indicating feature value (yellow: higher values; purple: lower values). Increased dot dispersion suggests greater influence on predictions. SHAP, shapley additive explanation.

This study constructed calibration curves and clinical decision curves based on an independent test set (Fig. 6). The logistic calibration curve exhibited a non-parallel deviation, with a calibration intercept of 1.064, indicating an overall systematic underestimation of predicted probabilities relative to the actual observed probabilities. The calibration slope was 1.358, suggesting that the model was overconfident in distinguishing between high-risk and low-risk groups. The model demonstrated good discriminative ability, with a Dxy value of 0.872, reflecting its effectiveness in risk stratification between cases and non-cases. The overall prediction error, as measured by the Brier score, was 0.113, which falls within an acceptable range and was primarily attributed to calibration bias. Furthermore, the nonparametric calibration curve showed a parallel deviation, indicating consistency between the risk ranking derived from the model and the observed risks. In summary, the model exhibited satisfactory discriminative performance and interpretability, however, its suboptimal calibration limits the reliability for direct application, and further validation or calibration refinement is warranted to enhance predictive accuracy.

Fig. 6.

Model calibration and clinical utility evaluation. (A) Calibration curve plotted on the testing set to evaluate the agreement between predicted and observed probabilities. The gray solid line represents the ideal reference where predicted probability equals observed probability. The black solid line indicates the logistic calibration fit, reflecting the overall calibration trend of the model. The dotted line depicts the nonparametric calibration curve, capturing local deviations from ideal calibration. (B) Clinical decision curve based on the testing set. The x-axis is the risk threshold, and the y-axis is the benefit. The blue line represents intervening only when the model predicts risk above the threshold, the red line represents treating all patients, and the green line represents treating none. Where the blue line lies above the other two, using the model provides the greatest net benefit.

The decision curve analysis revealed that across a broad threshold probability range of approximately 10% to 75%, the net clinical benefit of using this model was higher than that of the two empirical strategies: treat all and treat none. These results support the model’s utility as an effective tool for clinical decision-making.