Mammography-Based Radiomics for Prediction of Nodal Status and Disease Burden in Breast Cancer: A Temporally Validated Study

1. Introduction

Breast cancer (BC) is the most common cancer in women (excluding non-melanoma skin cancer), with 361,859 new cases and 92,749 deaths in the European Union in 2024 [1]. Mammography (MMG) is the most widely used breast imaging modality worldwide. In many countries, it is part of screening programs aimed at early detection of BC. In the United States, approximately 39 to 40 million MMG procedures are performed annually [2]. Given its reported sensitivity of approximately 72%, relatively low cost, and high availability, it remains the cornerstone of BC screening [3].

Radiomics focuses on extracting quantitative features from medical images. It is gaining momentum in oncology, as the field continually seeks novel and more robust predictive methods. The main promise of radiomics is that it captures what is imperceptible to the human eye, such as subtle variations in shape, texture, and signal intensity [4]. These features, when combined with clinical, histopathological, and genetic variables, can be used to develop predictive models.

Radiomics has been investigated in BC since the field’s inception. Initial reports appeared in 2015, followed by numerous studies [5]. Radiomics and multi-omics approaches have been successful in predicting multiple characteristics of BC, including hormone receptor and HER2 status, response to treatment, and risk of recurrence, among others [6,7,8]. However, the majority of BC radiomics and multi-omics studies use magnetic resonance imaging (MRI)-derived features. While MRI offers superior soft tissue contrast and three-dimensional evaluation of the breast, it is not routinely used as the primary diagnostic imaging modality in the general population.

In this study, we aimed to evaluate whether MMG-based radiomics, leveraging a more widely available modality, can predict the disease burden, concordance with radiologist-assessed BI-RADS category and the status of PIK3CA mutation in BC.

2. Materials and Methods

2.1. Patient Cohort

In this single-center retrospective study, we screened 136 BC cases from 134 female patients treated at Holy Cross Cancer Center, Kielce, Poland. Due to the retrospective nature of the study, with all data collected during the diagnostic and treatment process, the institutional review board approval was not required. The inclusion criteria were histopathological confirmation of BC and the presence of pre-treatment digital mammogram in the hospital’s archives. Exclusion criteria were lack of histopathological confirmation, lack of mammograms, unavailable radiologist’s report for the MMG, low quality of imaging or faulty genetic testing due to sample degradation. After the initial screening, 34 cases from 34 patients were excluded from the study. Utilized methods are reported per TRIPOD+AI guidelines [9]. Figure 1 demonstrates the inclusion process.

2.2. Endpoints

The primary endpoint was the evaluation of prognostic performance of MMG-derived radiomic features for prediction of the presence of nodal metastases in the pre-surgical diagnostic process (cN status). Secondary endpoints were the radiomics-based prediction of clinical tumor stage (cT), overall clinical disease stage, presence of PIK3CA mutation and concordance with radiologist-assessed BI-RADS 4 and 5 groups. The methodological workflow of this study is summarized in Figure 2.

2.3. Image Acquisition

Mammograms were taken using the following MMG hardware: Siemens Mammomat Inspiration, with kilovoltage peaks up to 35 kV, a tungsten anode and a 0.05 mm rhodium filter and Hologic Selenia Dimensions and Lorad Selenia systems with kilovoltage peaks up to 33 kV, molybdenum anodes and 0.03 mm rhodium filters. 76.7% of the mammograms were performed using the Siemens system, and the remaining 23.3% were taken using the Hologic systems. Open source, Image Biomarker Standardisation Initiative (IBSI)-compliant Z-Rad software (University Hospital Zurich, Switzerland) was used for image and mask preprocessing and radiomic features extraction [10]. All the images and masks were preprocessed, with resample resolution of 0.1 mm, linear image interpolation and nearest neighbor mask interpolation. The radiomic features were extracted from NIfTI files, with 2D averaged texture aggregation and fixed bin number discretization – 36 bins for mask1 and 22 for mask2 (based on IQR of pooled tumor and breast tissue intensities in the training subset, applied unchanged to test). 138 radiomic features were extracted for every segmentation. The full list of extracted radiomic features is provided in the Supplementary Materials.

a. 

Image segmentation

Mammograms of ‘tumor-side’ and ‘healthy’ breasts in craniocaudal (CC) and mediolateral oblique (MLO) projections were collected. Segmentations were made for the tumor region (mask1) and for the whole breast, excluding skin and nipple (mask2) for the ‘tumor-side’ breast, and a whole-breast segmentation for the ‘healthy’ breast. The segmentations were performed using 3D Slicer software and saved as a binary labelmap in NIfTI format for later processing [11]. All segmentations were performed manually. Segmentation of the tumor region was performed according to a written radiologist’s report in all cases.

a. 

Statistical analysis

For model validation, cases were divided into training (70%) and test (30%) sets with a temporal split based on the histopathological diagnosis date, with the most recent cases in the test set. This approach simulates a prospective deployment. The split was fixed throughout analyses. Random seed was fixed at 42 for reproducibility.

We evaluated 10 feature subsets to assess the extent to which they contribute to the predictive performance of MMG-based radiomics. The subsets included: tumor ROI-based features (all features and exclusively texture features, separately for CC, MLO views and combined), contralateral asymmetry features, a combination of tumor features and asymmetry features, whole breast mask features, as well as all features combined. For each of the described representations, we trained 6 linear models. The linear models included: L1-regularized logistic regression as selector with L2-regularized logistic regression or elastic net classifier; F-score univariate filter with L2-regularized logistic regression or elastic net classifier; mutual information univariate filter with L2-regularized logistic regression classifier, and an F-score univariate filter with a linear-kernel calibrated support vector machine classifier. The non-linear models (random forest, gradient boosting) were excluded during the initial analyses due to severe overfitting at the available sample size. For all the models, we employed median imputation for missing values, variance thresholding, feature standardization with the Z-score, and feature selection (limited to k = 2-5 features). Hyperparameter tuning was performed using a grid search within repeated stratified k-fold cross validation (CV; 5 folds x 3 repeats) on the training set. The winning model was the model with the highest mean CV area under the receiver operating characteristic curve (ROC AUC) for each endpoint. Both tumors from the same patient could appear in the same CV fold. Given the negligible number of bilateral cases (2/100), no patient-level grouping was applied.

The winning model was applied to the temporally split test set. The achieved discrimination was assessed using the AUC. Confidence intervals (CI) were calculated using bias-corrected, accelerated bootstrapping with 2000 resamples. Statistical significance was calculated using permutation testing (1000 permutations, yielding a minimum p-value of 0.001), and endpoints with permutation p-values of less than 0.05 were considered statistically significant. Benjamini-Hochberg correction with false discovery rate of 0.05 was applied across the five endpoints [12]. The winner model with its optimal hyperparameters was refitted on the entire dataset to obtain the final set of features and their coefficients.

In the paper, we report the generalization performance of the test set. We also assessed the calibration of the model by calculating the Brier score, intercept and slope. After that step, the model was compared to a clinical baseline model (consisting of age at diagnosis, tumor size and BI-RADS group, which were considered the baseline features available at the time of the MMG reporting). The clinical model was trained using L2-regularized logistic regression. The comparison was performed using the DeLong test for paired AUCs. BI-RADS group was excluded from the analysis for BI-RADS endpoint to prevent leakage. To address selection variability, the winning pipeline was refitted on each of the 15 CV folds for cN, cT and advanced stage endpoints, and feature selection frequency was reported. Features selected in ≥ 80% of folds are considered stable. To assess the stability of test set performance, the evaluation was repeated 200 times on random 80% subsamples of the test set. Median AUC and interquartile range were reported.

All analyses were performed in Python 3.13 with scikit-learn 1.8, statsmodels, SciPy, and matplotlib libraries [13,14,15,16]. Claude Code (artificial intelligence agent; Anthropic, San Francisco, USA) was used for code debugging.

3. Results

3.1. Cohort Baseline Characteristics and Sample Size

After applying the inclusion and exclusion criteria, 102 cases from 100 patients were included in the analysis. The median patient age in the studied group was 63 years (interquartile range, IQR: 49–70 years), while the median age at diagnosis was 59 years (IQR: 47-67 years). 42% cases were positive for PIK3CA mutation, 49% had clinically suspicious lymph nodes assessed using ultrasound or MRI (cN positivity), 65.6% were classified as BI-RADS 4 and 34.4% as BI-RADS 5 (BI-RADS 6 excluded from analyses as it stands for pathologically confirmed cancer), 41.6% were stage I or II and 58.4% were stage III or IV at diagnosis. 40.2% were clinically staged as cT1 or cT2, while 59.8% as cT3 or cT4. Baseline characteristics of the cohort, including the remaining hormonal, molecular and clinical features, with the split into training and test subsets are presented in Table 1.

We calculated the events-per-variable (EPV) index using the smaller class as the event count. At prespecified 3 maximum features, EPV was 11.3 for cN status, 10.3 for PIK3CA status, 9.7 for clinical stage, 8.3 for cT, and 7.0 for BI-RADS group. It therefore fell below TRIPOD-recommended 10 for a part of the endpoints at k=3 and for all endpoints at k=5. This issue is further addressed in the limitations section.

3.2. Primary Endpoint: Prediction of cN Status

The primary results of the analysis are presented in Table 2. MMG-based radiomics were effective for prediction of the clinical nodal status, clinical tumor stage and overall disease stage in this cohort. The model failed to discriminate BI-RADS 4 vs 5 groups and the presence of PIK3CA mutation in the studied group. For the primary endpoint, clinical nodal status prediction (cN0 vs cN1+), the best feature subset was texture features from the MLO projection for the tumor ROI, and the best statistical model was model-based feature selection with L1-regularized logistic regression followed by classification with L2-regularized logistic regression. The CV AUC was 0.684, with the test AUC being 0.726 (95% CI: 0.50–0.89), representative of an acceptable performance level at this sample size. Statistical significance was confirmed by permutation testing, with a p-value of 0.01 and Benjamini-Hochberg adjusted q-value of 0.017. Features selected for prediction of cN status were: the GLCM Information Correlation 1 (cm_info_corr1, β=−0.052), the GLDZM grey level variance (dzm_gl_var, β=+0.049), the GLDZM zone distance entropy (dzm_zd_entr, β=+0.137), and the NGLDM low dependence low grey level emphasis (ngl_ldlge, β=−0.103).

3.3. Prediction of cT Status and Overall Disease Stage

We observed stronger results for two of the secondary endpoints: the differentiation of early vs advanced clinical tumor status (cT1-2 vs 3-4) and clinical stage (I-II vs III-IV). For the cT status, the CV AUC was 0.824 and the temporally split test set AUC was 0.792 (95% CI: 0.58–0.93). Similarly, for the advanced stage prediction, CV AUC was 0.763 and test set AUC was 0.778 (CI: 0.57–0.92). Both endpoints were stable through the BH FDR correction. Radiomic features selected for the positive endpoints are presented together with their coefficients and prevalence across folds are shown in Table 3. Interestingly, for clinical stage prediction, the refit features were not consistently selected during CV (<13%), but intensity-based energy (stat_energy) was selected in 14/15 folds. This likely results from L1 feature selection instability.

3.4. Prediction of BI-RADS Group Concordance and PIK3CA Mutation

Remaining two of the secondary endpoints were clearly negative: differentiation between BI-RADS 4 and 5 groups (test AUC: 0.300) and prediction of PIK3CA status (test AUC: 0.432). A detailed table of calibration metrics for all models is provided in Table 4. ROC curves for all the endpoints are shown in Figure 3.

3.5. Comparison with the Clinical Baseline

The results were compared to a clinical model, which included the quantifiable features available at the point of MMG acquisition and reporting (age + BI-RADS + tumor size). The model was evaluated for cN and clinical stage only – cT was excluded because it is itself a derivative of tumor size. For cT, the clinical model included only age and BI-RADS category. The radiomic model outperformed the clinical baseline in all three positive endpoints (Table 5), however the differences were not significant. The study was not powered to detect AUC differences between the radiomic models and a clinical baseline at the magnitudes observed.

4. Discussion

The majority of radiomics research in BC is based on features derived from MR imaging. MRI provides much superior, three-dimensional visualization of the breast and tumor tissue. MRI-based radiomics predict complex outcomes, such as response to chemotherapy, presence of the triple-negative subtype, assessing receptor status and Ki-67 expression, among others [5,17,18,19]. However, the modality is much less frequently used in clinical practice and rarely employed as screening or an initial diagnostic method. By contrast, MMG is widely and easily available around the world, and most likely to be the first imaging performed when a suspicious breast mass is reported by the patient or identified by the physician during physical examination. Therefore, reliably extracting useful clinical data from a widely available imaging modality would provide a much-needed solution for risk stratification and triage of patients, with the prioritization of further diagnosis in those with higher suspicion of advanced or metastatic disease, especially in the age of BC overdiagnosis and possible overtreatment [20,21]. The body of literature, while present, is more limited for MMG than for other modalities [22,23].

In this study, we focused on the prediction of disease burden and PIK3CA status based on MMG-based radiomics. Among the studied endpoints, three – prediction of cN status, cT status and overall disease stage – showed statistically significant temporally validated discrimination after BH-FDR correction. The temporal split, fixed throughout the analyses, with the most recently diagnosed 30% of cases assigned to the test subgroup simulated prospective deployment.

Prediction of the cN status using MMG-based radiomics was chosen as a primary endpoint for its clinical relevance to the disease risk stratification and choice between surgical intervention and neoadjuvant chemotherapy. Clinical status, based on either ultrasound or MRI, was chosen instead of the surgical staging because the majority of the patients included in the study underwent neoadjuvant chemotherapy, often achieving postoperative partial or complete remission. Therefore, the surgical staging does not represent the initial nodal status at the point of MMG acquisition. The radiomic model achieved a temporally validated test-set AUC of 0.726. Four leading features for the predictive power were related to tumor texture, derived from the MLO projection. They encode the intratumoral heterogeneity, and we hypothesize that this mammographic heterogeneity could reflect the underlying tumor biology associated with nodal involvement, potentially including stromal disorganization and microcalcification patterns previously linked to advanced disease [24,25]. Further work, optimally also encompassing radiopathomics, is needed to explore this theory. Previous research, including the work by Tan et al. and Yang et al., among others, confirms the ability of MMG-based radiomics to predict lymph node metastases [26,27]. Although our test set AUC is lower than several of the previously reported MMG-based nodal status prediction models, this likely results from the more conservative temporal-split validation we applied. An external cohort and a larger calibration set will be required before the model can be considered clinically useful.

Prediction of cT and the overall clinical stage by MMG-based radiomics is, by its definition, partly explained by tumor size: tumor diameter, which can be measured directly on a mammogram, defines the cT category and is a major component of clinical stage. GLDZM large distance high grey level emphasis, the dominant feature for the cT model, rewards tumor regions that have many high-intensity zones located deep inside the lesion. Therefore, this feature carries additional signal intensity and tumor texture components. The clinical stage prediction model retained a related set of intra-tumoral intensity and zone-distance descriptors (global intensity peak, GLDZM small distance emphasis), which carry additional signal intensity and texture information. This pattern can serve as evidence that MMG-based radiomics can capture the biological phenotype of the disease and that some features track tumor aggressiveness in the same sense that the pathologist's nuclear grade and mitotic index do.

Two of the endpoints – the prediction of PIK3CA mutation and radiologist-assigned BI-RADS group – failed, with large CV-test AUC differences typical of overfitting and low EPV. Nevertheless, the negative results are retained deliberately, not to overestimate the potential yield of MMG-based radiomics. Interestingly, recent work by Yao et al. focused on the prediction of PIK3CA mutation in BC and found positive results (AUC 0.718) for the prediction of this genetic alteration, and even stronger results for ultrasound-based radiomics and a combined model [28]. This contradicts our findings and could be either a result of a larger sample (186 patients) or different statistical methodology. Further development of PIK3CA mutation presence prediction is key for earlier qualification for treatment with alpelisib in patients with advanced BC (29). Nevertheless, the majority of the research on MMG-based radiomics still focuses on fundamental cancer risk stratification [30].

We acknowledge the limitations of our study. First, the sample size is small given the number of variables produced by the radiomics software; EPV fell below the TRIPOD-recommended threshold, and the study was not powered to separate between-model AUC differences of the magnitude observed. The study is single-center, with no external validation cohort; however, the temporal split partially addresses this issue. Further limitations are the single-reader segmentation and suboptimal calibration of two positive endpoints: cN prediction, with slope of 3.762 suggesting overconfident extreme predictions, and cT prediction, with a slope of 0.444, typical of compressed predictions. A recalibration on a larger, external cohort would benefit the clinical applicability of the results. The cohort included mammograms acquired on two systems (Siemens and Hologic). This heterogeneity arguably enhances the generalizability of our findings beyond a single-scanner setting; however, the inter-vendor radiomic variability was not formally addressed. Our IBSI-compliant pipeline applies standardized preprocessing and discretization to all images to partially mitigate this issue.

We treat the pre-specified pipeline, and the blinded, leakage-proof temporal test set as significant methodological strengths. IBSI-compliant feature extraction, multiple-testing correction across endpoints and retention of negative endpoints also strengthen the results. A larger, externally validated study is required to confirm our findings, in particular regarding the cN status prediction, and to properly compare them to the clinical baseline.

5. Conclusions

Our study provides evidence that radiomic data extracted from mammography can discriminate between cN0 and cN-positive disease; that tumor texture features change with an increasing clinical tumor stage, and that radiomic features correspond to overall disease burden in BC. Conversely, our study found no correlation between radiomic features and PIK3CA mutational status in BC patients and showed no concordance between radiomic features and radiologist-assessed BI-RADS group. However, further - larger and externally validated - studies are needed to establish clinical utility and possibly examine other domains that can be predicted by MMG-derived radiomics.