Identification of Alternative Splicing Events and Mutated Genes in Breast Cancer Patients
To search for potential antigens in breast cancer, we analyzed alternative splicing events in TCGA patients and obtained 45421 alternative splicing events. Of these, cassette exon (ES) (17,702) had the most alternative clipping events, while matually exclusive exons (ME) (233) had the fewest (Figs. 1A-B). In addition, most of the time, ME events coincided with other alternative splicing events; Seven types of alternative splicing events occurred in genes COA1 and PTK2 (Fig. 1C). Next, we reduced the abnormal alternative splicing events to 2105 events using the thresh hold of P-value < 0.05 plus |logFC| > 1. As a consequence, 1267 alternative splicing events (988 genes) were upregulated, and 838 alternative splicing events (692 genes) were downregulated in breast cancer tissues versus normal adjacent tissues (Figs. 1D-E). Next, we analyzed the mutated genes in breast cancer samples and found 6227 frameshift mutations. The results revealed that the gene with the highest mutation rate was PIK3CA (Fig. 2A). Meanwhile, we also analyzed the altered genome fraction and mutation counts (Figs. 2B-C), and PIK3CA was still the most frequently altered gene in altered genome fraction and mutation counts (Figs. 2D-E).
Identification of Potential Antigens in Breast Cancer Patients
Recent studies have shown that transcripts containing transcoding mutations and abnormal splicing patterns may produce antigenic peptides in the presence of NMD disruption. In this case, abnormal transcripts were equally important as potential biomarkers for immune checkpoint inhibitors (ICI). Next, we took the intersection of genes with upregulated abnormal alternative splicing (988 genes) and genes with frameshift mutations (6227 genes) and agreed on 375 overlapped genes (potential antigen candidate genes) (Fig. 3A). These genes were found to be involved in BPs such as cell morphogenetic regulation, axonal formation, cellular junction formation, and epithelial cell proliferation, and their proteins were mostly found in adhesion spots, cell-matrix junctions, cell fronts, and intercell junctions. Actin binding, DNA-binding transcription factor binding, GTPase activity regulation, nucleoside triphosphatase activity regulation, and cadherin binding were relevant to the genes' MFs. Also, actin cytoskeleton modulation, axonal orientation, vascular smooth muscle contraction, and dilated cardiomyopathy were all shown to be overrepresented in the KEGG pathway study (Figs. 3B-E). 370 of the 375 overlapping genes were matched to RNA-Seq data from the TCGA breast cancer study, and 63 genes were differently expressed in tumor tissues versus normal tissues, including 10 upregulated and 53 downregulated genes (Figs. 3F-G and Supplementary Table S1).
The Confirmation of Five Potential Antigens for mRNA Aaccines targeting Breast Cancer
The action of NMD protects cells from producing probable abnormal transcripts, which would otherwise lead to the occurrence of tumors[14]. Since the NMD factor was linked to the degradation of abnormal splicing variants, we investigated the correlation between the expression of 370 intersecting genes and the NMD factors and presented the top 20 intersecting genes with the most significant differences. Notably, the expression of most of intersecting genes and NMD factors were identical (Supplementary Figures S1A-D). Finally, we analyzed the PSI values of 2545 alternative splicing events in the groups with different expression of NMD factors. As a result, the PSI values of the first 20 alternative splicing events were higher in the group with high expression of NMD factors (Supplementary Figures S1E-H).
Based on univariate Cox regression analysis, we discovered that 29 potential antigen candidate genes were associated with the OS of breast cancer patients, among which nine genes were protective factors and the other 20 genes were risk factors. The first six genes most correlated with the OS were SERPINA1, DLG3, CDHR3, PLP2, GREB1L, and IRF2 (Figs. 4A-C). After intersecting 29 survival-related genes and 370 possible antigen candidate genes, five genes including DST, AN06, LAMA3, FGD3, and NEDD9 were identified as putative antigens. Following that, using data from the TIMER database, the relationship between the expression of possible antigen genes and immune cell infiltration was investigated. The expressions of DST, ANO6, and NEDD9 were found to be adversely linked with tumor purity. LAMA3 expression was adversely linked with B cell infiltration. GENES DST, AN06, LAMA3, and NEDD expression were favorably linked with the infiltration of macrophage cells. Similarly, DST, AN06, and NEDD9 expression were substantially associated with dendritic cell infiltration (Figs. 4D-F, H). In addition to FGD3 (Fig. 4G), DST, AN06, LAMA3, and NEDD9 expression were relevant to the infiltration of antigen-presenting cells (APCs). Thus, the transcripts produced by these four genes might be directly processed by the APCs, delivered to T cells and recognized by B cells, then triggering an immune response. Consequently, there is hope for DST, AN06, LAMA3, and NEDD9 as mRNA vaccines for breast cancer.
Potential Prognostic Biomarkers that Predict Vaccine Response
Next, we hoped to find biomarkers that predict vaccine response. TCGA Breast cancer patients were clustered using WGCNA analysis based on the immune gene co-expression module, with a scale-free network soft threshold of 2. Finally, five modules were obtained, among which the gray modules did not aggregate together (Figs. 5A-B). Subsequently, the prognostic analysis was performed on the four modules (Fig. 5C). The results revealed that the Turquoise module had significant prognostic efficacy. To explore the genetic functions of this module, we performed KEGG enrichment analysis on genes in the Turquoise module and found that the module genes interacted with T cell receptors signaling pathway, cytokines receptors, cytokine signaling pathways, natural killer cell-mediated cytotoxicity effect, hematopoietic cell lineage, leishmania infection infection, produce IgA intestinal immune network, allogeneic graft rejection, graft versus host disease, and type 1 diabetes (Fig. 5D). As is well-known, the therapeutic potential of mRNA vaccines in cancer patients is largely determined by the infiltration and activation of T cells, as well as the immunosuppressive cells[8, 25, 26]. As T cell receptors signaling pathway was associated with Turquoise module, the genes in Turquoise module might be used as biomarkers to predict vaccination outcome. Thus, we selected hub genes in Turquoise module (correlation > 90%), and performed univariate Cox regression on the hub genes (n = 36) and got five genes as putative biomarkers, namely TRBC2, CD3D, CD27, CD3E, and TRBV28.
Furthermore, we also analyzed the potency of these five genes as prognostic genes through constructing a risk model and scored the samples. The samples with higher score than the median value were defined as high-risk patients, and those with lower score than the median value were defined as low-risk patients (Figs. 5E-F). We analyzed the survival of patients in the high and low-risk groups and tested the prognostic efficacy of the risk model using RiskScore. The results showed that the risk model accurately distinguished high- and low-risk patients (Fig. 5G). Finally, we analyzed the expression of five biomarker genes in the low-risk group and observed that all biomarker genes were higher expressed in the low-risk group (Fig. 5H).
Identification of Potential Immune Subtypes of BRCA
Immunotyping can be used to reflect the immunological state of tumors, assisting in identifying patients who are vaccine beneficiaries. As a result, we created an induced cluster by analyzing the expression profiles of immuno-prognostic genes in breast cancer tissues. In Figs. 6A-C, the stable cluster of immune-related genes was chosen based on the cumulative distribution and incremental region, and three immunological subtypes of breast cancer were defined, dubbed Cluster 1, Cluster 2, and Cluster 3, respectively. The three subtypes were then confirmed using principal component (PC) analysis, with the findings demonstrating that the three subtypes could be well distinguished (Fig. 6D). Immune gene expression was also revealed in Cluster 1, Cluster 2, and Cluster 3. Finally, we looked at the OS of Cluster 1, Cluster 2, and Cluster 3 patients and discovered that their prognosis differed significantly (Fig. 6E). Subsequently, we explored how patients from the three clusters were distributed across different breast cancer categories in the training cohort. The proportions of the three categories remained steady in patients with varied T stages, according to the findings (Fig. 6F). When we looked at triple-negative breast cancer (TNBC) on its own, we noticed that the proportion of Cluster 2 patients in stage 2 and stage 3 was substantially higher (Fig. 6G). This suggests that there may be a relationship between staging and immunological type in TNBC patients. Similarly, we used the validation data set to confirm the previous findings, and found the correlation results were consistent with the training set results (Figs. 6H-J). Despite the fact that the distribution pattern of Cluster 2 TNBC patients differed from that of the training data set, the distribution trend of patients from three clusters in various T stages matched that of the training test.
Cellular and Molecular Characteristics of the Immune Subtypes
The immunological condition of the patient influences the patient's reaction to an mRNA vaccine. Hence, we calculated the differences in immune score, stroma score, tumor purity, and cytolysis activity (CYT) among the three subtypes, and the results demonstrated that the four scores were considerably different among the three subtypes (Fig. 7A). Then we compared the differences of immune cell infiltration in the three clusters of patients and observed that the infiltration ratios of 18 types of the immune cells were statistically different between groups (Fig. 7B). Additionally, we double-checked the immune score, stroma score, and tumor purity (Fig. 7C) in the validation set (GSE25066) and also found clear distinctions among Cluster 1 versus Cluster 3, and Cluster 2 versus Cluster 3. Moreover, the infiltration ratios of 17 types of immunoinfiltrating cells were significantly different among the three groups (Fig. 7D).
After that, we looked at the link between this study's immunotyping and existing generalized immunotyping. Except for the identical proportion of C6 in the three immune subtypes, the proportion of C2 in Cluster 2 was significantly higher than that of Cluster 1 and Cluster 3, inferring that Cluster 2 may have a more significant relationship with known C2 clusters than Cluster 1 and Cluster 3 (Fig. 7E). C2 patients possesses IFN-γ-dominated feature, the highest lymphocytic infiltration, CD8 T-cell-related features, and the highest M1 macrophage content, which reflects the stronger antitumor immune response ability. Cluster 2 patients had a better prognosis in this study due to a larger immune cell infiltration. Furthermore, C4 only intersected with Cluster 1. C4 was an immune-suppressed microenvironment with macrophage-dominated, low lymphocytic infiltration, and high M2 macrophage content, which was consistent with the immunological and prognostic characteristics of Cluster 1 that we discovered in this investigation. The results show that our immunotyping is accurate and useful in addition to the existing generalized cancer subtypes. We also looked at the immunological activity scores of patients from various clusters. Significant differences in 24 types of immunological activity scores, such as T cells, CD4 T cells, and CD8 T cells, were discovered (Fig. 7F). Cluster 2 received the highest comprehensive score, indicating a highly active immune microenvironment and hence potential "hot tumors" that were not suitable for mRNA vaccine immunization. Cluster 1 and Cluster 3, on the other hand, had lower immune activity scores and were putative "cool tumors," making them more suitable for vaccination.
The Immune Landscape of BRCA
The mapping of immune landscape may be able to identify more refined subgroups within each cluster and thus be helpful in identifying vaccine-sensitive patients[22]. We constructed the immunomap of breast cancer using the expression profile of immune genes in each single patient (Fig. 8A). As illustrated in Fig. 8B, B cel1s, macrophages and T cells were substantially linked with PC1 and PC2. The point distribution of Cluster 1 and Cluster 3 was relatively discrete. Therefore, Cluster 1 and Cluster 3 were further divided into two subgroups (Cluster 1A, Cluster 1B, Cluster 3A, and Cluster 3B) according to the location of the immune cell population (Fig. 8C), and significant differences of immune cell enrichment scores were detected among further subdivided subgroups (Figs. 8F-G). In addition, we compared the outcomes of samples with extreme distribution in the immune landscape and found that patients in group 2 had higher OS (Figs. 8D-E). According to the above results, an immune landscape based on immune subtypes may identify the immune components of breast cancer patients and then predict responses to vaccines.