This single-center, retrospective study consecutively enrolled hospitalized patients with type 2 diabetes mellitus (T2DM) who underwent dual-energy X-ray absorptiometry (DXA) at Panyu Hospital of Chinese Medicine between January 2022 and December 2023. Patients who underwent DXA examination during hospitalization were eligible for inclusion. The inclusion criteria were as follows: (1) age ≥ 45 years; (2) completed standardized DXA assessment; and (3) availability of complete electronic medical records, including demographic information and routine laboratory test results. Exclusion criteria included: (1) secondary osteoporosis; (2) hematological disorders; (3) history of malignancy; (4) severe hepatic or renal insufficiency; (5) acute infectious disease; and (6) incomplete clinical data.Bone mineral density (BMD) was measured using a standardized DXA protocol, assessing the lumbar spine (L1–L4), left femoral neck, and total hip. According to the diagnostic criteria defined by the World Health Organization, participants were classified into the osteoporosis group (T-score ≤ –2.5 SD) and the non-osteoporosis group (T-score > –2.5 SD) ( Figure 1 ).

The study protocol adhered to the Declaration of Helsinki and was approved by the Institutional Review Board of Panyu Hospital of Guangzhou University of Chinese Medicine. Given the retrospective nature of this study, informed consent was waived.

Baseline demographic characteristics (age and sex) and standardized laboratory parameters were collected. Venous blood samples were collected from all fasting participants within 24 hours of hospital admission. Complete blood counts, including hemoglobin (HGB), neutrophil count (NEUT), red blood cell count (RBC), platelet count (PLT), lymphocyte count (LYMPH), and monocyte count (MONO), were analyzed using the Mindray BC-6800Plus hematology analyzer. Biochemical analyses included total cholesterol (TC), triglycerides (TG), high-density lipoprotein cholesterol (HDL-C), low-density lipoprotein cholesterol (LDL-C), alkaline phosphatase (ALP), alanine aminotransferase (ALT), albumin (ALB), uric acid (UA), fasting blood glucose (FBG), serum calcium (Ca), serum phosphate (Pi), and creatinine (Cr). All laboratory evaluations and DXA measurements were performed during the same hospitalization period.

Before model development, several preprocessing steps were performed to ensure data quality and consistency. Duplicate records and entries with apparent input errors were excluded. For variables with missing values less than 20%, imputation was performed using the mean for continuous variables and the mode for categorical variables. Variables with more than 20% missing data were excluded from model construction.

To capture complex interactions related to metabolic and inflammatory status, we derived seven composite indices from the original laboratory values: Monocyte-to-HDL ratio (MHR), Neutrophil-to-HDL ratio (NHR), Platelet-to-HDL ratio (PHR), Lymphocyte-to-HDL ratio (LHR), Triglyceride-Glucose index (TyG), Cholesterol-Glucose index (CHG), and Non-HDL-to-Neutrophil ratio (NHHR). These indices have been previously associated with chronic inflammation, insulin resistance, and cardiovascular risk—factors that are also implicated in osteoporosis, particularly among patients with type 2 diabetes mellitus.Seven novel indices were derived from original laboratory measurements using the following formulas, as previously described in the literature (7–13):

In this study, the caret package in R was used to randomly partition the dataset into a training set (70%) and a testing set (30%). The training set was utilized for model development, while the testing set was reserved for performance evaluation.

To identify potential predictors from baseline variables in the training set, three independent feature selection methods were employed: univariate analysis, least absolute shrinkage and selection operator (LASSO) regression, and the Boruta algorithm (14). Univariate analysis, a conventional statistical method, selected variables with a P-value < 0.05. LASSO regression identified predictors with non-zero coefficients, effectively addressing multicollinearity and reducing the risk of model overfitting due to excessive inter-variable correlation (15). In this study, LASSO regression with 10-fold cross-validation was performed to screen variables from high-dimensional data.The Boruta algorithm, a wrapper method based on feature importance, identifies relevant variables by comparing the Z-scores of actual features to those of permuted “shadow” features. In each iteration, Z-scores for the real features were computed using a random forest (RF) classifier, while those for the shadow features were obtained by randomly shuffling the original variables (16). Features consistently performing worse than the shadow features were iteratively removed.To ensure that there was no multicollinearity among the candidate variables, we calculated the variance inflation factor (VIF) scores. Variables with high collinearity (VIF > 10) were excluded to maintain model stability.

The intersection of variables selected by all three methods was used to develop prediction models. Ten supervised machine learning algorithms were implemented for model construction: logistic regression (LR), support vector machine (SVM), gradient boosting machine (GBM), neural network (NN), extreme gradient boosting (XGBoost), k-nearest neighbors (KNN), AdaBoost, and LightGBM.

The optimal model was determined through a comprehensive evaluation of the discriminative ability, calibration performance, and clinical applicability of all ten candidate models. The receiver operating characteristic (ROC) curve was constructed to visually represent the model’s discriminative power, with the area under the curve (AUC) serving as the primary quantitative metric. Additional performance indicators, including accuracy, sensitivity, specificity, precision, and F1 score, were also calculated to supplement the evaluation. Calibration and clinical utility were further assessed using decision curve analysis (DCA) and calibration plots, respectively.

To enhance model interpretability, Shapley Additive Explanations (SHAP) values were analyzed to quantify the contribution and importance of each feature in determining the final classification outcome. A higher SHAP value indicates greater influence on the model’s output prediction. We present a feature importance analysis based on SHAP values to interpret the results of the optimal model (17).

All statistical analyses were conducted using R software (version 4.4.3). Continuous variables with a normal distribution were presented as mean ± standard deviation (SD), while non-normally distributed continuous variables were expressed as median with interquartile range (IQR). Categorical variables were summarized as frequencies and percentages. Group differences in continuous variables were assessed using either the independent samples t-test or the Kruskal–Wallis test, as appropriate. Categorical variables were compared using the chi-square test or Fisher’s exact test. A two-sided p-value < 0.05 was considered statistically significant.

Lasso regression and Boruta feature selection were performed using the “glmnet” and “Boruta” R packages, respectively. Predictive model development and training were conducted with the “caret” package, with default hyperparameter tuning implemented via grid search. SHAP value analysis for model interpretation was carried out using the “shapviz” package.

A total of 1,560 patients were included in this study and randomly assigned to the training set (n = 1,093) and the testing set (n = 467) in a 7:3 ratio. No statistically significant differences were observed between the two sets across all clinical variables (P > 0.05), as detailed in Table 1 .

In the training cohort, patients with and without osteoporosis showed statistically significant differences in several clinical and biochemical variables, as detailed in Table 2 . These included age, sex, hemoglobin (HGB), neutrophil count (NEUT), lymphocyte count (LYMPH), red blood cell count (RBC), triglycerides (TG), high-density lipoprotein cholesterol (HDL-C), alkaline phosphatase (ALP), alanine aminotransferase (ALT), albumin (ALB), uric acid (UA), serum calcium (Ca), as well as derived indices such as neutrophil-to-HDL ratio (NHR), platelet-to-HDL ratio (PHR), lymphocyte-to-HDL ratio (LHR), triglyceride-glucose index (TyG), and non-HDL-to-neutrophil ratio (NHHR) (all P < 0.05). Detailed comparisons are provided in Table 2 . Conversely, no statistically significant differences were found in platelet count (PLT), monocyte count (MONO), total cholesterol (TC), low-density lipoprotein cholesterol (LDL-C), aspartate aminotransferase (AST), creatinine (Cr), serum phosphate (Pi), fasting blood glucose (FBG), monocyte-to-HDL ratio (MHR), and cholesterol-glucose index (CHG) between the two groups (all P > 0.05).

Eighteen variables identified through univariate analysis were subjected to further feature selection using LASSO regression and the Boruta algorithm. In the LASSO regression, the optimal lambda value was determined to be 0.024, which yielded 10 key predictive features. These included two demographic variables (age and sex) and eight laboratory indicators (ALP, UA, HGB, NEUT, ALT, LHR, Ca, and TG). Notably, the results of the Boruta algorithm fully corroborated those of LASSO regression, confirming the same 10 features ( Figure 2 ). The candidate variables also showed VIF values below the accepted threshold of 10, indicating no significant Multicollinearity.These variables were ultimately selected as the core predictors for model construction.

Eight machine learning algorithms were employed to develop models for predicting osteoporosis in patients with T2DM: logistic regression (LR), support vector machine (SVM), gradient boosting machine (GBM), neural network (NN), extreme gradient boosting (XGBoost), k-nearest neighbors (KNN), AdaBoost, and LightGBM. Among these, the LR model demonstrated the most robust and consistent performance in the testing set ( Table 3 ).

As illustrated in Figure 3 , the LR model achieved the highest discriminative performance, with an AUC of 0.812 (95% CI: 0.750–0.832). It also outperformed other models in terms of accuracy (0.762), sensitivity (0.809), and F1 score (0.771), while maintaining favorable specificity (0.761) and precision (0.778). In comparison, KNN (AUC = 0.734, 95% CI: 0.689–0.779), AdaBoost (AUC = 0.717, 95% CI: 0.674–0.761), and LightGBM (AUC = 0.723, 95% CI: 0.676–0.770) demonstrated inferior classification performance. Although the SVM model achieved slightly higher specificity (0.808) and precision (0.801), its lower sensitivity (0.650), F1 score (0.717), and AUC (0.773, 95% CI: 0.731–0.815) indicated limited generalizability. Moreover, the GBM, NN, and XGBoost models exhibited noticeably lower AUCs, accuracy, sensitivity, and F1 scores compared to the LR model, further supporting the superior generalization ability of the latter.

To elucidate the decision-making process of the LR model, we applied the SHapley Additive exPlanations (SHAP) framework for both global and local interpretation. The bar chart in Figure 4A displays the mean absolute SHAP values for each predictor, reflecting their overall contribution to the model.The SHAP summary plot ( Figure 4B ) ranked the importance of predictive features, with age, sex, ALP, UA, HGB, NEUT, ALT, LHR, Ca, and TG identified as the most influential variables in descending order. Figures 4C, D provide further insights into how individual SHAP values influence predictions. Age emerged as the most dominant predictor, and interaction plots ( Figures 5A–F ) illustrated the complex interrelationships between age and other variables.

To enhance the model’s clinical utility and streamline its application in practice, we optimized the input requirements by limiting them to the 10 key features, many of which are routinely obtained in standard blood tests. Furthermore, we developed a user-friendly, web-based calculator (available at:https://t2dm.shinyapps.io/t2dm-osteoporosis/) to facilitate real-time risk assessment of osteoporosis in patients with T2DM. This platform enables clinicians to enter patient data and receive immediate, individualized risk estimations, thereby supporting early intervention and personalized care strategies.