Level 3 RNA sequencing profiles of BC patients and accompanying clinical characteristic data were obtained from TCGA database. We excluded the samples which had absent or deficient information on age, TNM stage, overall survival (OS) time, or survival status. The inclusion criteria were as follows: 1) local invasive BC identified by using the pathological method; 2) OS time was longer than 1 month; and 3) both survival data and RNA expression profiles were completely available. Clinical features contained TNM stage, T stage, N stage, age, estrogen receptor (ER), progesterone receptor (PR), HER-2 status, tumor subtypes, OS, and OS time. Thus, 999 BC patients were divided into two groups (the high- and the low-expression groups) depending on the median value of YTHDF3 expression. The correlation between YTHDF3 expression and the clinicopathological factors was determined. Meanwhile, we analyzed the YTHDF3 differential expression in 27 diverse kinds of cancers.

The Tumor Immune Estimation Resource (TIMER) database provided a comprehensive analysis of immune infiltrates of various types of tumors (https://cistrome.shinyapps.io/timer). We performed the YTHDF3 expression in BC with the abundance of tumor-infiltrating immune cells (TIICs) through gene modules in TIMER. Furthermore, by using CIBERSORTx (http://cibersortx. Stanford. edu/), which is a deconvolution algorithm dependent on gene expression, we analyzed the immune microenvironment to explore the correlation between TIICs and YTHDF3 expression in BC patients.

In our report, 999 BC cases were divided into two groups (the high- and low-expression groups) according to the median value of YTHDF3 expression. We performed GSEA to determine the biological process that is enriched by YTHDF3 in BC patients. We selected “c2.cp.kegg.v6.2.symbols.gmt” as the gene set. The permutation number was set at 1,000 times, the minimum and maximum number sizes for gene analysis at 15 and 1,000, respectively. Gene sets with a nominal p-value < 0.05 or a false discovery rate (FDR) q < 0.25 were the significant threshold.

According to the analysis of univariate and multivariate CPHR, the YTHDF3-based model was conducted. The performance of the nomogram was quantified in reference to the ability of discrimination, calibration, and clinical utility. The score test was applied to estimate the goodness of fit of the model. The AUC of ROC curves was utilized to figure out the discrimination ability. The calibration ability was assessed via a calibration curve. In addition, DCA was used to illustrate the clinical utility of the model.

We employed “pRRophetic” R package to estimate the chemotherapy and the targeted therapy response combined with the Genomics of Drug Sensitivity in Cancer (GDSC, https://www.cancerrxgene.org) (Geeleher et al., 2014). Five commonly used chemotherapeutic agents vinorelbine (Adamo et al., 2019; Falvo et al., 2021; Wang et al., 2021), cisplatin (Hu et al., 2015; Jovanović et al., 2017; Tung et al., 2020), methotrexate (Cameron et al., 2017; Lu et al., 2020; Shakeran et al., 2021), doxorubicin (Mackey et al., 2016; Ganz et al., 2017; Schneeweiss et al., 2022), paclitaxel(Diéras et al., 2020; Fernandez-Martinez et al., 2020; Tolaney et al., 2021), and small-molecule inhibitor targeting EGFR lapatinib (Cortés et al., 2015; Powles et al., 2018; Johnston et al., 2021) were selected and kept the default values for all parameters. The response prediction relied on IC50 (the half-maximal inhibitory concentration) of each sample. The method for estimating the value of IC50 is ridge regression analysis. Through 10× cross-validation, the prediction accuracy was assessed.

The Wilcoxon rank-sum test was applied to compare the expression of YTHDF3 in BC and normal groups. We used the Kaplan–Meier (K-M) method and log-rank test to analyze the OS of two groups (the high- and low-YTHDF3 expression groups), respectively. CPHR analysis was conducted to screen the independent factors of OS. Ridge regression analysis was performed to estimate the value of IC 50. The association between YTHDF3 and TILs was determined via the TIMER and CIBERSORTx algorithm. R software (version 4.0.3) was applied to determine the correlation between the expression of the RNA-seq gene and clinical parameters in BC. A p value <0.05 meant a significant difference.

The clinical characteristics of 999 BC patients, including age, TNM stage, molecular subtypes, ER, PR, HER-2, and vital status, were extracted from TCGA database. In light of the inclusion and exclusion criteria, 999 patients were involved in later analysis. In our research, the patients’ age when diagnosed ranged from 26–89 years, and the medium of age was 57–91 years. The median follow-up time was 42.4 months, with median ranged from 1 to 286.8 months. The result was shown in Table 1.

In TCGA database, the YTHDF3 differential expression was analyzed between tumor and normal tissues in 27 diverse cancer types (Figure 1). The YTHDF3 expression was significantly highly upregulated in 1,098 BC patients versus 292 normal group individuals (p < 0.001), as well as in 21 of the 27 analyzed cancer types (Figure 1).

We conducted the CPHR analysis to figure out the correlation between OS and the YTHDF3 expression, as well as other clinical characteristics in BC patients. As Table 2 shows, univariate analysis suggested that several factors, including YTHDF3 (HR = 1.02, p = 0.003), age (HR = 1.032, p < 0.001), and TNM stage III are significantly associated with OS. Multivariate CPHR uncovered that YTHDF3 expression was an independent prognostic factor for OS. The expression distribution and profiles of YTHDF3 and survival status of BC patients are demonstrated in Figure 2A.

In order to determine the prognostic role of YTHDF3 in BC, the expression data were separated into high- and low-risk groups in light of the median value of YTHDF3 expression. The K-M method was conducted to describe the OS of two groups. The log-rank test of OS curves showed that elevated YTHDF3 expression was combined with shorter OS time in BC patients (p = 0.013).

The tumor-infiltrating lymphocytes (TILs) can independently predict the survival status (Azimi et al., 2012). Furthermore, TILs play an essential role in facilitating a favorable fate in BC patients (Stanton and Disis, 2016). Hence, we explored the association between immune infiltration and YTHDF3 expression in BC via TIMER. As shown in Figure 3A, YTHDF3 expression positively related to the levels of B cells (p = 5.02 × 10^–3), CD4⁺ T cells (p = 1.56 × 10^–3), CD8⁺ T cells (p = 1.11 × 10^–31), macrophages (p = 8.96 × 10^–24), neutrophils (p = 9.98 × 10^–13), and dendritic cells (p = 1.16 × 10^–7). These results revealed that YTHDF3 was a crucial part in the tumor immune infiltration of BC.

Meanwhile, we attempted to determine if the immune infiltrating microenvironment varied according to the high and low YTHDF3 expression levels. The CIBERSORTx algorithm (https://cibersortx.stanford.edu/) was utilized to identify different infiltrating levels of diverse immune cell subtypes in the high- and low-YTHDF3 expression groups (Figure 3B). Among the 22 subpopulations of immune cells, plasma cells, CD8⁺ T cells, memory B cells, resting memory CD4⁺ T cells, regulatory T cells, activated memory CD4⁺ T cells, activated NK cells, resting mast cells, and M2 macrophages were considerably influenced by YTHDF3 expression. Specifically, memory B cells (p = 0.019), plasma cells (p = 0.017), activated NK cells (p < 0.001), regulatory T cells (p < 0.001), and M2 macrophages (p = 0.004) apparently increased in low-YTHDF3 expression group. In contrast, CD8⁺ T cells (p < 0.001), resting memory CD4⁺ T cells (p < 0.001), resting mast cells (p = 0.003), and activated memory CD4⁺ T cells (p = 0.008) were significantly elevated in the high-YTHDF3 expression group.

The possible link among 22 subtypes of immune cells was also evaluated (Figure 3C). The heatmap demonstrated that the ratios of diverse TILs ranged from weakly to moderately correlated.

GSEA was performed to select the YTHDF3-related signaling pathway in BC. In light of the median YTHDF3 expression value, TCGA data were divided into the high- and low-YTHDF3 expression groups. A nominal p-value <0.05 and an FDR <0.25 were regarded as significant. The results of GSEA demonstrated that five KEGG items (selected based on the NES) were significantly differentially enriched in the high-YTHDF3 expression group, including ubiquitin-mediated proteolysis, inositol phosphate metabolism, epithelial cell signaling in Helicobacter pylori infection, and glycosylphosphatidylinositol GPI anchor biosynthesis, whereas the gene set involved in the low-YTHDF3 expression group was glycosaminoglycan biosynthesis chondroitin sulfate and ribosome. All snapshots of enrichment results are shown in Figure 4.

Age, TNM stage, BC subtype, and the YTHDF3 expression value were determined as independent variables through the multivariate CPHR analysis (Table 2). The YTHDF3-based nomogram is shown in Figure 5A. The result of the score test was 77.43, and the p value was 5e-13, which indicated good model fitness. The AUC of the model was 0.721 (95% CI: 0.654–0.789), which suggested good recognition ability of the model and prediction superiority over the TNM staging system (AUC: 0.72 vs. 0.3, p < 0.00028) (Figure 6B). The calibration curve showed high consistency between the predicted 5-year OS rate by model and the actual observed results (Figure 6A). Furthermore, with respect to clinical utility, the YTHDF3-based model was better than the TNM staging system on account of the wider range of threshold probability via DCA (Figure 6C).

Six commonly used BC therapeutic drugs (doxorubicin, vinorelbine, paclitaxel, methotrexate, cisplatin, and lapatinib) were selected as the candidates, which might show drug resistance in the high-risk group via the estimation of the IC₅₀ values. As displayed in Figure 7, the estimated IC50 values of doxorubicin, vinorelbine, paclitaxel, and methotrexate significantly increased in the high-risk group. There were no estimated IC50 value differences of cisplatin and lapatinib observed between the high- and low-risk groups.