This research builds upon our previous studies, extending the analysis to additional molecular pathways relevant to pyroptosis in breast cancer (Sirek et al., 2024a; Sirek et al., 2024b; Sirek et al., 2024c; Sirek et al., 2024a).

The study was designed to comprehensively evaluate pyroptosis-related molecular alterations in breast cancer across five molecular subtypes. Tumor tissue and corresponding control samples were obtained from patients diagnosed with luminal A, luminal B HER2-negative, luminal B HER2-positive, non-luminal HER2-positive, and triple-negative breast cancer (TNBC). For all collected tissue specimens, a multilevel molecular analysis was performed, including genome-wide mRNA and microRNA (miRNA) expression profiling using microarray technology. Selected differentially expressed genes were subsequently validated by quantitative reverse transcription polymerase chain reaction (qRT-PCR), while corresponding protein levels were assessed with enzyme-linked immunosorbent assay (ELISA).

In the second part of the study, an additional comparative group was included, consisting of patients diagnosed with breast fibroadenoma who were qualified for cryoablation. Peripheral blood samples were collected from these patients twice: immediately before the procedure and again at 7 days post-cryoablation. Gene expression in blood samples was quantified using qRT-PCR, and circulating protein levels were measured by ELISA.

The analyses performed in the fibroadenoma group focused exclusively on those genes and encoded proteins that, in the breast cancer microarray experiments, demonstrated consistent differential expression between tumor and control tissues regardless of cancer subtype. This allowed for the evaluation of whether early molecular responses observed in circulating blood following cryoablation reflect key pyroptosis-related pathways identified in breast cancer tissue.

Total RNA was isolated from tissue samples using TRIzol reagent (Invitrogen, Carlsbad, CA, United States; Cat. No. 15596026) following the manufacturer’s protocol. Commercially available PAXgene Blood RNA kit (Qiagen, Valencia, CA, United States, Cat No. 762174) was used to extract the total RNA from the whole blood as per the manufacturer’s protocol. Purification of RNA was subsequently performed with the RNeasy Mini Kit (QIAGEN, Hilden, Germany; Cat. No. 74104). To eliminate residual genomic DNA, samples were treated with DNase I (Fermentas International Inc., Burlington, ON, Canada; Cat. No. 18047019).

RNA integrity was verified using 1% agarose gel electrophoresis stained with 0.5 mg/mL ethidium bromide, enabling visualization of rRNA bands. RNA concentration and purity were quantified spectrophotometrically by measuring absorbance at 260 nm.

Genes associated with pyroptosis were identified using the GeneCards human gene database (accessed 10 November 2025) by searching the keyword “pyroptosis.” This search yielded 877 mRNAs that met the inclusion criteria.

Expression analyses of these genes in tumor versus control samples were performed using the HG-U133_A2 microarray platform (Affymetrix, Santa Clara, CA, United States) together with the GeneChip™ 3′IVT PLUS Reagent Kit (Affymetrix; Cat. No. 902416). The analytical workflow followed previously described and manufacturer-recommended procedures. Among the 22,277 probes present on the microarray, 65 were specific to the pyroptosis.

The protocol included synthesis of double-stranded complementary DNA (cDNA), in vitro transcription to amplify RNA (aRNA), fragmentation, and hybridization to the arrays. Fluorescence detection was performed using the Affymetrix GeneArray Scanner 3000 7G, with data acquisition via GeneChip® Command Console® Software. This approach enabled a detailed assessment of expression alterations across pyroptosis-associated transcripts.

To explore miRNAs potentially regulating pyroptosis-associated gene expression, samples were analyzed using the GeneChip miRNA 2.0 Array (Affymetrix). All procedures adhered strictly to the manufacturer’s instructions to ensure reproducibility and analytical reliability.

Differentially expressed miRNAs were identified by comparing tumor and control tissues. Potential miRNA–mRNA regulatory interactions were predicted using TargetScan (http://www.targetscan.org/) (Agarwal et al., 2015) and miRanda/mirDB (http://mirdb.org) (Chen and Wang, 2020). These databases provide computational models predicting miRNA binding based on sequence complementarity and evolutionary conservation (Chen and Wang, 2020; Liu and Wang, 2019).

For the purposes of this study, targets with prediction scores >80 were considered highly reliable. Predictions <60 were treated with caution and interpreted only in conjunction with additional evidence. Integrating data from multiple prediction tools strengthened confidence in inferred miRNA–mRNA interactions and helped delineate regulatory networks associated with pyroptosis-related signaling.

Selected genes identified through microarray analysis were validated using quantitative reverse transcription polymerase chain reaction (qRT-PCR). Experiments were conducted with the SensiFast SYBR No-ROX One-Step Kit (Bioline, London, UK), following the manufacturer’s protocol.

Gene expression levels were calculated using the 2^−ΔΔCT method, where a fold change of 1 represented the control group, values >1 indicated upregulation, and values <1 indicated downregulation. β-actin (ACTB) served as the endogenous control. Primer sequences used for amplification are provided in Table 1.

To quantify the overall activation of pyroptosis-related signaling across breast cancer molecular subtypes, a Pyroptosis Index (PI) was constructed based on the expression profiles of nine core pyroptosis-associated genes that were consistently and significantly dysregulated in tumor tissue compared with matched controls, irrespective of subtype. These genes were identified from the microarray dataset as the most robust subtype-independent markers and included: CXCL8, BAX, CASP1, CASP9, TP53, MMP9, BCL2, CDKN1A, and CDKN1B.

Based on their established biological functions, genes were categorized into two functional groups reflecting opposing roles in the regulation of pyroptotic and inflammatory cell-death pathways. The pro-pyroptotic group comprised CXCL8, BAX, CASP1, CASP9, TP53, and MMP9, which promote inflammasome activation, caspase signaling, and inflammatory cell death. The anti-pyroptotic (pyroptosis-suppressive) group included BCL2, CDKN1A, and CDKN1B, genes involved in inhibition of programmed cell death and maintenance of cell-cycle control.

To account for inter-gene variability in expression magnitude, fold-change values were first log₂-transformed using a log₂ (1 + FC) transformation to stabilize variance and reduce the influence of extreme values. Subsequently, for each gene, transformed expression values were standardized across the study cohort using Z-score normalization (mean = 0, standard deviation = 1). These standardized values were then used for index construction.

For each sample, the PI was calculated as the difference between the mean Z-score of pro-pyroptotic genes and the mean Z-score of anti-pyroptotic genes, according to the formula: PI = mean Z pro ‐ pyroptotic   genes − mean Z anti ‐ pyroptotic   genes

Higher PI values indicate stronger activation of pyroptosis-associated and inflammatory cell-death pathways, whereas lower or negative values reflect reduced pyroptotic signaling and relative dominance of survival or cell-cycle inhibitory mechanisms.

The degree of inflammasome pathway activation was quantified using an Inflammasome Activation Score (IAS), calculated exclusively from experimentally measured microarray expression data. The IAS was constructed based on a curated panel of ten inflammasome-associated genes that were directly represented on the microarray platform and passed quality-control criteria. These genes reflect key functional components of inflammasome signaling and included canonical sensors, adaptor proteins, inflammasome-regulated cytokines, signaling regulators, and downstream effector caspases: IL1B, IL18, NLRP3, PYCARD (ASC), TLR9, RIPK1, TNF, STING1, JAK3, and CASP1.

For each gene, log₂ fold-change (log₂FC) values were calculated relative to matched control tissue using the same preprocessing and normalization pipeline applied to all transcriptomic analyses in this study.

For each sample, the IAS was calculated as the arithmetic mean of the log₂ fold-change values of the ten inflammasome-related genes, according to the formula: IAS = mean log ⁡ 2 FC inflammasome   genes

Higher IAS values indicate stronger activation of inflammasome signaling, increased IL-1β and IL-18–associated inflammatory activity, and an overall intensified pro-inflammatory tumor microenvironment.

Venous blood samples were obtained from all participants by standard antecubital venipuncture using sterile, single-use equipment. Approximately 5 mL of peripheral blood was collected into serum-separating tubes (SST) containing a clot activator. The samples were allowed to clot undisturbed for 30–45 min at room temperature. Subsequently, the tubes were centrifuged at 1,500–2,000 × g for 10–15 min to separate the serum fraction. The resulting serum was carefully aliquoted into sterile microtubes to avoid repeated freeze–thaw cycles. All samples were stored at −80 °C until analysis. Quantification of target proteins was performed using commercially available ELISA kits according to the manufacturer’s protocol.

Protein concentrations of the selected targets were measured using commercially available ELISA kits, following the manufacturer’s instructions for each assay. The following kits (all from MyBioSource, Inc., San Diego, CA, United States) were employed: Human Interleukin-8 (IL-8) ELISA Kit (MBS703104); Human BCL-2–Associated X Protein (BAX) ELISA Kit (MBS701787); Human Caspase-1/ICE ELISA Kit (MBS163158); Human Caspase-9 ELISA Kit (MBS704900); Human P53 ELISA Kit (MBS041908); Human Antioncogene p21/CDKN1A ELISA Kit (MBS731903); and Human Cyclin-Dependent Kinase Inhibitor 1B (CDKN1B) ELISA Kit (MBS2122089); Human Matrix Metalloproteinease-9 (92 kDa Type IV Collagenase) ELISA Kit (MBS8415331).

Statistical analyses were performed using Statistica 13.0 PL (StatSoft, Kraków, Poland) and Transcriptome Analysis Console (Affymetrix). Normality of distribution was checked using the Shapiro–Wilk test (p < 0.05). Depending on data characteristics, group comparisons were conducted using analysis of variance (ANOVA) with Benjamini–Hochberg correction and Tukey’s post hoc test, or the Student’s t-test.

Overall survival analyses were performed with the Kaplan–Meier Plotter (http://kmplot.com/) (Győrffy, 2024a; Győrffy, 2024b). Gene–gene and protein–protein association networks were examined using the STRING database (version 11.0; accessed 20 November 2025). STRING network enrichment was evaluated by the Log10 (observed/expected) strength parameter and the false discovery rate (FDR), adjusted using the Benjamini–Hochberg method (Szklarczyk et al., 2023).

Sample size estimation was performed using an online sample-size calculator (Kalkulator doboru próby, 2022), assuming a 95% confidence interval and considering that approximately 19,620 women were diagnosed with breast cancer in Poland in 2019 (Krajowy Rejestr Nowotworów, 2019), the calculated minimum sample size required was 377 participants.

For contextual comparison, the distribution of breast cancer subtypes based on national epidemiological data (Dai et al., 2015; Jassem et al., 2020) was contrasted with their representation in our study population. Reported proportions in the literature include luminal A (23.7%), luminal B HER2– (38.8%), luminal B HER2+ (14%), HER2+ (11.2%), and TNBC (12.3%). These values correspond closely to the subtype distribution observed in our cohort.

Out of the 877 mRNAs associated with pyroptosis, ANOVA identified 48 transcripts that significantly differentiated tumor tissue from control samples (p < 0.05, log₂FC >3 or <−3). Subsequent Tukey’s post hoc test revealed that 9 mRNAs consistently distinguished neoplastic tissue from non-tumor controls regardless of the breast cancer subtype, indicating their subtype-independent discriminatory potential.

A graphical summary of these findings is presented in Figure 1, where a Venn diagram illustrates the overlap of differentially expressed mRNAs across the five breast cancer subtypes and highlights the subset of universally altered transcripts.

Microarray analysis confirmed substantial differences in the expression of nine pyroptosis-associated genes across all breast cancer subtypes compared with control tissue (Table 2). CXCL8 showed consistent and strong upregulation in all subtypes, with the highest expression in TNBC (8.18–8.21-fold) and non-luminal HER2+ tumors (7.21–7.65-fold). BCL2 was markedly overexpressed in luminal tumors (5.19–5.41-fold in Luminal A; 3.19–3.43-fold in Luminal B HER2−; 4.01–4.17-fold in Luminal B HER2+) but demonstrated significant downregulation in non-luminal HER2+ (−3.19 to −4.19) and TNBC (−4.12 to −4.22). The pro-apoptotic gene BAX displayed progressive upregulation from luminal to highly aggressive subtypes, reaching the highest levels in TNBC (9.18–9.87-fold) and non-luminal HER2+ cancers (9.12–9.18-fold). CASP1 and CASP9 showed a similar pattern of stepwise activation, with modest increases in luminal tumors (3.1–4.7-fold) and markedly elevated expression in non-luminal HER2+ and TNBC, where CASP1 reached 7.01–7.32-fold and CASP9 5.10–5.18-fold. TP53 was moderately upregulated in luminal subtypes (3.11–3.41-fold) and strongly increased in HER2+ (5.91–6.43-fold) and TNBC (7.89–7.92-fold).

Cell-cycle regulators showed divergent patterns: CDKN1A (p21) was elevated in luminal tumors (3.45–4.98-fold) but suppressed in aggressive subtypes, particularly in TNBC (−8.19). CDKN1B (p27) followed a similar trend, remaining upregulated in luminal cancers (3.11–3.99-fold) and strongly downregulated in both HER2+ (−3.14 to −6.16) and TNBC (−17.55).

Finally, MMP9 exhibited modest upregulation in luminal cancers (3.19–3.87-fold) but was dramatically increased in non-luminal HER2+ (9.19-fold) and especially TNBC (12.91-fold), consistent with its known role in invasiveness.

qRT-PCR analysis confirmed the differential expression patterns identified in the microarray experiment. CXCL8, BAX, CASP1, CASP9, TP53, and MMP9 were progressively upregulated across all breast cancer subtypes, with the strongest increases observed in non-luminal HER2+ and TNBC tumors. BCL2 showed marked overexpression in luminal subtypes but was downregulated in HER2-enriched and TNBC samples. In contrast, CDKN1A and CDKN1B exhibited subtype-specific suppression, most prominently in TNBC, consistent with loss of cell-cycle control in aggressive tumors. Overall, qRT-PCR results closely mirrored the microarray findings, confirming subtype-dependent molecular alterations (Figure 2).

Differential expression of miRNAs was determined directly from microarray data by comparing tumor and control tissues across breast cancer subtypes. Subsequently, putative miRNA–mRNA regulatory interactions were identified using in silico target prediction algorithms, and only high-confidence pairs with target scores between 89 and 97 were selected for further analysis (Table 3). Importantly, the predicted interactions were evaluated independently of expression directionality, as miRNA-mediated regulation does not necessarily result in a linear inverse relationship at the transcript level.

Prediction analysis identified several high-confidence miRNA–mRNA pairs corresponding to genes that were consistently differentially expressed across breast cancer subtypes. CXCL8, predicted to be regulated by hsa-miR-140-3p, was markedly upregulated in all molecular subtypes, with the highest expression observed in TNBC (7.65 ± 0.19). Similarly, BAX and CASP9, paired with hsa-miR-1843 and hsa-miR-124-3p, respectively, showed progressive upregulation, reaching maximal expression levels in aggressive tumor subtypes. TP53, linked to hsa-miR-300, exhibited a comparable subtype-dependent increase.

In contrast, CDKN1A and CDKN1B, predicted targets of hsa-miR-608 and hsa-miR-30 family members (hsa-miR-30d-3p and hsa-miR-30a-3p), demonstrated divergent expression patterns. While both genes remained upregulated in luminal subtypes, they were markedly downregulated in non-luminal HER2-positive and TNBC tumors, with CDKN1B showing particularly strong suppression in TNBC (−6.19 ± 0.76). These findings indicate that predicted miRNA–mRNA interactions should be interpreted as biologically plausible regulatory relationships rather than direct expression correlations.

Table 4 summarizes the estimated concentrations of nine pyroptosis-related proteins measured by ELISA in control tissue and across five breast cancer subtypes. All proteins showed significant differences compared with controls (p < 0.05). CXCL8, BAX, CASP1, CASP9, TP53, and MMP9 displayed a progressive increase from luminal tumors to the most aggressive subtypes, reaching the highest concentrations in non-luminal HER2+ and TNBC. In contrast, BCL2, CDKN1A (p21), and CDKN1B (p27) were elevated in luminal cancers but markedly reduced in non-luminal HER2+ and TNBC, reflecting subtype-specific differences in apoptotic and cell-cycle regulation.

The transcriptional response showed a distinct early activation pattern. CXCL8, BAX, CASP1, CASP9, and MMP9 demonstrated a rapid rise already at T1, reaching peak expression at T2 (8–12 h), followed by a gradual decline between T3–T6. TP53 increased moderately, with a significant elevation at T3. BCL2 showed only minimal and transient upregulation at T2–T3, while CDKN1A and CDKN1B exhibited mild, non-significant fluctuations without pronounced temporal changes. Overall, the data indicate an early inflammatory and apoptotic transcriptional peak with subsequent normalization within 1–3 months (Table 5).

Protein levels paralleled the mRNA dynamics, with the highest concentrations for CXCL8, BAX, CASP1, CASP9, P53, CDKN1A, CDKN1B, and MMP9 observed at T2, and in several cases also at T1 and T3. A progressive reduction was noted thereafter, with values approaching baseline by T5–T6. BCL2 protein remained relatively stable, exhibiting only minor elevations at early time points (Table 6).

The STRING-based PPI analysis performed on the nine subtype-independent, differentially expressed mRNAs revealed a highly interconnected molecular network (Figure 3). The network consisted of 9 nodes representing the encoded proteins and 31 edges, substantially exceeding the expected 10 interactions for a network of this size. The average node degree was 6.89, indicating that most proteins were linked to nearly seven interaction partners, reflecting strong functional connectivity. The average local clustering coefficient of 0.883 further supported the presence of densely interconnected clusters within the network. Importantly, the PPI enrichment p-value (6.31 × 10⁻⁸) demonstrated that the observed interactions were significantly more frequent than would occur by chance, confirming that these proteins are functionally related and likely engage in shared biological processes. The central positioning of TP53, BCL2, BAX, CASP1, and CASP9 highlights their coordinated involvement in regulated cell death pathways, including apoptosis–pyroptosis crosstalk (Figure 3).

The PI demonstrated a clear and progressive increase from luminal to more aggressive breast cancer subtypes, indicating a gradual shift toward enhanced inflammatory and pyroptosis-associated cell-death activity (Table 7). Luminal A tumors exhibited the lowest PI value (−0.67), consistent with minimal activation of pyroptotic signaling pathways. Luminal B HER2– tumors showed a near-neutral PI (0.13), whereas luminal B HER2+ tumors displayed a moderate increase (PI = 3.33). In contrast, markedly elevated PI values were observed in aggressive subtypes, with non-luminal HER2-positive tumors reaching a PI of 11.65 and TNBC exhibiting the highest PI (18.46), reflecting pronounced activation of pyroptosis-related molecular programs.

The IAS values followed a pattern closely paralleling that of the PI, with progressively increasing values across breast cancer subtypes of rising biological aggressiveness (Table 8). Luminal A tumors showed the lowest IAS (3.35), indicative of relatively modest inflammasome activity. Incremental increases were observed in luminal B HER2– (IAS = 3.92) and luminal B HER2+ tumors (IAS = 5.15). Substantially higher IAS values were detected in non-luminal HER2-positive tumors (7.00) and TNBC (8.41), reflecting intensified activation of inflammasome-associated signaling pathways.

This gradient corresponded with the stepwise upregulation of key inflammasome-related genes, including IL1B, IL18, NLRP3, CASP1, and additional inflammatory regulators. Detailed log₂ fold-change values for individual inflammasome-associated transcripts across molecular subtypes are provided in Supplementary Table S1.

Analysis was performed for BCL2, BAX, CASP1, CASP9, TP53, CDKN1A, CDKN1B, and MMP9 with a follow-up threshold of 60 months. The database did not contain data for the CXCL8 gene, therefore it was omitted from this analysis (Figures 4–8).

In luminal A subtype, high BCL2 expression with low TP53 and CDKN1B expression were associated with worse OS (Figure 5).

Low BAX and TP53 expression were associated with poorer OS in HER2– luminal B cancer.

In the HER2+ luminal B subtype, the change in expression of the studied genes is probably not associated with worse OS.

In turn, in non-luminal HER2+ cancer, low levels of CDKN1A and MMP9 promote worse OS.

In TNBC, low levels of CASP9, CDKN1B, and MMP9 were associated with poorer OS.