Positive margin is significantly associated with higher tumor stage and lumpectomy

The probability of attaining positive margins after surgery was observed to be significantly (p ≤ 0.05) associated with higher tumor stage, larger tumor size and chest wall involvement (T4), positive lymph nodes (N2, N3), and distant metastasis (M1), based on univariable logistic regression models and Fisher’s exact test (Table 1). The type of first surgery to remove tumor also influenced margin status with lumpectomy (as reference) having significantly higher chance of obtaining positive margins than mastectomy (Simple Mastectomy: p = 0.002, Odds Ratio (OR) = 0.30, Confidence Interval (CI) = 0.13 − 0.62; Modified Radical Mastectomy: p < 0.001, OR = 0.30, CI = 0.15 − 0.57). Among PAM50 subtypes, Luminal A subtype (as reference) was observed to be significantly contributing towards positive margin in the univariable regression model compared to the Basal-like (Basal) subtype (p = 0.05, OR = 0.44, CI = 0.18 − 0.94). Her2-enriched (Her2) subtype was associated with positive margins (OR = 1.39) although it was not significant (p = 0.397). The results of Fisher’s exact test were consistent with the logistic regression results except for PAM50 subtype which did not show any association with margin status.

Table 1 Summary of clinical characteristics of TCGA-BRCA data (n = 951) and their association with margin status

The significant factors in the univariable regression model (Stage, PAM50, TNM: T = Tumor size, N = Lymph Node status, M = Metastasis, Type of first surgery) were used in the multivariable model with margin status as response variable. Tumor stage, size, and lymph node status, which were highly significant in the univariable model, were no longer significant in the multivariable model (Supplementary Table S1). Further evaluation using various multivariable models proved that Stage and TNM were confounding (Supplementary Table S2); hence, only Stage was used in the final multivariable model (Table 2). The final regression model, in agreement with the univariable model, showed that patients diagnosed at higher tumor stage (Stage III: p < 0.001, OR = 4.85, CI = 2.09 − 12.41; Stage IV: p < 0.001, OR = 80.83, CI = 18.65 − 411.45) were significantly associated with positive margins. Similarly, in case of type of surgery for tumor removal, the multivariable regression model reemphasized that lumpectomy (as reference) was significantly associated with positive margin compared to simple mastectomy (p = 0.002, OR = 0.27, CI = 0.12 − 0.59) and modified radical mastectomy (p < 0.001, OR = 0.17, CI = 0.08 − 0.35). For the PAM50 subtypes, Luminal A (as reference) was significantly associated with positive margins compared to basal subtype (p = 0.042, OR = 0.41, CI = 0.16 − 0.91).

Table 2 Multivariable logistic regression analysis on the association of margin status with significant clinical features from univariable analysisEffect of margin status and other factors on disease progression

Of the 951 cases included in our study, one was excluded from survival analysis due to missing follow-up information. The univariable survival models using 950 cases for margin status showed that positive margins were significantly associated with worse survival with both PFI (p < 0.001) and OS (p = 0.006) as endpoints (Fig. 1). In addition to margin status, stage, TNM, PAM50 subtype, and hormone receptor (Estrogen Receptor (ER), Progesterone Receptor (PR)) status were significantly associated with disease progression (Table 3). While examining the survival models based on histology, mucinous carcinoma was found to have significant survival difference compared to ductal carcinoma (Table 3). However, as the sample size (n =15) and the number of events (n =3) were low for mucinous carcinoma (Supplementary Fig.S1), these results were regarded as unreliable. It is worth noting that the type of first surgery, though significantly associated with margin status, does not significantly impact survival.

Fig. 1

The Kaplan–Meier (K–M) curves for cumulative survival in years for margin status for two end points: progression-free interval (PFI) (a) and overall survival (OS) (b). P value, Hazard ratio (HR), and the number of events ‘/’ number of cases are given in the legends of plots

Table 3 Univariable survival analysis to assess the effect of each clinicopathologic factor on disease progression (Progression Free Interval, PFI)

In order to assess the combined effect of margin status and other factors that were significant in the univariable model on survival, multivariable survival analysis was performed. TNM, though significant in the univariable model, was excluded in the multivariable models since tumor stage is derived from TNM and the inclusion of both features in the same model was observed to be confounding in the previous logistic regression model. Surprisingly, PAM50 and ER status were not significant in this model (Supplementary Table S3). Further exploration using different multivariable models (Supplementary Tables S4–S5) indicated that hormone receptor status and PAM50 were confounding to each other; hence only PAM50 was retained in the final model (Table 4). Higher tumor stages (III and IV), and the Basal and Her2 subtypes were significant (p ≤ 0.05) in contribution to disease progression in the final model, while margin status was not significant (p = 0.135, HR = 1.54, CI = 0.88 − 2.70). The bi-variable survival models (Table 5) demonstrated that margin status remained highly significant when PAM50 or either of the hormone receptor (ER/PR) status was added to the model whereas in the model with tumor stage, margin status was only close to significance (p = 0.067).

Table 4 Final Cox proportional hazards regression model for multivariable survival analysisTable 5 Assessment of impact of each significant factor from univariable models on margin status using bi-variable survival analysis with Progression Free Interval (PFI) as endpointAssociation of gene expression with margin status identified 29 DEGs

To address the sample imbalance between positive and negative margins, a matched dataset (n = 142; Supplementary Table S6) was extracted from our cohort to perform unbiased molecular analyses. Principal component analysis (PCA) of matched samples using 2000 highly varying genes did not clearly cluster the samples by margin status but clustered them instead by PAM50 subtypes (Supplementary Fig. S2). Differential expression analysis between positive and negative margin cases discovered 53 upregulated and 50 downregulated DEGs and the subsequent LASSO regression selected 29 DEGs for the prediction of margin (Supplementary Table S7). The unsupervised clustering for these 29 genes demonstrated largely subtype-driven clusters (Fig. 2). We also observed two main level clusters that have different positive margin enrichment (~ 59% for left cluster, ~ 41% for right cluster, Fisher’s exact p value = 0.044). This show the genes to some degree can separate the positive margin from negative margin. Leave-One-Out Cross-Validation (LOOCV)-based prediction models with the 29 genes showed an accuracy of 0.7.

Fig. 2

Unsupervised clustering for 29 significant genes derived using LASSO regression from Deseq2 analysis for TCGA RNA-Seq data

Among the 29 genes, 16 were upregulated and 13 were downregulated in positive margin cases. It included 17 protein-coding genes, 4 pseudogenes (AC084880.1, BEND3P1, CPHL1P, AP002001.2), and 8 long non-coding RNA (LncRNA) genes (AC004947.1, AC008663.2, AC099329.2, LINC01344, SLC26A4-AS1, AF015262.1, AC114296.1, LINC00589).

Pathway analysis identified 8 differentially expressed pathways (Table 6) between positive and negative margin cases. The 7 upregulated pathways include three cell proliferation associated pathways (E2F_TARGETS, G2M_CHECKPOINT, MYC_TARGETS_V1); two cell signaling-related pathways (ESTROGEN_RESPONSE_LATE, ESTROGEN_RESPONSE_EARLY); and two immune-related pathways (INTERFERON_ALPHA_RESPONSE, TNF_SIGNALING_VIA_NFKB). The only downregulated pathway was associated with progression and metastasis (EPITHELIAL_MESENCHYMAL_TRANSITION).

Table 6 Significant pathways observed in TCGA RNA-Seq data (n = 142) using GSEA Preranked test