The primary tumor (T) and adjacent non-tumorous (N) tissue samples were obtained from 25 patients with TNBC undergoing surgery at Jiangxi Cancer Hospital (Nanchang, China). All paired tumor and normal tissue samples were independently identified by two pathologists. These samples were examined by quantitative real-time polymerase chain reaction (qRT-PCR) and western blotting analyses. All samples were stored at –80 °C until analysis. The tissue specimens were collected with the consent of patients from July 2019 to September 2020. The Ethics Committee of the Jiangxi Cancer Hospital approved the study design. All patients signed an informed consent form.

Three gene microarray datasets (GSE38959 [7], GSE62931 [8] and GSE64790 [9]) of expression profiles of TNBC tissues or TNBC cells and non-TNBC samples or healthy mammary ductal cells or healthy breast tissues were obtained from the GEO database (https://www.ncbi.nlm.nih.gov/geo/). High-throughput RNA sequencing data and clinicopathological data of breast cancer patients were downloaded from TCGA [10]. We confirmed the ER, PR and HER2 expression status of all breast cancer patients in TCGA based on the immunohistochemical (IHC). Information on patients with TNBC was extracted according to the classification system proposed by Voduc et al. [11]. These patients were enrolled in the TCGA-TNBC cohort.

The data were divided into TNBC and non-TNBC subsets. DEGs between the TNBC and non-TNBC samples were filtered by the Limma package in R (version 3.5.0, R Foundation for Statistical Computing, Vienna, Austria). To screen TNBC-related DEGs, the following criteria were used: false discovery rate (FDR) 0.05 and $|$log2 fold change (FC)$|$ $>$1.5. The VennDiagram package in R was used to overlap DEGs obtained from three GEO datasets. Finally, the overlapping genes were defined as DEGs.

To clarify the biological processes in which AFF3 and all DEGs, GO [12] and KEGG [13] analyses were performed using the clusterProfiler R package. Significant DEGs were identified based on the following criteria: p 0.05 and FDR q-values 0.05.

GSEA is a computational method that determines whether an a priori defined set of genes shows statistically significant, concordant differences between two biological states. In order to analyze correlations among all DEGs, clusterProfiler package was used to perform GSEA [14, 15].

GSEA firstly generated an ordered list of all DEGs according to their correlation with AFF3 expression and gene set permutations were performed 1000 times for each analysis. The criteria for significant correlations were as follows: p 0.05; FDR q-values 0.25.

Based on the TNBC-related DEGs, the L1000 fireworks display (L1000FWD) database [16] was utilized to predict prospective drugs that could attenuate or enhance the biological status of TNBC. The DEGs were submitted to the L1000FWD database for potential small-molecular drugs for TNBC. The closer the similarity score is to –1, the higher the efficacy of the drug against TNBC.

Survival analysis was conducted using survival and survminer packages. The survival duration of 152 patients with TNBC for whom detailed survival data were available was 0–9.61 years. The Kaplan–Meier method was used to draw the survival curve. Statistical significance was assessed by the log-rank test and p 0.05 was considered significant.

To further determine the effect of gene expression and clinical characteristics in TNBC patients, univariate Cox regression analysis was used to calculate the association between gene expression and clinical characteristics and patient’s overall survival (OS) in TCGA-TNBC cohort. Afterwards, a multivariate analysis was used to assess the independent prognostic factor for TNBC patient survival. The survival package in R was used to perform univariate and multivariate Cox regression.

The CIBERSORT algorithm [17], which can analyze the composition of immune cells in samples based on RNA high-throughput sequencing data, was applied to assess tumor-infiltrating immune cell (TICs) in tumor tissues of TCGA-TNBC cohort. The permutation (perm) was set at 1000 to obtain more stable results.

Radioimmunoprecipitation lysis buffer was used to lyse TNBC samples and adjacent non-tumorous tissues. The supernatant was collected after centrifugation for 10 min at 12,000 g. The bicinchoninic acid (BCA) protein quantification kit (Beijing Cwbiotech Co., Ltd., Beijing, China) measured the total protein content. After mixing the protein with the loading buffer, it was boiled for 10 min in a water bath. Protein lysates (20 µg) was separated by 10% sodium dodecyl sulfate polyacrylamide gel electrophoresis (SDS-PAGE) and transferred onto a polyvinyl difluoride (PVDF) membrane (10600023, GE Healthcare Life Sciences, Woburn, MA, USA). Nonspecific binding protein on PVDF was blocked by 10% Bovine Serum Albumin (BSA, GC305006, Servicebio, Wuhan, Hubei, China). Then, the membranes were incubated with antibodies against AFF3 (1:1000, ab106231) and $\beta$-tubulin (1:500; ab6040) (Abcam, Cambridge, MA, USA). After washing the PVDF membrane several times, the PVDF membrane was incubated with a corresponding secondary antibody (1:4000, ab6721). After washing the PVDF membrane several times again, bands on the PVDF membrane were detected on a Bio-Rad ChemiDoc XRS system.

Total RNA was extracted and reverse-transcribed into complementary DNA (cDNA) using the all-in-one first-strand cDNA synthesis super mix kit (TransGen Biotech, Beijing, China). qRT-PCR analysis was performed using the fast green qRT-PCR supermix (TransGen Biotech, Beijing, China). All samples were divided into tumor tissue (T) or non-tumorous tissue (N) to detect the mRNA expression of AFF3. The 2${}^{-△△\text{CT}}$ method was used to calculated relative gene expression. The corresponding primer sequences were as follows: AFF3: 5${}^{\prime }$-ACTCAACAGGATGATGGC-3${}^{\prime }$ (forward) and 5${}^{\prime }$-TGCCTAAAGTGTTCTGGATC-3${}^{\prime }$ (reverse); glyceraldehyde-3-phosphate dehydrogenase (GADPH): 5${}^{\prime }$-GGTGAAGGTCGGAGTCAACG-3${}^{\prime }$ (forward) and 5${}^{\prime }$-CAAAGTTGTCATGGATGHACC-3${}^{\prime }$ (reverse).

All statistical analyses were performed using the R software (version 3.5.0, the Vienna University of Economics and Business, Vienna, Austria). Clinicopathological and immune infiltration data between different groups were tested by Wilcoxon test. The Correlation between AFF3 expression and TICs was examined using Spearman correlation analysis. All statistical tests were two-sided, and the level of significance was set at p 0.05.

As shown in Table 1, we selected three GEO datasets for analysis. The number of upregulated and downregulated genes in different datasets were as follows: GSE38959, 939 upregulated genes and 509 downregulated genes (Fig. 1A); GSE62931, 356 upregulated genes and 483 downregulated genes (Fig. 1B); GSE64790, 509 upregulated genes and 630 downregulated genes (Fig. 1C). Fig. 1D shows a Venn diagram with 75 overlapping DEGs, comprising 38 upregulated and 37 downregulated genes (Supplementary Table 1).

Fig. 1.

Identification of differentially expressed genes (DEGs) in triple-negative breast cancer (TNBC). (A–C) Volcano map of three genes expression profiles in GEO datasets, GSE38959 (A), GSE62931 (B), and GSE64790 (C). Red and blue colors indicate upregulated and downregulated genes in tumor tissues, respectively. (D) Venn diagram showing the common DEGs in three datasets.

GO and KEGG functional enrichment analyses were performed with ClusterProfiler to examine the functional implications of 75 DEGs between TNBC and non-TNBC samples. The biological process (BP) included mitotic nuclear division, nuclear division, organelle fission, chromosome segregation and kinetochore organization. Cellular components (CC) were condensed chromosome kinetochore, condensed chromosome, condensed nuclear chromosome kinetochore, condensed chromosome, centromeric region, kinetochore; KEGG enrichment suggested progesterone-mediated oocyte maturation, cell cycle, and oocyte meiosis (Fig. 2, Supplementary Table 2).

Fig. 2.

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses of DEGs between TNBC and non-TNBC tissues. (A) Enriched GO terms in the “biological process” category. (B) Enriched GO terms in the “Cellular Components” category. (C) KEGG pathway annotations. The size of the circle represents the number of genes contained in a category. The higher the number of genes, the larger the circle. The color represents the adjusted p-value in enrichment analysis. The smaller the p-value, the higher the red color intensity. The larger the p-value, the higher the blue color intensity.

To screen for drugs for TNBC treatment, upregulated and downregulated DEGs were separately uploaded into the L1000FWD database. Small-molecules with anticancer effects on TNBC progression were selected with a similarity score of zero and an adjusted p-value of 0.01. The five most similar small-molecule drugs were idarubicin, teniposide, homosalate, palbociclib, and tremulacin (Table 2).

The correlation between gene expression and clinical characteristics in TCGA-TNBC dataset was examined. Gene expression matrices and clinicopathological data of 152 TNBC patients were obtained. Univariate Cox regression analysis revealed that 13 DEGs (NUF2, TPX2, BUB1, EZH2, ASPM, EXO1, DEPDC1, AFF3, H1-1, PSAT1, NDC80, CENPF, and KIF18B), race, T stage, and N stage were correlated with prognosis of TNBC patients (p 0.1 was considered statistically significant, Supplementary Table 3). Multivariate analysis revealed that among the 13 genes and 3 clinicopathological variables mentioned above, only AFF3 downregulation (hazard ratio (HR): 0.22; 95% confidence interval (CI): 0.07–0.72; p = 0.012) was an independent risk factor for OS in patients with TNBC (Table 3).

The transcription levels of AFF3 in TCGA-TNBC cohort were analyzed. The expression of AFF3 in TNBC tumor tissues was significantly lower than that in healthy mammary tissues (p = 4.40 $\times$ 10${}^{-27}$, Fig. 3A). Compared with those in the paired non-tumorous tissues, the AFF3 mRNA levels were downregulated in the TNBC tumor tissues. Next, the mRNA and protein expression levels of AFF3 in 25 paired samples of TNBC collected from patients at Jiangxi Cancer Hospital were examined using qRT-PCR (Fig. 3C,D) and western blotting (Fig. 3E). Compared with those in the non-tumorous tissues (N), the mRNA (with an average fold change of 2.68, Fig. 3C) and protein levels (Fig. 3F) of AFF3 were downregulated in the TNBC tumor tissues (T).

Fig. 3.

Expression of AFF3 in TNBC and non-tumorous tissues. (A) The relative mRNA expression of AFF3 in 152 TNBC tissues and 99 non-tumorous tissues in TCGA dataset. Significant differences between the two groups were evaluated using the Wilcoxon rank sum test. (B) The relative mRNA expression of AFF3 in 10 TNBC tissues and paired-adjacent tissues from TCGA dataset. Significant differences between the two groups were analyzed using the Wilcoxon signed-rank test. (C) The AFF3 mRNA levels in 25 paired samples of TNBC are represented as a histogram. N/T expression value $\ge$2 suggests a significantly lower expression, which $\le$1/2 is a significantly higher expression, and values between 1/2 and 2 show no significant change. N, non-tumorous tissue; T, tumor tissue. (D) The AFF3 mRNA expression levels in paired N and T from 25 TNBC patients. Significant differences between the two groups were evaluated using the Wilcoxon signed-rank test. (E) The AFF3 protein levels in 25 paired N and T. $\beta$-tubulin was an internal control. (F) The quantified AFF3 protein level in 25 paired N and T. Significant differences between the two groups were evaluated using the Wilcoxon signed-rank test.

Next, the correlation between clinicopathological and AFF3 expression was further investigated in TCGA-TNBC cohort. All patients in the TCGA-TNBC cohort were categorized into AFF3-high and AFF3-low groups, based on the median AFF3 expression levels. The AFF3 mRNA expression level was significantly correlated with T stage (p = 0.025, Fig. 4A), N stage (p = 0.008, Fig. 4B), pathologic stage (p = 0.046, Fig. 4D), race (p = 0.023, Fig. 4E). In TCGA-TNBC cohort, the median OS in the AFF3-high group was higher than that in the AFF3-low group (p = 0.023, Fig. 4F).

Fig. 4.

Correlation of AFF3 mRNA expression with clinicopathological characteristics. (A–C) The correlation of AFF3 mRNA expression with T stage (A), N stage (B), and M stage (C). Significant differences between the two groups examined using the Wilcoxon Mann-Whitney test. (D,E) The correlation of AFF3 mRNA expression with pathologic stage (D) and race (E). Significant differences between multiple groups were examined using the Kruskal-Wallis rank sum test. (F) Survival analysis of patients with TNBC in the AFF3-high and AFF3-low groups in TCGA-TNBC cohort. Significant differences between multiple groups were examined using the log-rank test. The number of samples as follows: T1 + T2: 133, T3 + T4: 19; N0 + N1: 131, N2 + N3: 15; M0: 129, M1: 2; stage I: 27, stage II: 96, stage III: 24, stage IV: 2; Asian: 8, White: 84, African: 53.

To further investigate the function of AFF3, the DEGs between AFF3-high and AFF3-low groups were examined using the Limma package. The criteria for selecting DEGs were as follows: FDR 0.05 and $|$log2 FC$|$ $\ge$2. In total, 182 up-regulated genes, and 67 downregulated genes were obtained. The correlation between the expression of DEGs and samples was shown using a heatmap (Fig. 5A, Supplementary Table 4). GO and KEGG enrichment analyses revealed that DEGs were significantly correlated with some biological processes and signaling pathways, such as cornification, cornified envelopes, and oxidoreductase activity (Fig. 5B,C, Supplementary Table 5). Cytochrome P450 enriched steroid hormone biosynthesis, drug metabolism, and xenobiotic metabolism were annotated in KEGG enrichment analysis (Fig. 5D, Supplementary Table 5).

Fig. 5.

GO/KEGG enrichment analysis for DEGs in the AFF3-high and AFF3-low groups. (A) Heatmap of DEGs between AFF3-high and AFF3-low groups. The row of the heatmap represents the gene symbol of the 249 DEGs and the column represents ID number of the samples in TCGA-TNBC cohort. Gene symbols and identification numbers are not shown in the graph. (B) Enriched GO terms in the “biological process” category. (C) Enriched GO terms in the “cellular component” and “molecular function” category. (D) KEGG pathway annotations. Only some notable and leading gene sets are displayed in the graph.

GSEA was performed to investigate the differential activation of AFF3-related signaling pathways in TNBC. GSEA identified some critical signaling pathways in MSigDB collection (c2.cp.v7.2.symbols.gmt and c7.all.v7.2.symbols.gmt). Based on the normalized enrichment score (NES) and p-values, some significant signaling pathways are shown (Fig. 6 and Supplementary Table 6). For the C2 (one of the human collections of the molecular signatures database) collection, neuroactive ligand receptor interaction, steroid hormone biosynthesis and phase II conjugation of compounds (Fig. 6A, Supplementary Table 6) were enriched in AFF3-high group. Meanwhile, the formation of the cornified envelope, keratinization, translation, respiratory electron transport, and complex I biogenesis were enriched in AFF3-low group (Fig. 6B, Supplementary Table 6). For the C7 collection, some different immune-related signaling pathways were enriched in AFF3-high group (Fig. 6C, Supplementary Table 7) or AFF3-low group (Fig. 6D, Supplementary Table 7). These results further suggested that AFF3 had an important role in the immune microenvironment of TNBC.

Fig. 6.

GSEA identifies AFF3-related signaling pathways. (A,B) Enriched signaling pathways in the C2 collection in the AFF3-high (A) or in AFF3-low expression groups (B). (C,D) Enriched signaling pathways in the C7 collection in AFF3-high (C) and AFF3-low groups (D). GSEA, Gene Set Enrichment Analysis.

GSEA revealed that AFF3 may be involved in the immune function in TNBC. Next, the distribution of immune cells in TNBC was analyzed with CIBERSORT algorithm. M0, M1, and M2 macrophage subsets accounted for most of the infiltrating immune cells (Fig. 7A). Violin plots were used to compare immunity subset distributions between the AFF3-low and AFF3-high groups (Fig. 7B). The proportions of naive B cells, neutrophils, and M2 (a subtype of macrophage) macrophages were significantly different between the AFF3-low and AFF3-high groups (Fig. 7B). The proportions of naive B cells, mast cells resting and T cells CD4 (Cluster of differentiation 4) memory resting were positively correlated with AFF3 expression (Fig. 7C–E). In contrast, the proportions of neutrophils and follicular helper T cells were negatively correlated with AFF3 expression. Therefore, the level of AFF3 regulates the proportion of TICs and immune activity in TME (Tumor microenvironment).

Fig. 7.

Correlation of tumor-infiltrating immune cell (TIC) proportion with AFF3 expression. (A) Bar plot shows the distribution of TICs in each tumor sample in TCGA-TNBC cohort. Column are the ID numbers of the samples, which are not shown in the plot. (B) The violin plot shows the distribution of TICs between the AFF3-high and AFF3-low groups in TCGA-TNBC cohort. (C–G) The scatter plot shows the distribution of the 5 TICs that were significantly correlated with AFF3 expression.