Abstract

Background: Cardiovascular assessment in patients undergoing maintenance hemodialysis is complicated by elevations and variability in cardiac biomarkers. We aimed to develop and temporally validate an interpretable machine learning model for classifying a study-defined composite of elevated B-type natriuretic peptide (BNP) or high-sensitivity cardiac troponin I (hs-cTnI) at the same clinical assessment episode. Methods: This study included 960 adults receiving maintenance hemodialysis in China. Patients treated from January 2021 to December 2023 comprised the development-period cohort, and those treated from January 2024 to December 2025 comprised the temporal validation cohort. The development-period cohort was divided into training and same-period holdout sets. Least absolute shrinkage and selection operator (LASSO) regression selected predictors, and six models were developed. Three top-performing models were combined in a stacking ensemble. Discrimination, grouped calibration, Brier score, decision-curve analysis, and Shapley additive explanations (SHAP) were evaluated. Results: LASSO retained 13 predictors. The stacking ensemble achieved ROC-AUC point estimates of 0.848, 0.863, and 0.835 in training out-of-fold evaluation, same-period holdout evaluation, and temporal validation, respectively. Corresponding precision-recall AUC point estimates were 0.77, 0.76, and 0.70, and Brier scores were 0.15, 0.15, and 0.16. In temporal validation, sensitivity, specificity, and accuracy point estimates were 0.75, 0.76, and 0.75. SHAP identified prognostic nutritional index (PNI), log-transformed D-dimer, age, red cell distribution width (RDW), and beta-2 microglobulin as leading predictors. Small between-model differences should be interpreted cautiously because confidence intervals were not estimated. Conclusions: The stacking ensemble showed comparable point-estimate performance across same-period holdout and temporal validation data and may support further investigation of contemporaneous biomarker-burden classification. Multicenter external validation, standardized sampling, uncertainty quantification, and prospective evaluation are required before clinical implementation.

Keywords

Maintenance Hemodialysis, Cardiac Biomarkers, Machine Learning, Stacking Ensemble, Temporal Validation, Contemporaneous Classification, Explainable Artificial Intelligence

Share and Cite:

1. Introduction

Cardiovascular disease is a major contributor to morbidity and mortality among people with kidney failure treated with dialysis. Conventional in-center hemodialysis adds recurrent hemodynamic fluctuations and rapid changes in volume and solute concentrations to an already complex cardiovascular risk profile [1]. These features make cardiovascular assessment difficult, particularly when symptoms are nonspecific and conventional risk factors do not fully characterize subclinical myocardial stress or injury.

B-type natriuretic peptide (BNP) and high-sensitivity cardiac troponin I (hs-cTnI) represent different aspects of cardiovascular stress. BNP is released in response to myocardial wall stress, whereas hs-cTnI reflects myocardial injury; both have been investigated for cardiovascular risk assessment in chronic kidney disease [2]. Their interpretation in patients receiving dialysis is not straightforward. High-sensitivity troponin concentrations are frequently elevated in clinically stable dialysis patients and show substantial within-person variability, so a single measurement may not reliably distinguish chronic elevation from acute ischemic injury [3]. In a systematic review and meta-analysis of 61 studies in patients with end-stage kidney disease, higher BNP or N-terminal pro-BNP concentrations were associated with greater risks of cardiovascular and all-cause mortality, although the estimated risks varied across thresholds and some estimates were imprecise [4]. Classifying a study-defined state of biomarker elevation is therefore distinct from using either marker alone to diagnose acute myocardial infarction or heart failure.

Machine learning offers a way to model heterogeneous clinical data and potentially nonlinear relations among demographic characteristics, comorbidities, vital signs, and laboratory measurements. Previous studies have applied random forests and gradient-boosting methods to mortality prediction in hemodialysis, including an incident cohort of 1571 patients and the international Monitoring Dialysis Outcomes data set of 95,142 patients [5] [6]. These studies demonstrate feasibility while also showing that model performance and selected predictors depend on the algorithm, population, and assessment target. Evaluation of a new classification target should therefore extend beyond a single discrimination statistic to include calibration, temporal validation, and transparent reporting. Explainable artificial intelligence methods, including Shapley additive explanations (SHAP), can describe the direction and magnitude of variable contributions at the population and individual levels [7]. The TRIPOD + AI statement further emphasizes complete reporting of model development and evaluation regardless of whether regression or machine-learning methods are used [8].

Accordingly, an interpretable classification framework may offer a structured way to integrate multidimensional clinical information when assessing cardiac biomarker elevation in hemodialysis. In this retrospective, single-center model-development and temporal-validation study, we aimed to classify a study-defined composite of elevated cardiac biomarkers at the index assessment episode and to evaluate model performance in a later-period cohort. The intended use is research-stage identification of patients whose routinely recorded clinical profile is consistent with contemporaneous high biomarker burden when biomarker results are unavailable or pending; the model does not predict a future event and is not intended to diagnose heart failure or acute myocardial infarction. We further sought to determine whether explainable classification could support transparent characterization of model output at both the population and individual levels.

2. Methods

2.1. Ethics and Study Design

This study was reviewed and approved by the Ethics Committee of the participating hospital (Approval No. YJ202570) and was conducted in accordance with the Declaration of Helsinki and relevant institutional requirements. Because this retrospective observational study used routinely collected and de-identified clinical data, the requirement for written informed consent was waived in accordance with the decision of the ethics committee. All data were accessed only by authorized investigators, and strict measures were implemented to protect patient privacy and confidentiality.

This was a single-center model-development and temporal-validation study conducted at a tertiary hospital in Jingzhou, Hubei Province, China. Patients treated from January 2021 through December 2023 constituted the development-period cohort, whereas patients treated from January 2024 through December 2025 constituted the temporal validation cohort. The development-period cohort was used for model training and same-period holdout evaluation, and all fitted models were subsequently applied without refitting or retuning to the later-period cohort. Because the later-period data originated from the same hospital, this cohort represents temporal validation rather than external validation in an independent setting. All observations represented actual clinical records; no synthetic data generation, bootstrap-based sample expansion, artificial perturbation of clinical variables, model-guided outcome assignment, or performance-targeted cohort construction was performed.

2.2. Data Sources and Study Population

Demographic characteristics, comorbidities, vital signs, anthropometric measurements, and laboratory findings were obtained from the hospital electronic medical record system, laboratory information system, and hemodialysis management database. Consecutive adult patients receiving maintenance hemodialysis during the study periods were screened. Patients were eligible if they were aged 18 years or older, had received maintenance hemodialysis for at least 3 months, had available BNP and hs-cTnI measurements during the index clinical assessment, and had sufficient information for outcome definition, feature selection, and model development. Patients receiving temporary dialysis for acute kidney injury, those without a definable primary outcome, those lacking variables required for the final model, and those with duplicate or internally inconsistent records that could not be resolved were excluded. When more than one eligible assessment was available for the same patient during a study period, only the first qualifying assessment was retained. Patient identifiers were cross-checked between study periods to prevent the same patient from contributing observations to both cohorts. After application of the eligibility criteria, 600 patients treated from January 2021 through December 2023 were included in the development-period cohort, and 360 patients treated from January 2024 through December 2025 were included in the temporal validation cohort. Each patient contributed one patient-level analytical record.

2.3. Data Collection and Outcome Definition

Clinical information was extracted according to a prespecified data dictionary, and identical variable definitions, coding procedures, measurement units, and extraction rules were applied to the development and temporal validation cohorts. Predictor values and index cardiac-biomarker measurements were obtained from the same clinical assessment episode. When multiple predictor measurements were available, the result obtained closest to and before, or on the same day as, the index BNP and hs-cTnI assessment was retained. Blood samples were obtained as part of routine clinical care, but the retrospective source systems did not consistently encode whether each specimen was collected before dialysis, during dialysis, or after dialysis, nor the immediately preceding ultrafiltration volume. Sampling time relative to dialysis therefore could not be standardized retrospectively and may have influenced BNP and hs-cTnI concentrations. Body mass index (BMI) was calculated as weight in kilograms divided by height in meters squared. The prognostic nutritional index (PNI) was calculated as serum albumin in g/L plus five times the absolute lymphocyte count in 10^9/L, and the neutrophil-to-lymphocyte ratio (NLR) was calculated as the absolute neutrophil count divided by the absolute lymphocyte count. D-dimer concentrations were natural-log transformed because of their right-skewed distribution. The primary outcome was a study-defined contemporaneous composite elevation of cardiac biomarkers. A patient was classified as having the outcome when BNP was at least 1045.13 pg/mL or hs-cTnI was at least 47.62 ng/L. These thresholds corresponded to the 75th percentiles of the observed BNP and hs-cTnI distributions in the complete development-period cohort and were calculated using the type 8 quantile estimator before the cohort was randomly partitioned. The outcome definition was subsequently applied unchanged to the training set, same-period holdout set, and temporal validation cohort. BNP and hs-cTnI were used exclusively to define the outcome and were excluded from feature selection and model development. Because this was a study-specific contemporaneous biomarker endpoint, these thresholds should not be interpreted independently as diagnostic cutoffs for heart failure or acute myocardial infarction.

2.4. Data Preprocessing, Partitioning, and Feature Selection

The completeness and consistency of the data were assessed before model development. Because no missing values were present in the variables included in the analytical datasets, neither single nor multiple imputation was performed. Data-quality checks were conducted to identify duplicate observations, impossible values, inconsistent measurement units, and data-entry errors. No patient was excluded solely because of an extreme but clinically plausible value, and no Winsorization or model-driven outlier removal was performed. After the outcome had been defined, the development-period cohort was divided once into a training set and a same-period holdout set at a ratio of 7:3 using outcome-stratified random sampling and random seed 123. Because the biomarker percentile thresholds had been estimated from the complete development-period cohort before this split, the holdout set was not fully independent of outcome-definition development and is described as same-period holdout evaluation rather than strict internal validation. Least absolute shrinkage and selection operator (LASSO) logistic regression was conducted exclusively in the training set to select predictors for model development. The penalty parameter lambda was determined through 10-fold cross-validation using the one-standard-error rule (lambda.1se), and predictors with nonzero coefficients at the selected lambda value were retained. The selected predictor set was subsequently fixed and used in all base models and the stacking ensemble, with no further feature selection performed in the holdout or temporal validation cohorts. Spearman correlation coefficients were calculated among the retained continuous predictors, and predictor pairs with an absolute correlation coefficient greater than 0.50 were flagged for review rather than automatically excluded.

2.5. Model Development and Stacking Ensemble

Six base models were developed using the predictors retained by LASSO: logistic regression, classification and regression tree (CART), elastic net, random forest, radial support vector machine, and gradient boosting machine. Outcome-stratified 10-fold cross-validation was conducted within the training set, and identical fold assignments were used for all base models. Hyperparameters were selected by maximizing the mean cross-validated area under the receiver operating characteristic curve (ROC-AUC), and out-of-fold predicted probabilities were retained for every patient in the training set. Centering and scaling were incorporated into the modeling pipelines for logistic regression, elastic net, and the radial support vector machine, whereas tree-based models used predictors on their retained analytical scales. Logistic regression was fitted without a tuning grid. For elastic net, alpha values of 0, 0.5, and 1 and 10 logarithmically spaced lambda values ranging from 0.001 to 1 were evaluated. CART evaluated 10 complexity-parameter values ranging from 0.001 to 0.05. Random forest used 600 trees and evaluated mtry values of 3, 5, 8, and 11 and minimum node sizes of 5 and 12 under Gini splitting. For the radial support vector machine, automatic tuning was performed through the caret training interface using a tuning length of six. Gradient boosting evaluated 75, 125, and 200 trees; interaction depths of 1, 2, and 3; shrinkage values of 0.03 and 0.06; and a minimum of 10 observations per terminal node.

The three base models with the highest mean cross-validated ROC-AUC values in the training set were selected for inclusion in the stacking ensemble. Their out-of-fold predicted probabilities were used as inputs to an L1-penalized logistic meta-model. To obtain training-set stacking predictions, each outer fold was predicted using a meta-model fitted on the remaining folds, with the penalty parameter selected through an additional five-fold cross-validation procedure. A final meta-model was subsequently fitted using the complete out-of-fold prediction matrix, with its penalty parameter selected through 10-fold cross-validation. The final meta-model was used to combine the corresponding base-model predictions in the same-period holdout and temporal validation cohorts. No feature selection, preprocessing parameter, model coefficient, hyperparameter, stacking weight, or classification threshold was re-estimated using the temporal validation cohort.

2.6. Model Performance Evaluation

ROC-AUC was the primary measure of discrimination and was reported separately for the training out-of-fold, same-period holdout, and temporal validation data sets. Precision-recall area under the curve was calculated by trapezoidal integration of precision over recall. Additional classification measures included accuracy, balanced accuracy, sensitivity, specificity, positive predictive value, negative predictive value, F1 score, Matthews correlation coefficient, Jaccard index, and detection prevalence. For each model, the classification threshold was selected by maximizing the Youden index using its training-set out-of-fold predictions and was then applied unchanged to the same-period holdout and temporal validation cohorts. Calibration was described by Brier score and by grouping predicted probabilities into 10 approximately equal-sized groups and comparing mean predicted probability with observed event proportion. Calibration-in-the-large/intercepts, calibration slopes, and confidence intervals were not estimated in the original analysis; therefore, the reported calibration results should be interpreted as descriptive rather than as a complete assessment of systematic miscalibration. Likewise, ROC-AUC, precision-recall AUC, Brier score, and threshold-based metrics are reported as point estimates without confidence intervals, so small differences between models should not be interpreted as evidence of superiority. Decision curve analysis was conducted over threshold probabilities from 0.01 to 0.80 in increments of 0.01. Net benefit was calculated as the true-positive proportion minus the false-positive proportion weighted by the odds of the threshold probability and was compared with the treat-all and treat-none strategies. Because no clinical intervention was assigned according to model predictions, decision curve analysis was interpreted as an assessment of potential decision-analytic value rather than direct evidence of improved clinical outcomes.

2.7. SHAP Analysis Methods

SHapley Additive exPlanations (SHAP) were used to interpret the stacking ensemble. SHAP values were estimated on the predicted-probability scale using a Monte Carlo approximation with 60 repetitions and adjustment enabled. The training predictor matrix was used as the background dataset, and patients in the same-period holdout set were used as the cases to be explained. Global feature importance was defined as the mean absolute SHAP value for each predictor. SHAP summary plots displayed the magnitude and direction of predictor contributions, whereas dependence plots illustrated the relationships between predictor values and their corresponding SHAP contributions.

2.8. Statistical Analysis

Continuous variables were summarized as medians and interquartile ranges (IQRs) and were compared between the development-period and temporal validation cohorts using the two-sided Wilcoxon rank-sum test. Categorical variables were summarized as frequencies and percentages and were compared using Pearson’s chi-square test or Fisher’s exact test, as appropriate. Baseline comparisons were descriptive and were not used for feature selection, model tuning, or modification of the fitted models. All statistical analyses were performed using R version 4.4.3. All tests were two-sided, and a P value < 0.05 was considered statistically significant.

3. Results

3.1. Participants

A total of 960 patients undergoing maintenance hemodialysis were included: 600 in the January 2021-December 2023 development-period cohort and 360 in the January 2024-December 2025 temporal validation cohort. The development-period cohort was divided into a training set (n = 420) and a same-period holdout set (n = 180). Overall, the median age was 64.0 years (IQR 54.0 - 72.0), and 590 of 960 patients (61.5%) were male. High cardiac biomarker burden was present in 340 of 960 patients (35.4%), including 213 of 600 (35.5%) in the development-period cohort and 127 of 360 (35.3%) in the temporal validation cohort (P = 0.944; Table 1). Compared with the development-period cohort, the temporal validation cohort had higher red cell distribution width (median 14.59%, IQR 13.66 - 15.52 vs 14.30%, IQR 13.50 - 15.27; P = 0.006), neutrophil count (median 5.37, IQR 3.60 - 7.83 vs 4.80, IQR 3.36 - 6.62 × 10^9/L; P = 0.008), beta-2 microglobulin concentration (median 21.25, IQR 17.21 - 26.18 vs 19.25, IQR 15.84 - 23.95 mg/L; P < 0.001), and neutrophil-to-lymphocyte ratio (median 6.46, IQR 4.27 - 11.27 vs 5.77, IQR 3.79 - 9.67; P = 0.033). BNP and high-sensitivity cardiac troponin I distributions did not differ between cohorts (P = 0.971 and P = 0.322, respectively; Table 1).

Table 1. Baseline characteristics of the development-period and temporal validation cohorts.

Values are median [25th percentile; 75th percentile] or n (%). Continuous variables were compared using the two-sided Wilcoxon rank-sum test; categorical variables were compared using Pearson’s chi-square test. P values compare the January 2021-December 2023 development-period cohort with the January 2024-December 2025 same-hospital temporal validation cohort; this later-period cohort is not an external validation cohort. BMI, body mass index; SBP, systolic blood pressure; DBP, diastolic blood pressure; FBG, fasting blood glucose; WBC, white blood cell; PLT, platelet; HB, hemoglobin; MCV, mean corpuscular volume; RDW, red cell distribution width; MPV, mean platelet volume; NEU, neutrophil; LYM, lymphocyte; MONO, monocyte. TP, total protein; ALB, albumin; GLB, globulin; ALT, alanine aminotransferase; AST, aspartate aminotransferase; GGT, gamma-glutamyl transferase; ALP, alkaline phosphatase; Scr, serum creatinine; UA, uric acid; PTH, parathyroid hormone; B2M, beta-2 microglobulin; PNI, prognostic nutritional index; NLR, neutrophil-to-lymphocyte ratio. TC, total cholesterol; TG, triglyceride; LDL-C, low-density lipoprotein cholesterol; HDL-C, high-density lipoprotein cholesterol; hs-cTnI, high-sensitivity cardiac troponin I; BNP, B-type natriuretic peptide; PT, prothrombin time; INR, international normalized ratio; FIB, fibrinogen; APTT, activated partial thromboplastin time; TT, thrombin time; FDP, fibrin degradation products.

3.2. Feature Selection and Correlation Assessment

Ten-fold cross-validated least absolute shrinkage and selection operator logistic regression performed in the training set retained 13 predictors with nonzero coefficients at the lambda.1se value: age, male sex, body mass index, systolic blood pressure, diabetes mellitus, coronary heart disease, hemoglobin, red cell distribution width, prognostic nutritional index, neutrophil-to-lymphocyte ratio, log-transformed D-dimer, beta-2 microglobulin, and phosphate (Figure S1). These 13 predictors constituted the fixed predictor set for all subsequent model development and evaluation analyses. Spearman correlation analysis was subsequently performed among the 10 continuous predictors. No predictor pair exceeded the prespecified absolute correlation threshold of 0.50. The strongest inverse correlation was observed between the prognostic nutritional index and log-transformed D-dimer (rho = −0.42), whereas the strongest positive correlation was observed between beta-2 microglobulin and phosphate (rho = 0.34). Therefore, no predictor was removed after the correlation assessment, and all 13 LASSO-selected predictors were carried forward into the six base models and the stacking analysis (Figure S2).

3.3. Model Development

Using the same 13 predictors retained by LASSO, 6 base models were developed: logistic regression, classification and regression tree, elastic net, random forest, radial support vector machine, and gradient boosting machine. Training-set out-of-fold evaluation identified elastic net and logistic regression as the strongest individual models for discrimination, with ROC-AUCs of 0.852 and 0.851 and precision-recall AUCs of 0.777 and 0.756, respectively. The radial support vector machine achieved a ROC-AUC of 0.835 and a precision-recall AUC of 0.744, followed by gradient boosting (0.828 and 0.749), random forest (0.821 and 0.732), and CART (0.744 and 0.646; Figures 1(a)-(e)). The unweighted overall performance scores ranged from 0.68 for CART to 0.77 for logistic regression (Figure 1(c)). Decision-curve analysis showed positive net benefit relative to the treat-none strategy across a broad range of low-to-moderate threshold probabilities, although the magnitude and persistence of benefit varied among models (Figure 1(f)). Grouped calibration plots provided a descriptive visual assessment of agreement between predicted and observed risks (Figure 1(g)). Cross-validation results were consistent with the out-of-fold estimates. Mean ROC-AUCs across the 10 folds were 0.86 (SD 0.06) for logistic regression, 0.86 (SD 0.07) for elastic net, 0.85 (SD 0.07) for the support vector machine, 0.84 (SD 0.06) for gradient boosting, 0.83 (SD 0.07) for random forest, and 0.75 (SD 0.05) for CART (Figure 1(h)). The corresponding precision-recall AUC results across the 10 cross-validation folds are shown in Figure 1(i). Logistic regression, elastic net, and the support vector machine consequently provided the 3 out-of-fold probability inputs for the stacking ensemble.

Figure 1. Development-period evaluation of six machine learning models. (a) Performance profiles across evaluation metrics; (b) Heat map of metric values; (c) Unweighted overall performance scores; (d) Training out-of-fold ROC curves; (e) Training out-of-fold precision-recall curves; (f) Decision-curve analysis; (g) Grouped calibration plots; (h) ROC-AUC across 10 cross-validation folds; (I) PR-AUC across 10 cross-validation folds.

3.4. Same-Period Holdout and Temporal Validation of Individual Machine Learning Models

In the same-period holdout set, logistic regression and elastic net had the highest discrimination point estimates among the individual models, with ROC-AUCs of 0.87 and 0.86, respectively. The support vector machine and gradient boosting each achieved a ROC-AUC of 0.84, compared with 0.80 for random forest and 0.77 for CART. In the temporal validation cohort, ROC-AUC point estimates were 0.83 for logistic regression, elastic net, and gradient boosting; 0.82 for the support vector machine; 0.81 for random forest; and 0.74 for CART (Figure 2(a)). The corresponding precision-recall curves for the same-period holdout and temporal validation cohorts are shown in Figure 2(b). Decision curves indicated that most models provided positive net benefit over treat-none across low-to-moderate thresholds in both evaluation cohorts, with net benefit generally diminishing at higher thresholds (Figure 2(c)). Descriptive calibration performance varied by model and validation period. Logistic regression had the lowest reported Brier score in same-period holdout evaluation (0.149), whereas gradient boosting had the lowest Brier score in temporal validation (0.162). The corresponding holdout and temporal validation Brier scores were 0.185 and 0.221 for CART, 0.189 and 0.190 for elastic net, 0.177 and 0.170 for random forest, and 0.161 and 0.173 for the support vector machine (Figure 2(d)). Because uncertainty intervals were not estimated, the small between-model differences in these point estimates should not be interpreted as evidence of superiority.

Figure 2. Performance of six machine learning models in the same-period holdout and temporal validation cohorts. (a) ROC curves; (b) Precision-recall curves; (c) Decision-curve analysis; (d) Grouped calibration plots with Brier scores.

3.5. Stacking Ensemble Development and Performance

The three base models with the highest mean cross-validated ROC-AUCs were logistic regression, elastic net, and the radial support vector machine; these models were selected as components of the stacking ensemble. Each component model was developed using the same 13 LASSO-selected predictors. Their out-of-fold predicted probabilities, rather than the original predictor values, were used as inputs to an L1-penalized logistic meta-model. The ensemble achieved ROC-AUC point estimates of 0.848 in training out-of-fold evaluation, 0.863 in same-period holdout evaluation, and 0.835 in temporal validation. The corresponding precision-recall AUC point estimates were 0.77, 0.76, and 0.70, and the Brier scores were 0.15, 0.15, and 0.16, respectively (Figure 3(a)). Thus, the ensemble remained comparable by point estimate with the strongest individual models; the small observed differences cannot establish algorithmic superiority without uncertainty intervals.

Figure 3. Performance and interpretation of the stacking ensemble. (a) Performance profiles and heat maps for the base models and stacking ensemble across training out-of-fold evaluation, same-period holdout evaluation, and temporal validation; (b) ROC curves for the stacking ensemble across the three data sets; (c) Mean absolute SHAP values; (d) SHAP summary and categorical-feature dependence plots; (e) Continuous-feature SHAP dependence plots.

At the prespecified classification threshold, the stacking model had sensitivity and specificity point estimates of 0.68 and 0.85 in training out-of-fold evaluation, 0.60 and 0.86 in same-period holdout evaluation, and 0.75 and 0.76 in temporal validation. Accuracy point estimates were 0.79, 0.77, and 0.75 across these data sets, respectively (Figure 3(a)). The ROC curves for the stacking ensemble across the three data sets are shown in Figure 3(b).

3.6. Global and Individual SHAP Explanations

In the same-period holdout set, global SHAP analysis quantified the contributions of the 13 LASSO-selected predictors to the stacking ensemble. In descending order of mean absolute SHAP value, the predictors were prognostic nutritional index, log-transformed D-dimer, age, red cell distribution width, beta-2 microglobulin, coronary heart disease, male sex, body mass index, neutrophil-to-lymphocyte ratio, hemoglobin, diabetes mellitus, systolic blood pressure, and phosphate (Figure 3(c) and Figure 3(d)). Lower prognostic nutritional index values and higher values of log-transformed D-dimer, age, red cell distribution width, beta-2 microglobulin, and body mass index were associated with larger positive contributions to predicted risk across their observed ranges (Figure 3(e)). Coronary heart disease, male sex, and diabetes mellitus also contributed positively to model predictions, whereas systolic blood pressure and phosphate had comparatively small mean absolute SHAP values.

The single-patient explanation illustrated how these population-level patterns combined within an individual prediction. For the example patient, the predicted probability was 0.450 compared with an adjusted background expectation of 0.342. Red cell distribution width of 17.90%, log-transformed D-dimer of 1.92, male sex, body mass index of 24.01 kg/m², and a neutrophil-to-lymphocyte ratio of 8.19 increased the prediction, whereas a prognostic nutritional index of 44.65, age of 43 years, absence of coronary heart disease, beta-2 microglobulin of 16.86 mg/L, and hemoglobin of 83.00 g/L decreased it (Figure S3).

4. Discussion

In this single-center retrospective study of 960 adults receiving maintenance hemodialysis, we developed and temporally evaluated an interpretable model for a study-defined contemporaneous composite of high cardiac biomarker burden. Training-set LASSO retained 13 routinely available demographic, clinical, and laboratory predictors spanning nutritional, inflammatory, hematologic, coagulation-related, and uremic domains. The stacking ensemble had ROC-AUC point estimates of 0.848 in training out-of-fold evaluation and 0.835 in the later-period cohort; its temporal-validation precision-recall AUC was 0.70 and its Brier score was 0.16. SHAP analysis identified the prognostic nutritional index (PNI), log-transformed D-dimer, age, red cell distribution width (RDW), and beta-2 microglobulin as the leading contributors to model predictions. These findings support the feasibility of integrating routine data to characterize contemporaneous biomarker burden, but they do not establish future-risk prediction, a diagnostic test for heart failure or acute myocardial infarction, a causal model of cardiovascular injury, or an intervention that improves outcomes.

The ensemble should be interpreted as a relatively stable aggregation of complementary models rather than as evidence that stacking was uniformly superior. Logistic regression achieved a slightly higher ROC-AUC point estimate than stacking in same-period holdout evaluation (0.87 vs 0.863), and several individual models had temporal-validation ROC-AUCs between 0.81 and 0.83. The observed differences in point estimates are too small to support a claim of decisive algorithmic superiority, particularly because confidence intervals were not estimated. The strong performance of logistic regression and elastic net suggests that much of the available signal was captured by additive or regularized relations, whereas stacking may have helped preserve performance by combining their probability estimates with those from a radial support vector machine. Recent hemodialysis studies have likewise reported useful mortality prediction with random forests, deep learning, and other machine-learning approaches, but their outcomes, prediction horizons, feature sets, and validation designs differ substantially from those used here; their reported ROC-AUCs therefore should not be compared directly with ours [9]-[11].

Temporal validation was informative because the later-period cohort differed from the development-period cohort in RDW, neutrophil count, beta-2 microglobulin, and neutrophil-to-lymphocyte ratio (NLR), yet the outcome prevalence remained similar and discrimination declined only modestly by point estimate. This provides a limited test of robustness to calendar-time change within the same hospital, not external validation across independent settings. Methodological guidance regards temporal validation as weaker evidence of transportability than evaluation in an independent setting and recommends examining calibration and clinical utility in addition to discrimination [12]. The relevance of this distinction is illustrated by a recent eight-center external validation of a Bayesian-network mortality tool in hemodialysis, in which performance and sensitivity differed between incident and prevalent patients [13].

The outcome in this study was deliberately defined as biomarker burden rather than a clinical diagnosis or future event. This distinction is particularly important in dialysis. In 200 asymptomatic patients receiving hemodialysis, elevated hs-cTnI and hs-cTnT were common and the two troponins showed different changes after dialysis [14]. In a prospective cohort of 178 patients, ultrafiltration volume and dialysis treatment time influenced intradialytic hs-cTnI concentrations [15], whereas monthly measurements in another cohort demonstrated marked between-patient variation in chronic troponin baselines [16]. A randomized crossover study further showed that the trajectories of hs-cTnI, hs-cTnT, and N-terminal pro-B-type natriuretic peptide differed by dialysis modality and sampling time [17], and a smaller crossover trial found different pre-post patterns of hs-cTnI and natriuretic peptide concentrations with routine versus cooled dialysate [18]. These studies support cautious interpretation of a single absolute biomarker value. Because dialysis-relative sampling time and immediately preceding ultrafiltration volume were not consistently encoded in our retrospective source systems, residual measurement heterogeneity may have affected outcome classification and model performance. The present endpoint should not be used to rule in or rule out an acute cardiovascular condition.

BNP and hs-cTnI nevertheless provide complementary information about myocardial wall stress and injury. A review focused on dialysis emphasized that natriuretic peptide concentrations reflect a combination of cardiovascular disease, volume status, reduced renal clearance, and dialysis-related factors, while suggesting that serial measurements may be more informative than isolated values [19]. In 3182 participants with nondialysis chronic kidney disease, BNP and hs-cTnT contributed to a multibiomarker score associated with incident heart failure, although adding the biomarker score to clinical factors produced only modest gains in discrimination [20]. That study used a different population and a longitudinal clinical outcome, so it does not validate our percentile-based thresholds. It does, however, support the biological rationale for considering markers of stretch and injury together while maintaining clinical context.

PNI had the largest mean absolute SHAP value, and lower PNI contributed to higher predicted probability. Because PNI combines serum albumin and lymphocyte count, it summarizes nutritional and inflammatory information rather than a single pathway. In a cohort of 101,616 patients initiating hemodialysis, higher PNI was associated with progressively lower mortality and provided better one-year mortality prediction than albumin or lymphocyte count alone [21]. A 2024 meta-analysis of 29 dialysis cohorts also found that malnutrition was associated with higher mortality, with a similar direction in the hemodialysis subgroup [22]. These mortality data do not establish that low PNI causes cardiac biomarker elevation, but they make the prominence of PNI in our model clinically coherent and identify nutritional assessment as an appropriate component of subsequent clinical review. The contributions of D-dimer, NLR, and RDW should also be interpreted as associations with a broader systemic state rather than as disease-specific signals. In a direct study of 167 maintenance hemodialysis patients, D-dimer exceeded 500 micrograms/L in 75% overall and remained positive in 52% of patients without an additional acute illness or predisposing chronic disease, limiting its specificity for thromboembolic exclusion in this population [23]. A recent meta-analysis of 19 chronic hemodialysis studies associated elevated NLR with both all-cause and cardiovascular mortality [24], and a cohort study of 181 maintenance hemodialysis patients associated higher RDW with all-cause mortality while finding no significant relation with cardiovascular mortality [25]. These prior outcomes are different from the present composite, and the SHAP findings neither diagnose thrombosis nor justify anticoagulation. They instead suggest that coagulation-related, inflammatory, and erythropoietic abnormalities jointly inform the model’s estimate of biomarker burden. Beta-2 microglobulin ranked fifth in the global SHAP analysis. In 5332 participants from the Dialysis Outcomes and Practice Patterns Study, beta-2 microglobulin remained positively associated with mortality in the high-flux hemodialysis era and was also related to dialysis vintage, residual urine volume, and inflammation [26]. Its contribution here may therefore represent accumulated uremic and dialysis-related burden rather than a specific cardiac mechanism. Phosphate and systolic blood pressure had relatively small mean absolute SHAP values. Low importance in this model does not imply low clinical importance because feature importance is conditional on the selected predictors, their distributions, the fitted model, and the contemporaneous endpoint. Longer-term phosphorus exposure above 4.5 mg/dL was associated with cardiovascular mortality in 17,414 hemodialysis patients [27]; phosphate variability was associated with mortality in a cohort of 302,613 incident hemodialysis patients [28]; and a contemporary Australian and New Zealand registry analysis found a U-shaped relation between serum phosphate and mortality [29]. A 2026 systematic review and meta-analysis of dialysis cohorts similarly found that high phosphate was associated with higher cardiovascular mortality [30]. These longitudinal findings also show why a single phosphate measurement may have limited predictive prominence for the cross-sectional endpoint used here.

SHAP adds transparency by showing how each recorded feature shifted a prediction relative to a background expectation. This can help clinicians audit whether a prediction is driven by plausible recorded information and can expose unexpected dependence on a feature. However, explainable artificial intelligence frameworks distinguish population-associative explanations from mechanistic explanations. The observed SHAP direction and magnitude therefore should not be interpreted as an intervention target, a causal effect, or proof that changing a predictor will change BNP, hs-cTnI, or clinical outcomes.

The next step is prospective, multicenter evaluation with standardized biomarker sampling, an explicitly defined clinical action, and assessment of workflow and patient consequences. The size of an external validation cohort should be chosen to estimate calibration, discrimination, and net benefit with adequate precision rather than by an arbitrary event count [31]. If the model is embedded in clinical decision support, early-stage evaluation should also assess safety, human factors, and actual clinical performance before any large-scale implementation or impact trial [32].

This study has several strengths but also important limitations. First, feature selection and model fitting were confined to the training data, the 13-predictor set was fixed before evaluation, and the stacking model used out-of-fold predictions. Performance was evaluated in a later-period cohort without refitting using multiple point-estimate measures, together with global and individual SHAP analyses. Second, the retrospective single-center design and same-hospital temporal cohort limit generalizability and do not constitute independent external validation. Third, the percentile-based composite was study-specific and contemporaneous rather than diagnostic or prognostic. Its thresholds were estimated from the complete development-period cohort before the training/holdout split; this allowed holdout information to influence outcome-definition development, so the same-period holdout evaluation was not fully independent. Re-estimating thresholds from the training set would change some outcome labels and would require rebuilding and reevaluating the complete modeling pipeline, which was not performed in the present revision. The later-period temporal cohort remained uninvolved in threshold estimation. Fourth, dialysis-relative sampling time and ultrafiltration immediately before sampling were not consistently available, and longitudinal dialysis-related variables were unavailable; these factors may influence biomarker concentrations. Fifth, calibration intercepts and slopes and confidence intervals for performance measures were not estimated, so grouped calibration plots, Brier scores, and small between-model differences should be interpreted descriptively. The complete-case dataset also did not permit evaluation of model performance under missingness. Finally, retrospective decision-curve analysis and the absence of subgroup fairness or prospective impact evaluation preclude claims of clinical superiority or implementation readiness.

5. Conclusion

A stacking ensemble based on 13 LASSO-selected routine clinical variables showed comparable point-estimate discrimination in same-period holdout and temporal validation for a study-defined contemporaneous composite of high BNP or hs-cTnI burden in patients receiving maintenance hemodialysis. PNI, D-dimer, age, RDW, and beta-2 microglobulin contributed most strongly to predictions, linking the model to nutritional, inflammatory, coagulation-related, and uremic domains. The findings support further investigation of an interpretable classification approach, but not future-risk prediction or immediate diagnostic or therapeutic use. Multicenter external validation, standardized longitudinal sampling, clinically anchored thresholds, complete calibration assessment, uncertainty quantification, and prospective impact evaluation are required before the model can inform care.

Ethics Approval and Consent to Participate

This retrospective observational study was conducted in accordance with the Declaration of Helsinki and relevant institutional requirements. The study protocol was reviewed and approved by the Ethics Committee of The First Affiliated Hospital of Yangtze University (Approval No. YJ202570). Because the study used routinely collected, de-identified clinical data and involved no prospective intervention or direct participant contact, the requirement for written informed consent was waived by the ethics committee. All data were accessible only to authorized investigators and were handled in accordance with institutional requirements for patient privacy and confidentiality.

Availability of Data and Materials

The de-identified patient-level data supporting the findings of this study are not publicly available because they contain sensitive clinical information and are subject to ethical and institutional data-protection restrictions. Reasonable requests for access may be directed to the corresponding author and will be considered subject to approval by The First Affiliated Hospital of Yangtze University and its ethics committee, and completion of an appropriate data-use agreement.

Code Availability

The analysis code used to develop and evaluate the models is available from the corresponding author upon reasonable request.

Author Contributions

All authors contributed to the study conception and design. Writing-original draft preparation: Tong Wang; Writing-review and editing: Jingchun Yang; Conceptualization: Xiaoxiao Zhang; Methodology: Tong Wang; Formal analysis and investigation: Tong Wang; Resources: Jing Wang; Supervision: Lin Li. All authors commented on previous versions of the manuscript. All authors read and approved the final manuscript.

Supplement

Figure S1. LASSO selection of model predictors. (a) Predictor coefficient trajectories across log-transformed λ values. (b) Ten-fold cross-validated mean squared error used to select the penalty parameter.

Figure S2. Spearman correlation matrix of the continuous LASSO-selected predictors. Cell values represent Spearman correlation coefficients, and no predictor pair exceeded the prespecified absolute correlation threshold of 0.50.

Figure S3. Individual SHAP explanation for an example patient. Contributions of individual features that increased or decreased the predicted probability relative to the adjusted background expectation.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References