This retrospective study included patients who underwent mechanical heart valve replacement surgery at Beijing Anzhen Hospital, Capital Medical University, and subsequently achieved a stable warfarin dose. Comprehensive clinical data and genetic information was collected from patients who underwent surgery between January 2021 and January 2024.

The inclusion criteria were: (1) northern Chinese population; (2) patients aged $\ge$18 years on postoperative warfarin anticoagulation therapy, with INR monitoring; (3) achieved warfarin anticoagulation stability, as defined by warfarin usage for $\ge$3 months, with a consistent dose over the last three consecutive follow-ups (interval $\ge$7 days), and INR maintained between 1.5–2.5. A dose is considered to be a stable maintenance dose once anticoagulation stability is achieved.

The exclusion criteria were: (1) the patient had severe liver and kidney dysfunction; (2) abnormal preoperative INR values (0.8 or $>$1.2); (3) inability to determine genotype; (4) concomitant administration of nonsteroidal anti-inflammatory drugs (NSAIDs); (5) hematologic disorders (e.g., hemophilia, thrombocytopenia, and platelet function disorders); (6) undergoing other concurrent cardiac surgeries (except for the Cox-Maze procedure); (7) patient experienced thromboembolic events, bleeding, or other anticoagulation-related complications, or died during hospitalization or follow-up. Patients who met any of the above criteria were excluded from the study.

The study protocol received ethical approval from the Institutional Review Board of Beijing Anzhen Hospital (Ethics Approval Number: 2020067X), and was conducted in accordance with the Declaration of Helsinki. All participants voluntarily provided written informed consent prior to enrollment.

To ensure the reliability of the study results, the effect size was set to medium (Cohen’s d = 0.5), the significance level to 0.05 ($\alpha$ = 0.05), and the power to 0.8 (1 – $\beta$ = 0.08). Based on these parameters, the sample size was calculated using G*Power software (version 3.1; Heinrich-Heine-Universität Düsseldorf, Düsseldorf, North Rhine-Westphalia, Germany). The study required at least 332 valid samples to achieve sufficient power in regression analysis tasks. A total of 413 subjects were actually included in the study, ensuring an adequate sample size.

Patient clinical data was collected through face-to-face interviews, regular phone calls, and review of the medical records from our hospital information system. These data included demographic characteristics, comorbidities, concomitant medications, preoperative echocardiographic parameters, cardiac function classification (New York Heart Association (NYHA) classification), surgical approaches, stable maintenance dose of warfarin therapy, and steady-state INR. The demographic characteristics included gender, age, height, weight, body surface area (BSA), body mass index (BMI), smoking history, and alcohol consumption history. Comorbidities were defined as chronic conditions or diseases diagnosed before the patient underwent heart valve replacement surgery. These included hypertension, diabetes, cerebrovascular disease, coronary artery disease and atrial fibrillation (AF). The presence of each comorbidity was confirmed by assessing the patient’s medical history and relevant diagnostic test results. Concomitant medications were defined as any medication being taken by the patient alongside the prescribed warfarin regimen, before or after surgery. This included digoxin, amiodarone, statin, and angiotensin converting enzyme inhibitors (ACEI) that could potentially interact with warfarin. A detailed list of these medications was obtained from the patient’s medical records. Preoperative echocardiographic parameters included ejection fraction (EF), left ventricular end-diastolic diameter (LVEDD), left ventricular end-systolic diameter (LVESD), and left atrial diameter (LAD).

Blood samples were drawn from each patient on the morning of the day before surgery and under fasting conditions. The levels of B-type natriuretic peptide (BNP), alanine transaminase (ALT), aspartate aminotransferase (AST), creatinine, triglycerides (TG), total cholesterol (TC), high density lipoprotein (HDL), low density lipoprotein (LDL), INR, fibrinogen, prothrombin time (PT) and other biochemical parameters were determined in the clinical laboratory of Beijing Anzhen Hospital, Capital Medical University, Beijing, China.

A total of 5 mL of venous blood was collected from eligible patients the day before surgery using an EDTA vacuum tube (BD Vacutainer, Franklin Lakes, NJ, USA) and stored at 4 °C before DNA extraction. Genomic DNA was extracted using the Blood DNA System DNA isolation kit (Catalog number: CW2320, CWBIO, Beijing, China) according to the manufacturer’s instructions and stored at 4 °C for subsequent use. DNA quality and purity were assessed by agarose gel electrophoresis and by measuring optical absorbance at A260/A280.

Genotyping of 6 SNPs was performed using the Illumina SNP GoldenGate Assay (Catalog number: WG-2500A-1001, Illumina, San Diego, CA, USA) according to the manufacturer’s instructions. These SNPs were CYP2C9 (rs1957910), VKORC1 (rs9923231), CYP4F2 (rs2108622), APOE (rs7412), CYP1A2 (rs2069514) and CYP3A4 (rs28371759). Briefly, 250 ng of genomic DNA was amplified at 37 °C for 20 h, followed by fragmentation and precipitation. The dried pellet was then resuspended and hybridized to the beadchips. These were subsequently incubated at 48 °C for 20 h, washed, and a single-base extension step performed. The beadchips were then stained, washed, coated, and dried. Finally, signal-intensity data were generated by an Illumina BeadArray Reader (Catalog number: BeadArray-1000, Illumina, San Diego, CA, USA). Twenty percent of the samples were selected at random for duplicate genotyping, with 99.8% concordance observed. Inconsistent data were excluded from the final analysis.

The collected data were entered into a Microsoft Excel spreadsheet and saved as a .csv file. During data preprocessing, missing values were first analyzed and the percentage of missing data determined for each variable to ensure completeness. If the missing data for a variable was 5%, the average of multiple estimates was used to fill in. If the missing data for a variable was $\ge$10%, the variable was removed or filled in using the mean value, depending on the distribution of data. If an individual patient record had $>$30% missing values, the data for that patient was deleted. Outliers were then detected using box plots and the 3$\sigma$ rule (i.e., values greater than $\pm$3 times the standard deviation of the mean were considered outliers). If the outliers were input errors, they were manually corrected or excluded. If they were extreme values, the data were excluded as appropriate.

All clinical and genetic data were included as candidate variables. Analysis of covariance (ANCOVA) was used to evaluate the effect of clinical and genetic variables on the target variable while controlling for potential confounders, selecting variables with statistical significance based on p-values (p 0.05) and partial $\eta$² effect sizes ($\eta$² $\ge$0.002). Initially, candidate variables were screened by their p-values (p 0.05), ensuring that only those significantly related to the target variable were retained. Among the remaining variables, partial $\eta$² was then calculated to assess the effect size. Variables with a very small contribution (i.e., $\eta$² 0.002) were excluded, as they were deemed to have minimal impact on model performance. Multiple iterations of this process were performed, with each iteration involving the recalculation of partial $\eta$² values and p-values. In each round, the variables contributing least to model performance were excluded, and the remaining filtered variables were used as input for the ML algorithm, with the stabilized warfarin dose as the output variable. By combining these statistical tools with effect size measures, more scientifically meaningful features were selected to improve model performance and interpretability.

The processed dataset was used for model development and validation, as shown in Fig. 1. The patients included in the study were randomized into the training set and validation set at an 8:2 ratio. Based on previous literature [21, 22, 23], we selected both traditional linear models, such as the general linear model (GLM), and 10 commonly used ML algorithms to build the models, allowing for comparison between traditional and modern approaches. Specifically, these models included support vector machine with a linear kernel (SVM Linear), support vector machine with radial basis function kernel (SVM Radial), recursive partitioning and regression trees (RPART), gradient boosting machine (GBM), random forest (RF), generalized linear model with elastic net regularization (Glmnet), extreme gradient boosting with a linear booster (XGB Linear), kernel k-nearest neighbors (KKNN), convolutional neural network (CNN), and extreme gradient boosting (XGB). Considering the small size of the training and validation sets, we used 10-fold cross-validation (10-fold CV) for model evaluation and grid search for parameter optimization to ensure the robustness and accuracy of the model. In order to reduce the risk of overfitting, 10-fold CV was performed on the training set. In this process, the training set was first randomly divided into 10 subsets, with the model using 9 subsets as training data and the remaining subset as validation data for performance evaluation in each round. To further enhance the predictive performance of the model, grid search was used to optimize the model hyperparameters. The grid search method traverses all possible parameter combinations in the preset hyperparameter space, uses cross-validation to evaluate the performance of each set of hyperparameters, and finally selects the best hyperparameter combination. Specifically, hyperparameter optimization was performed for models such as SVM, RF, and XGB by adjusting parameters such as the C value, gamma value, number of trees, maximum depth, and learning rate to improve the model’s generalization and predictive accuracy. The best-performing model was selected after final evaluation, and its performance was validated using the validation set to ensure it had good generalization ability and clinical application potential. The parameter settings for the 10 ML algorithms are detailed in Supplementary Table 1.

Fig. 1.

Flow chart used for modeling. ANCOVA, analysis of covariance; GLM, general linear model; MAE, mean absolute error; ML, machine learning; CV, cross-validation.

To assess the predictive accuracy and performance of the model, the primary evaluation metrics used were mean absolute error (MAE) and R-squared (R²). R² is the coefficient of determination (also known as goodness of fit) and represents the extent to which the regression model explains the variation in the dependent variable, or how well the model fits the observed data. The value of R² ranges from 0 to 1. Generally, the closer R² is to 1, the better the model fit, indicating a higher degree of explanation of the dependent variable by the independent variables. A higher R² value implies that a larger percentage of the total variability is explained by the model. When the observed values are closer to the regression line, the fit is more accurate [24].

The MAE metric is used to evaluate the accuracy of a predictive model by calculating the average of the absolute differences between the actual values and the predicted values. Lower MAE and higher R² indicate higher predictive accuracy, and thus better model performance.

Additionally, the ideal prediction percentage was used to evaluate the clinical utility of the model, defined as the percentage of predicted doses within $\pm$20% of the actual dose. An increase in the ideal prediction percentage indicates greater clinical utility of the model.

A dose subgroup analysis was conducted to account for heterogeneity in clinical practice. The conventional warfarin dose in China is 2.5 mg/day. Based on clinical recommendations and existing studies [16, 17, 18], dose thresholds were defined as the 25th and 75th percentiles of the conventional dose, i.e., high-dose ($\ge$3.125 mg/day), medium-dose (1.875–3.125 mg/day), and low-dose ($\le$1.875 mg/day). This categorization enables a more precise evaluation of the model’s efficacy across different dose levels, thereby optimizing clinical decision support.

Continuous variables are presented as means $\pm$ standard deviations (SD) or median (interquartile range (IQR)), while categorical variables are presented as frequencies (n) and percentages (%). The normality of the continuous variables was assessed using the Shapiro-Wilk test. For comparisons of variance between groups, the F-test was applied to assess the homogeneity of variances. To compare continuous variables between two independent groups, the Student’s t-test was used when the data followed a normal distribution, while the Mann-Whitney U test was applied for non-normally distributed variables. For categorical variables, intergroup differences were assessed using the Chi-square test, with Fisher’s exact test applied when expected frequencies in any cell of the contingency table was less than 5.

Multiple ML techniques were utilized to model and analyze the data, with a variety of R packages employed to perform these analyses. These included the stats package (version 4.2.2; R Foundation for Statistical Computing, Vienna, Austria; https://www.R-project.org/) for basic statistical functions, e1071 (version 1.7-14; https://CRAN.R-project.org/package=e1071) for SVM implementations, RPART (version 4.1.23; https://CRAN.R-project.org/package=rpart) for recursive partitioning and decision tree algorithms, GBM (version 2.1.9; https://CRAN.R-project.org/package=gbm) for generalized boosting machines, RF (version 4.7-1.1; https://CRAN.R-project.org/package=randomForest) for ensemble learning, Glmnet (version 4.1-8; https://CRAN.R-project.org/package=glmnet) for penalized regression models, XGBoost (version 1.7.7.1; https://CRAN.R-project.org/package=xgboost) for gradient boosting machines, Class (version 7.3-22; https://CRAN.R-project.org/package=class) for classification algorithms, and Keras (version 2.15.0; https://CRAN.R-project.org/package=keras) for deep learning models. All statistical and machine learning analyses were performed using R (version 4.2.2; R Foundation for Statistical Computing, Vienna, Austria; https://www.R-project.org/). For all statistical tests, a significance threshold of p-value 0.05 was considered to indicate statistical significance, and all tests were conducted with a two-sided approach.

Following execution of the inclusion and exclusion criteria, a total of 413 patients were enrolled in the study, including 227 males (54.96%). The median age of participants was 53 (46–59) years, with an average height of 165.86 $\pm$ 8.37 cm, median weight of 65 (59–75) kg, average BSA of 1.75 $\pm$ 0.18 m², and median BMI of 24.02 (22.04–27.03) kg/m². The median stable daily warfarin dose was 3 (2.25–3.75) mg, and the median steady-state INR was 1.97 (1.78–2.20). A total of 332 patients were randomly selected as the training set (80%), while the remaining 81 patients were used as the validation set (20%). Baseline characteristics of the study participants are presented in Table 1. No significant differences in demographic baseline data or clinical and genetic information were found between the training and validation cohorts (all p $>$ 0.05).

ANCOVA was used to select variables from the potential independent factors. Based on these criteria, 13 primary input variables were chosen from the initial factors for subsequent construction of the model (Table 2), namely sex, age, BSA, height, weight, smoke, digoxin, AF, TC, LDL, Cox-Maze, VKORC1, and CYP2C9. Among the variables included in the model, AF had the largest partial $\eta$² value (0.016), sex and TC had the smallest partial $\eta$² values (0.002), while the remaining variables had partial $\eta$² values ranging from 0.003 to 0.012.

To assess the predictive accuracy and performance of the models, we first evaluated the R² of 10 ML models and of the GLM model. In the training set, the SVM Radial model and the RF model have the highest R² of 0.98 and 0.95, respectively, while the SVM Linear model and the GLM model have the lowest R² of 0.31 and 0.33, respectively. The remaining models have R² ranging between 0.48 and 0.72 (Fig. 2A). In the validation set, the SVM Radial model had the highest R² of 0.98, demonstrating high stability and generalization ability, and indicating it has significant advantages in terms of handling complex data patterns and lower susceptibility to overfitting (Fig. 2B). However, RF had the lowest R² (0.22) in the validation set (Fig. 2B), performing noticeably worse than in the training set. This result suggests potential overfitting, where the model learns too many details and noise from the training data that do not apply to new, unseen data. The R² for the remaining models with the validation set ranged from 0.34 to 0.71 (Fig. 2B).

Fig. 2.

Comparison of R² for the training and validation sets for different algorithms. (A) Box plot of R² values for the training set. (B) Box plot of R² values for the validation set. The horizontal axis represents the different models and the vertical axis represents the R² values. The box plot shows the median, interquartile range and 95% confidence interval (CI) of the R² values. GLM, general linear model; SVM Linear, support vector machine with linear kernel; SVM Radial, support vector machine with radial basis function kernel; RPART, recursive partitioning and regression trees; GBM, gradient boosting machine; RF, random forest; Glmnet, generalized linear model with elastic net regularization; XGB Linear, extreme gradient boosting with linear booster; KKNN, kernel k-nearest neighbors; CNN, convolutional neural network; XGB, extreme gradient boosting.

Next, we calculated the MAE of all models. The MAE of SVM Radial was the lowest at 0.14 mg/day (95% CI: 0.11–0.17), indicating it was a stable model. In contrast, the MAE in the validation set was highest in RF at 0.93 mg/day (95% CI: 0.88–0.98) (Table 3). This indicates RF was an unstable model, probably due to the small sample size and overfitting.

Calculation of the ideal prediction percentage revealed the SVM Radial and RF models performed best. Moreover, the results obtained with the validation set (Fig. 3B) for these two models were as good as with the training set (Fig. 3A). The lowest ideal prediction percentage observed in the training set was GLM, while in the validation set it was SVM Linear.

Fig. 3.

Comparison of ideal prediction percentage for different algorithms. (A) Ideal prediction percentage in the training set. (B) Ideal prediction percentage in the validation set. The horizontal axis shows the different models, while the vertical axis shows the proportion of samples where the predicted value is within $\pm$20% of the true value. GLM, general linear model; SVM Linear, support vector machine with linear kernel; SVM Radial, support vector machine with radial basis function kernel; RPART, recursive partitioning and regression trees; GBM, gradient boosting machine; RF, random forest; Glmnet, generalized linear model with elastic net regularization; XGB Linear, extreme gradient boosting with linear booster; KKNN, kernel k-nearest neighbors; CNN, convolutional neural network; XGB, extreme gradient boosting.

In summary, the SVM Radial model demonstrated the highest R², the best ideal prediction percentage, and the smallest MAE. We therefore selected SVM Radial as the optimal model for predicting stable warfarin dose in patients following surgery for mechanical heart valve replacement.

To further evaluate the overall predictive performance of the SVM Radial model, we plotted the predicted dose versus actual dose fit assessment in both the training and validation datasets (Fig. 4). In the training and validation sets, R² was 0.98, indicating a strong correlation between the predicted values and the actual observations, and demonstrating the excellent predictive capabilities of the SVM Radial model.

Fig. 4.

Plot of predicted versus actual dose fit assessment for the training and validation sets. (A) Plot of predicted versus actual dose fit assessment for the training set. (B) Plot of predicted versus actual dose fit assessment for the validation set. The horizontal axis shows the predicted dose, and the vertical axis shows the actual dose. The solid red line indicates the ideal fit (y = x), and the black dots represent the median of the actual values within each prediction interval, while the error bars reflect the variability of the data within that interval.

In the dose subgroup analysis, the SVM Radial model demonstrated high accuracy for predicting warfarin doses. For the total sample, this model successfully predicted the ideal dose in 98.06% of cases. A detailed analysis of different dose levels revealed a prediction success rate of 95.00% for the 40 subjects requiring a low daily dose ($\le$1.875 mg/day), 98.78% for the 245 subjects in the medium-dose group (1.875–3.125 mg/day), and 97.66% for 128 subjects in the high-dose group ($\ge$3.125 mg/day) (Table 4).

The SVM Radial model accurately predicted the dose for 93.83% of the 81 subjects in the validation cohort. In the different dose groups, the prediction accuracy was 85.71% for 7 subjects in the low-dose group ($\le$1.875 mg/day), 95.92% for 49 subjects in the medium-dose group (1.875–3.125 mg/day), and 92.00% for 25 subjects in the high-dose group ($\ge$3.125 mg/day). These results demonstrate the superior performance of the SVM Radial model across various dose levels and highlight its potential value in clinical applications. In particular, this model has excellent predictive ability within the medium-dose range, providing guidance for precise dose adjustment in clinical practice.