Integration of Patellar Radiomics and Clinical Risk Factors via Machine Learning: A Proof-of-concept Study for Opportunistic Osteoporosis Screening

Osteoporosis is a systemic skeletal disease characterized by low bone mass and microarchitectural deterioration that increases fracture risk, and it remains substantially underdiagnosed and undertreated worldwide (1, 2). Estrogen withdrawal after menopause accelerates osteoclastic bone resorption and trabecular deterioration (3), while age-related bone loss further increases skeletal fragility (4).

Dual-energy X-ray absorptiometry (DXA) is the reference standard for measuring bone mineral density (BMD) and operationally defines osteoporosis as a T-score of -2.5 or lower, with osteopenia between -1.0 and -2.5 (5). Current screening recommendations target women aged 65 years or older and younger postmenopausal women at increased fracture risk, generally identified using clinical risk assessment (6, 7). Nevertheless, limited DXA access and underutilization leave many at-risk patients unscreened. This gap has motivated opportunistic screening, in which BMD-related information is extracted from imaging examinations obtained for other clinical indications without additional radiation exposure; opportunistic computed tomography (CT) screening may also be cost-effective (8).

Radiomics converts medical images into quantitative descriptors of intensity and texture that can be linked to biological characteristics through machine learning (9). Reproducibility across software, acquisition settings, and preprocessing remains a recognized challenge, and standardization initiatives have therefore become central to radiomics research (10). The approach conceptually parallels trabecular bone score, a gray-level textural index derived from lumbar DXA that reflects trabecular microarchitecture (11). Opportunistic BMD assessment has been demonstrated using deep learning on routine CT (12), proximal-femur CT radiomics (13), lumbar magnetic resonance apparent diffusion coefficient maps (14), clival CT radiomics (15), and volumetric Hounsfield unit analysis (16).

Plain radiographs are more widely available than CT or magnetic resonance imaging. A small exploratory study reported an association between patellar radiographic texture and lumbar BMD (17), and recent work indicates that knee (18) and foot or ankle radiographs (19) can support opportunistic osteoporosis screening. However, most radiograph-based studies have relied on deep-learning models with limited interpretability, and the contribution of simple clinical variables to a transparent patellar radiomics model has not been systematically examined.

We therefore evaluated whether patellar radiomic features from routine lateral knee radiographs, alone and in combination with age, body mass index (BMI), and estrogen-free years, could discriminate DXA-defined normal BMD, osteopenia, and osteoporosis using an interpretable Naive Bayes classifier. A secondary exploratory analysis evaluated the binary classification of osteoporosis versus non-osteoporosis.

Study Design and Patient Selection

This retrospective, single-center study was conducted in accordance with the Declaration of Helsinki. The institutional picture archiving and communication system and the hospital information system were searched to identify patients who had both a lateral knee radiograph and a DXA examination performed within 6 months of each other. A total of 146 potentially eligible patients were identified.

Thirty-seven patients were excluded: 14 were male; seven were younger than 18 years; eight had an implant, prosthesis, or a history of major fracture in the evaluated region; five had missing or incomplete radiographic or DXA data; two had inadequate image quality or anatomy unsuitable for standardized region-of-interest (ROI) placement; and one had advanced osteoarthritis. The final study cohort comprised 109 adult women (Figure 1).

Age, BMI, and estrogen-free years were retrieved from the hospital information system. Estrogen-free years were defined as the number of completed years since menopause. A value of 0 was assigned to women who were premenopausal or who had entered menopause less than 1 year previously; these two states could not be reliably distinguished retrospectively. Five participants had a value of 0 (normal BMD, n=2; osteopenia, n=3; osteoporosis, n=0). When menopausal age could not be determined from the medical records, the variable was treated as missing. DXA reports were reviewed, and the lowest T-score documented across skeletal sites assessed during the clinical examination was used as the reference standard. Participants were grouped as normal (T-score>-1.0), Participants were grouped as normal (T-score ≥ -1.0), osteopenia (-2.5 < T-score < -1.0), or osteoporosis (T-score ≤ -2.5) (5). The analytic dataset available for this retrospective study did not consistently retain data on the DXA scanner model, site-level quality-control records, or vertebra-level adjudication of degenerative/invalid lumbar measurements; the study, therefore, relied on the finalized clinical DXA report rather than on retrospective re-analysis of the raw DXA examination, and this limitation is acknowledged below.

Image Acquisition and Standardized ROI Placement

Lateral knee radiographs were obtained as part of routine clinical care, exported in DICOM format, and imported into 3D Slicer version 5.10.0 (www.slicer.org) (20). A custom Python script configured the Segment Editor paint tool with a fixed absolute diameter of 15 mm and set the segmentation reference geometry based on the source image. A musculoskeletal radiologist with 10 years of experience divided the patella into three equal craniocaudal thirds and placed the circular ROI using a single click at the junction between the middle and lower thirds (Figure 2). The predefined location was selected to predominantly sample trabecular bone while avoiding cortical margins, articular osteophytes, and the patellar tendon insertion. When more than one eligible lateral radiograph was available, the examination closest to the DXA date was selected. Because this was a retrospective routine-practice dataset, a single standardized exposure protocol could not be prospectively imposed; the per-examination radiographic system and exposure details were not consistently retained in the analytic dataset used for this study.

Radiomic Feature Extraction and Selection

Quantitative feature extraction was performed using PyRadiomics version 3.1.0 (9). Single-slice NRRD images and their corresponding masks were explicitly converted to two-dimensional SimpleITK image objects before extraction. Whole-image pixel intensities were normalized using the PyRadiomics Z-score normalization procedure (normalize=True; default normalizeScale=1), and gray levels within each ROI were discretized using a fixed bin count of 64. No spatial resampling was applied; the original DICOM in-plane pixel spacing was retained. Accordingly, texture features remained potentially sensitive to differences in detector sampling and acquisition, which was considered in the interpretation of the results.

Because the same geometrically fixed 15-mm circular ROI was used for every participant, shape features were disabled. Exactly 930 features were extracted per ROI. These comprised 18 first-order features and 75 texture features: 24 gray-level co-occurrence matrix features, 16 gray-level run-length matrix features, 16 gray-level size-zone matrix features, 14 gray-level dependence matrix (GLDM) features, and five neighboring gray-tone difference matrix features. The 93 features were calculated from the original image, four level-one, two-dimensional Coiflet-1 wavelet decompositions (LL, LH, HL, and HH), and five Laplacian-of-Gaussian-filtered images with sigma values of 0.5, 1.0, 1.5, 2.0, and 3.0 mm. Diagnostic metadata generated by PyRadiomics were excluded. The string “3D” retained in LoG-derived output names reflects the fixed PyRadiomics nomenclature and does not indicate a three-dimensional texture calculation in this workflow.

Dimensionality reduction was performed using the fast correlation-based filter (FCBF), which prioritizes features associated with the outcome while reducing redundancy (21). After re-examination of the archived analysis workflow, the entire modeling analysis was repeated using a strictly fold-contained feature selection procedure. FCBF was therefore fitted exclusively on each training fold and re-fitted independently in every fold and repetition. To reduce model complexity, given the limited number of osteoporosis cases, the principal revised analysis prespecified a maximum of 10 radiomic features. Sensitivity analyses used maximum feature counts of 5, 20, and 30. Continuous variables used for FCBF scoring were discretized by equal-frequency binning, consistent with the original Orange workflow. Clinical variables were not subjected to radiomic feature selection and were appended to the selected radiomic feature set when applicable.

Machine Learning and Statistical Analysis

A Naive Bayes classifier was retained to preserve the original modeling framework and its transparent probabilistic structure. The primary analysis was a three-class classification problem. Five feature sets were evaluated: (1) clinical-only (age, BMI, and estrogen-free years), (2) radiomics-only, (3) BMI+radiomics, (4) age+BMI+radiomics, and (5) full clinical+radiomics (age, BMI, estrogen-free years, and radiomics). A secondary exploratory binary analysis combined the normal and osteopenia groups into a non-osteoporosis class and compared clinical-only, radiomics-only, and combined clinical-radiomics models.

Internal validation was conducted using 20 repetitions of patient-level stratified 10-fold cross-validation with different fold assignments. For every training/held-out split, all data-dependent operations were confined to the training fold: FCBF feature selection was re-fitted, missing clinical values were imputed using training-fold means, standardization parameters were estimated from the training data, and the Naive Bayes model was fitted only after these steps. The fitted transformations and selected features were then applied unchanged to the corresponding held-out fold. Thus, no held-out observation contributed to imputation, scaling, feature selection, or model fitting (Figure 3).

Each participant generated one held-out predicted probability for each of the 20 repetitions. For the principal patient-level evaluation, the 20 held-out probabilities for each participant were averaged to obtain one aggregated out-of-fold prediction. Class labels were assigned to the class with the highest posterior probability; for binary classification, this corresponds to a posterior probability threshold of 0.5. In the three-class analysis, osteoporosis discrimination was summarized in a one-versus-rest framework. Overall multiclass performance was evaluated using accuracy, balanced accuracy, macro-F1 score, and linearly weighted Cohen’s kappa. A majority-class benchmark, assigning every participant to osteopenia, was calculated explicitly.

For the secondary binary analysis, area under the curve (AUC), accuracy, sensitivity, specificity, precision, F1 score, and Matthews correlation coefficient were calculated. Ninety-five percent CIs for the principal AUC, sensitivity, and specificity estimates were obtained using 5,000 patient-level bootstrap resamples of the aggregated out-of-fold predictions. The difference in AUC between the full clinical+radiomics and radiomics-only models was evaluated by paired patient-level bootstrap resampling, using the same resampled participants for both models. To provide a conservative sensitivity analysis addressing ambiguity in menopausal status, all five participants whose estrogen-free years were coded as 0 were excluded, and the binary analysis was repeated (n=104). The revised analysis was executed in Python 3.10, using Orange Data Mining version 3.36.0 to run the official FCBF implementation and NaiveBayesLearner. Custom code orchestrated repeated stratified cross-validation, training-fold imputation and standardization, out-of-fold probability aggregation, bootstrap confidence intervals, and performance summaries. The pinned numerical environment included NumPy 1.23.5, SciPy 1.10.1, scikit-learn 1.1.3, and pandas 1.5.3.

Study Population

Of the 109 women included, 22 (20.2%) had normal BMD, 66 (60.6%) had osteopenia, and 21 (19.3%) had osteoporosis. Mean age was 65.7±10.4 years and mean BMI was 29.8±5.5 kg/m². Estrogen-free years were available for 98 participants and missing for 11 participants (10.1%). Baseline characteristics are summarized in Table 1.

Radiomic Feature Selection and Stability

Exactly 930 radiomic features were extracted from each standardized patellar ROI. In the revised analysis, FCBF feature selection was repeated independently within each training fold, with a maximum of 10 radiomic features retained for the principal models. Across the 200 training folds generated by 20 repetitions of 10-fold cross-validation, feature selection showed substantial variability. In the three-class analysis, 220 different radiomic features were selected at least once; the most frequent feature, log-sigma-0-5-mm-3D_gldm_DependenceEntropy, was selected in 129/200 folds (64.5%). In the binary analysis, 191 different features were selected at least once; the most frequent feature, log-sigma-0-5-mm-3D_glszm_ZoneVariance, was selected in 180/200 folds (90.0%). No feature was selected in every training fold in either analysis, arguing against a single stable radiomic signature in this cohort.

Sensitivity analyses using feature caps of 5, 20, and 30 showed that increasing the number of retained features did not result in a monotonic improvement in discrimination (Table 2). In the binary radiomics-only analysis, the AUCs were 0.598, 0.653, 0.657, and 0.616 with maximum feature counts of 5, 10, 20, and 30, respectively. In the primary three-class analysis, the osteoporosis-class AUCs for the corresponding radiomics-only models were 0.625, 0.660, 0.661, and 0.655.

Primary Three-class Analysis

In the primary three-class analysis, discrimination for the osteoporosis class was modest across all internally validated models (Table 3). The clinical-only model yielded an osteoporosis-class AUC of 0.570 (95% CI: 0.417-0.714). The radiomics-only model increased the AUC numerically to 0.660 (95% CI: 0.532-0.773), whereas the addition of BMI did not improve discrimination (AUC: 0.647; 95% CI: 0.518-0.763). The age+BMI+radiomics model yielded an AUC of 0.666 (95% CI: 0.541-0.784), and the full clinical+radiomics model yielded an AUC of 0.668 (95% CI: 0.543-0.784), with a sensitivity of 28.6% and a specificity of 86.4%.

Overall, three-class performance remained limited. The full clinical+radiomics model achieved an overall accuracy of 49.5%, balanced accuracy of 41.8%, macro-F1 score of 0.417, and linearly weighted Cohen’s kappa of 0.145. By comparison, the majority-class reference classifier, which assigned every participant to the osteopenia group, achieved a higher overall accuracy of 60.6% but a lower balanced accuracy of 33.3%, a macro-F1 score of 0.251, and a weighted kappa of 0. These results indicate that the model provided some class-balanced information beyond the majority-class rule but did not improve the accuracy of exact three-category classification.

Class-wise Performance

Class-wise evaluation of the full clinical + radiomics model further demonstrated the limitations of exact three-category classification (Table 4 and Figures 4 and 5). Discrimination was highest for osteoporosis (AUC: 0.668; 95% CI: 0.543-0.784), whereas discrimination for normal BMD and osteopenia was close to chance level (AUCs: 0.549 and 0.500, respectively). The model correctly classified 8 of 22 participants as normal, 40 of 66 as having osteopenia, and 6 of 21 as having osteoporosis.

Secondary Binary Analysis

When normal BMD and osteopenia were combined into a non-osteoporosis category, the clinical-only model showed limited discrimination, with an AUC of 0.565 (95% CI: 0.414-0.708). The radiomics-only model achieved an AUC of 0.653 (95% CI: 0.526-0.774), and the full clinical+radiomics model yielded an identical AUC of 0.653 (95% CI: 0.525-0.773) (Table 5 and Figure 6). Using the model-derived classification rule, the radiomics-only model identified 5 of 21 osteoporosis cases, corresponding to a sensitivity of 23.8% and a specificity of 78.4%. The full clinical+radiomics model identified 6 of 21 cases of osteoporosis, with a sensitivity of 28.6% and a specificity of 76.1%. The paired patient-level bootstrap comparison showed no evidence of incremental discrimination after adding the clinical variables to radiomics (Delta AUC: 0.000; 95% CI: -0.050 to 0.048).

Sensitivity Analysis Related to Menopausal Status

Because a value of 0 for estrogen-free years could represent either premenopausal status or less than 1 year since menopause, a conservative sensitivity analysis excluded all five participants with that value. The resulting cohort comprised 104 women (normal BMD, n=20; osteopenia, n=63; osteoporosis, n=21). The binary findings were materially unchanged. The radiomics-only model yielded an AUC of 0.665 (95% CI: 0.543-0.785), whereas the full clinical+radiomics model yielded an AUC of 0.661 (95% CI: 0.533-0.785). The paired difference remained negligible (Delta AUC: -0.004; 95% CI: -0.044 to 0.034). Thus, the principal binary findings were not driven by participants with zero-coded estrogen-free years.

This exploratory proof-of-concept study evaluated whether quantitative radiomic features extracted from a standardized patellar region on routine lateral knee radiographs contain information related to DXA-defined BMD status. After repeating the analysis using a leakage-free repeated internal validation framework with fold-specific feature selection, three main findings emerged. First, patellar radiomics showed modest discrimination of osteoporosis, with an AUC of 0.660 in the primary three-class setting and 0.653 in the secondary osteoporosis-versus-non-osteoporosis analysis. Second, the addition of age, BMI, and estrogen-free years did not improve discrimination beyond radiomics for the binary task. Third, exact three-category classification remained limited. Although the combined model performed better than the majority-class reference on class-balanced measures, its overall accuracy remained below the 60.6% obtained by assigning all participants to the osteopenia category. These findings support a measurable, albeit limited, radiographic texture signal and emphasize the need for cautious interpretation in a small, internally validated cohort.

The primary three-class task was particularly challenging. The full clinical+radiomics model yielded an osteoporosis-class AUC of 0.668; however, its overall accuracy was 49.5%, compared with 60.6% for the majority-class baseline. The apparently discordant results across metrics are informative. Overall accuracy was strongly influenced by the predominance of osteopenia, whereas balanced accuracy, macro-F1, and weighted Cohen’s kappa placed greater weight on performance across all three categories. The combined model improved balanced accuracy from 33.3% for the majority-class reference to 41.8% and macro-F1 from 0.251 to 0.417, suggesting that it captured some information beyond class prevalence. Nevertheless, these values remain insufficient for reliable individual-level classification. Osteopenia is an intrinsically broad, intermediate BMD category, and its overlap with adjacent normal and osteoporotic groups may make precise separation into three classes more difficult. Measurement variability in the reference standard, class imbalance, and limited sample size may have contributed.

The secondary binary analysis provided a more direct assessment of osteoporosis discrimination. Radiomics alone yielded an AUC of 0.653, compared with 0.565 for clinical variables alone. Bone architecture adapts to mechanical loading (22), and age, BMI, and age at menopause are associated with osteoporotic fracture risk (23). Adding age, BMI, and estrogen-free years did not increase the AUC (0.653), and the paired bootstrap estimate of the difference was essentially null (Delta AUC: 0.000; 95% CI: -0.050 to 0.048). Although sensitivity increased from 23.8% to 28.6% after the addition of the clinical variables, the increase occurred without improvement in overall discrimination and was accompanied by a reduction in specificity from 78.4% to 76.1%. Accordingly, the present data do not support a reproducible incremental contribution of the evaluated clinical variables beyond the radiomic features. In particular, the low sensitivity observed in the model-derived classification rule limits the suitability of the present models for osteoporosis screening or DXA referral in clinical practice.

An additional finding was the instability of individual radiomic features. Despite restricting the principal analysis to a maximum of 10 features, 220 distinct features were selected at least once across the 200 training folds in the three-class analysis, and 191 were selected at least once in the binary analysis; no single feature was retained in every training fold. The most frequently selected features were retained in 64.5% of three-class folds and 90.0% of binary folds. This variability is consistent with the high dimensionality of the original 930-feature space relative to the number of participants and osteoporosis cases, and indicates that the current data do not support a fixed patellar radiomic signature. Importantly, increasing the maximum number of selected features did not consistently improve performance. In the binary task, radiomics-only AUCs were 0.598, 0.653, 0.657, and 0.616 when the feature cap was set at 5, 10, 20, and 30, respectively. Thus, the results favor a parsimonious modeling strategy while also illustrating substantial uncertainty in feature-level interpretation.

The present findings should be interpreted within the broader literature on opportunistic osteoporosis assessment. Quantitative CT, CT radiomics, and deep-learning approaches applied to routine imaging have demonstrated the potential to derive bone-health information from examinations obtained for other indications (12-16,24). Radiograph-based investigations of the patella, knee, foot, and ankle have similarly suggested that conventional radiographs contain information associated with BMD status (17-19). Clinical covariates have also been incorporated into lumbar radiograph-based models (25). A systematic review of deep-learning models using plain radiographs reported promising diagnostic performance but also emphasized the limited number of available studies and the need for further clinical validation (26). More recently, multi-site validation of hip and pelvic radiographs has highlighted the importance of testing radiographic osteoporosis models across scanners, patient subgroups, and clinical environments (27). Direct numerical comparison with the present results is difficult because these studies differ substantially in anatomical site, outcome definition, model architecture, cohort size, and validation strategy.

The standardized 15-mm patellar ROI remains a practical feature of the present approach. It provides a simple and anatomically predefined sampling strategy and avoids the segmentation complexity of whole-bone or end-to-end models. Nevertheless, gray levels on plain radiographs are not calibrated measurements of physical density. Detector characteristics, exposure parameters, image processing, pixel spacing, and positioning can influence the apparent texture from which radiomic features are derived. Image normalization and gray-level discretization may reduce some sources of intensity variation, but cannot fully eliminate acquisition dependence. The absence of spatial resampling in the present analysis may have introduced additional variability in texture measurements. These considerations are especially important, given the observed feature-selection instability, and underscore the need for multi-vendor and external validation before biological meaning is assigned to individual features or clinical use is considered.

Menopausal status coding also required careful interpretation. Estrogen-free years were unavailable for 11 participants, and a value of 0 represented either premenopausal status or a period of less than 1 year since menopause. Because these states could not be reliably distinguished retrospectively, we performed a conservative sensitivity analysis, excluding all five participants with a value of 0. Binary discrimination was materially unchanged: the radiomics-only and full models yielded AUCs of 0.665 and 0.661, respectively, with no evidence of incremental benefit from the clinical variables (Delta AUC: -0.004; 95% CI: -0.044 to 0.034). Thus, the principal binary findings were not driven by participants with potentially premenopausal or very early postmenopausal status. Nevertheless, prospective studies should record menopausal status explicitly and use age-appropriate DXA interpretation.

Study Limitations

This study has several limitations. First, it was retrospective and single-center, with only 109 participants and 21 osteoporosis cases, resulting in substantial uncertainty in model performance estimates despite repeated stratified cross-validation. Second, 930 radiomic features were initially available, creating a high-dimensional setting relative to the sample size. Although feature selection was confined to each training fold and the principal model was deliberately restricted to a maximum of 10 radiomic features, feature selection stability remained limited. Third, all performance estimates were derived from internal validation; an independent external test cohort was not available. The results should therefore be considered exploratory rather than evidence of clinical generalizability.

Fourth, the three-class models showed limited performance and did not exceed the majority-class benchmark in overall accuracy. The low sensitivity of the binary models under the model-derived classification rule also precludes use as a clinical screening or DXA-referral system in their present form. Fifth, a single experienced radiologist performed ROI placement; intra- and interobserver reproducibility were not assessed. Although the fixed-size, anatomically predefined ROI was intended to simplify placement, its reproducibility requires direct evaluation. Sixth, radiographs were obtained during routine clinical practice; no spatial resampling was performed; variation in detector characteristics, exposure, positioning, post-processing, and pixel spacing may have influenced texture features. Detailed per-examination acquisition parameters were not consistently retained in the analytic dataset.

Seventh, data on estrogen-free years were unavailable for 11 participants, and zero-coded values could not distinguish premenopausal women from those within their first postmenopausal year. The conservative exclusion sensitivity analysis reduced concern that these participants drove the main binary results but did not replace prospectively recorded menopausal status. Eighth, the reference classification was based on the lowest T-score documented in the finalized clinical DXA report. Raw DXA scans, scanner-specific quality-control records, and site-level adjudication of degenerative or invalid lumbar measurements were not available in the analytic dataset for retrospective re-review; this also limits interpretation for potentially premenopausal participants, for whom Z-scores are generally preferred. The present study evaluated a predefined Naive Bayes-based proof-of-concept framework rather than conducting an extensive comparison of machine-learning algorithms. Future work should prioritize larger multicenter cohorts, standardized acquisition and preprocessing, reproducibility assessment, independent external validation, and prospective evaluation before consideration of clinical implementation.

Patellar radiomic features extracted from routine lateral knee radiographs showed modest discrimination for DXA-defined osteoporosis under leakage-free, repeated internal validation. Clinical variables did not demonstrate an incremental AUC benefit beyond radiomics in the binary analysis, and the exact three-class BMD classification remained limited. The substantial variability of selected radiomic features and uncertainty in performance estimates indicate that the present findings should be regarded as exploratory evidence for a radiographic texture signal rather than a clinically applicable screening model. Larger externally validated studies are required to determine whether this signal can be developed into a robust opportunistic approach for osteoporosis assessment.