The RNA-sequencing profile for 375 primary GC, 32 normal gastric tissues, and their corresponding clinicopathological characteristics (age, grade, stage, and sex) and survival time were download from The Cancer Genome Atlas-Stomach Adenocarcinoma (TCGA-STAD) database (https://portal.gdc.cancer.gov/). Similarly, the profile of 100 normal gastric tissues and 300 GC tumor were download from GEO database (GSE66229).

Anti-PD-1 mAb (Nivolumab) was obtained from Bristol-Myers Squibb Company, New York, New York, USA. FGF 401 (Selleck Chemicals, Houston, Texas, USA), FGFR1(D8E4, #9740), FGFR4(D3B12, #8562), Erk1/2(137F5, #4695), Phospho- Erk1/2 (Thr202/Tyr204) (D13.14.4E, #4370), Akt (C67E7, #4691), Phospho-Akt (Ser473) (D9E, #4060), E-Cadherin (4A2, #14472), MMP-9 (D6O3H, #13667) antibody (Cell Signaling Technology, Danvers, Massachusetts, USA), anti-fibronectin (ab32419), anti-MMP2 (ab97779) (Abcam, Cambridge, UK), and β-actin (Sigma-Aldrich, St Louis, Missouri, USA) were used as obtained.

Twenty-seven pairs of GC tissues and their adjacent normal tissues were obtained from General Surgery Department of the First Affiliated Hospital at Xi’an Jiaotong University, China, in 2020 and 2021. Immediately after resection, GC and normal adjacent tissues specimens were stored at -80°C. This study was approved by the Ethics Committee of the First Affiliated Hospital of Xi’an Jiaotong University, and conducted in accordance with the Declaration of Helsinki principles. Informed consent was obtained from all patients.

All human GC cell lines were obtained from the Cell Resource Center, Peking Union Medical College (Peking, China) and cultured in Roswell Park Memorial Institute Medium (RPMI) or Dulbecco’s Modified Eagle Medium (DMEM) with 10% Fetal Bovine Serum (Invitrogen) at 37°C and 5% CO₂. The 293T cells were obtained (Clontech Laboratories, San Francisco, California, USA) and cultured in DMEM with high glucose. pHAGE-FGFR1 and pHAGE-FGFR4 were gifts from Gordon Mills & Kenneth Scott, Addgene plasmid #116740, #116743. Retroviral infections were performed as previously described (19, 20). Transfected cells were established under puromycin selection, until stable cell lines were generated as N87-EV, N87-FGFR1, and N87-FGFR4.

Gene Expression Profiling Interactive Analysis database (GEPIA, gepia.cancer-pku.cn) provided vast integrated genome expression profiles in GC and normal samples (21). The mRNA expression levels of FGFR in GC tissues were compared to normal tissues. We further examined the expression of FGFRs in multiple tumor samples through ONCOMINE (www.oncomine.org) and Tumor Immune Estimation Resource (TIMER) databases (https://cistrome.shinyapps.io/timer/) to validate the former result. In addition, we used the Cancer Genome Atlas (TCGA) and Genotype-Tissue Expression (GTEX) data to evaluate the expression of FGFRs members in tumor and adjacent normal tissues to explore the expression in GC patients of different clinicopathological characteristics. The data was analyzed through the “ggplot2” R package using R software from Bell Laboratories (formerly AT&T, now Lucent Technologies, New Jersey, USA) by John Chambers and colleagues.

Kaplan–Meier plot (www.kmplot.com) was utilized to examine the prognosis of FGFRs expression in STAD. The survival outcomes were OS and progress-free survival (PFS) (22). In addition, the association were evaluated between mRNA levels of FGFRs and clinical characteristics. The univariate and multivariate analysis were applied using the Cox regression model. A p-value<0.05 was considered significant.

The LinkedOmics database was used to analyze the correlated genes with FGFRs through the “limma” R package to extract the expression profiles of the correlated gene in TCGA-STAD data sets. Then, we drew the correlation heatmap of the top 10 genes through TCGA based on FGFRs different expression groups (Pearson Correlation Coefficient, p<0.05). The Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were performed using the R package “clusterProfiler” to assess the relative biological functions and molecular pathways regulated by the correlated genes. A p-value<0.05 was set as significant for screen criteria.

The cBioPortal database (http://cbioportal.org) was used to explore the FGFRs gene alterations in GC (23). Then the TIMER somatic copy number alternation (SCNA) module compares tumor infiltration levels in tumors with different SCNA for FGFR family members through a two-sided Wilcoxon rank-sum test.

TIMER database is a publicly available source for evaluating the impact of different immune cells in diverse cancers (24). The correlation between FGFRs and the 24 different kinds of immune cell infiltrates in gastric cancer samples, including CD8+ T-cells, CD4+ T-cells, neutrophils, macrophages, and natural killer cells, was achieved by the “ssGSEA” method through the GSVA package with R computing. TIMER was utilized to assess the correlation between FGFRs expression and tumor purified immune cell infiltration.

Functional experiments (cell proliferation, transwell invasion assay, and scratch assay), RT-qPCR, and immunoblotting were performed as described in our previous articles (25–26). FGFR1, FGFR2, FGFR3, and FGFR4 primers (OriGene, Rockville, Maryland, USA) were used to detect mRNA levels by real-time quantitative reverse transcription polymerase chain reaction (RT- PCR). All the experiments were repeated three times independently.

Immunodeficient NOD Cg-Prkdc^scidIl2rg^tm1Sug/JicCrl (NOG) mice (Weitong Lihua Experimental Animal Co., Ltd, Beijing, China) were used to receive human immune cells and established a humanized immune system. Peripheral blood mononuclear cells (PBMCs) from healthy donors were isolated and 2×10⁷ cells were administered intraperitoneally into the mice. BALB/c nude mice were purchased from the Laboratory Animal Center of Xi’an Jiaotong University. N87-EV, N87-FGFR1, and N87-FGFR4 cells (5 ×10⁶) were administered subcutaneously into mice. After 10 days, the subcutaneous tumors had grown to a size that could be measured (approximately 90 mm³). Then FGF401was diluted in 6% 1M HCl/100 mM citrate buffer pH 2.5 (citrate salt) and administered intragastrically at 30 mg/kg daily. Mouse was injected with Anti-PD-1 mAb (nivolumab) at 200 mg/mouse twice a week for 4 weeks. Every 2 days, the tumor volume and body weight were measured and calculated as (length × width²)/2. Mice were euthanized when Tumor volume (TV) exceeded 1,500 mm³ or when the mouse’s weight decreased by < 70% of that at Day 1. Survival (in days) was defined as the time for each mouse from Day 0 until the day of death or euthanization. Survival curve was built for each group through GraphPad Prism 7.0 (GraphPad Software, San Diego, California, USA). The mice experiments were performed following the protocols approved by the by the Ethics Committee of the First Affiliated Hospital of Xi’an Jiaotong University.

All results were statistically analyzed and summarized as a mean and standard error of measurement (SEM) by GraphPad Prism Software. The two-tailed Student’s t-test was applied to assess the differences between the groups. Spearman’s correlation test was performed to determine the correlation between groups. p ≤ 0.05 was considered statistically different.

To determine differences in FGFR family expression between tumor and normal tissues, we analyzed the mRNA level both in various datasets and our GC samples. First, the differential expression of FGFR family members was evaluated with pan-cancer analysis with Oncomine and TIMER databases, which revealed that the expression of FGFR family in lung cancer, breast cancer, liver cancer, colorectal cancer, head-neck cancer, GC, and sarcoma was higher than that in normal tissues ( Figure S1 ). Then the transcriptional levels of FGFR family in GC and their adjacent normal tissues were evaluated in the GEPIA database ( Figure 1A ). According to the result, FGFR2 and FGFR4 expression were considerably increased in tumors than that in normal tissue. But there was no significant difference between the tumors and normal tissues in FGFR3 and FGFR1. Next, we downloaded the original files from TCGA-STAD and GTEx database to analyze FGFR1-4 expression. In this database, the expression levels of FGFR2, FGFR3, and FGFR4 were remarkably higher in STAD than that in normal tissues, while FGFR1 was significantly lower ( Figure 1B ). Furthermore, the analysis of paired samples from the TCGA database suggested that only FGFR4 expression was significantly elevated in tumor samples ( Figure 1C ). To further ensure the accuracy of results we analyzed above, we examined the mRNA level in 27 pairs of GC tissues and their adjacent normal tissues by qRT-PCR. The results indicated that the RNA level of FGFR2 (p<0.05) and FGFR4 (p<0.001) was considerably higher in tumor tissues compared to the matched adjacent non-tumor tissues ( Figure 1D ). No significant difference was found between tumor and adjacent normal tissues in the mRNA level of FGFR1 and FGFR3. Moreover, we further validated the expression of FGFRs in another GSE66229 database. The expression of FGFR3 and FGFR4 were significantly higher in tumor than in normal samples in GC (p<0.001) ( Figure 1E ).

To further assess the prognostic and diagnostic value of FGFRs in STAD, we analyzed the impact of FGFR family expression on survival using the Kaplan–Meier plot. High expression mRNA levels of the whole FGFR family were significantly associated with OS and PFS in GC patients ( Figures 2A, B ). The result from GSE66229 database showed the same trend. The expression of FGFR1 and FGFR4 were also correlated with the poorer overall survival in GC. But FGFR2 and FGFR3 ( Figure 2C ). Furthermore, we analyzed the correlation between FGFRs and GC patients’ characteristics and found FGFR1 expression was substantially associated with the T stage, whereas other clinical features were not markedly different with high expression of FGFR ( Figure 2D , Figure S2 ). FGFR1, FGFR 3, and FGFR 4 expression varied in different histological types, among which the difference between diffuse and tubular tumors was most pronounced ( Figure 2E ). Then, we assessed the correlation between FGFR family expression and the clinicopathological features with univariate and multivariate Cox regression analysis using TCGA data. As shown in Table 1 , N stage, distant metastasis, age, and FGFR4 were independent prognostic factors for GC patients, while FGFR1, FGFR2, and FGFR3 had no significant association with an increased risk of GC and poor OS and PFS in univariate or multivariate analyses. Additionally, the ROC analysis revealed that the FGFRs expression levels had high diagnostic potential for GC, among which FGFR4 was the most accurate (FGFR1: AUC: 0.615, FGFR2: AUC: 0.533, FGFR3: AUC: 0.615, FGFR4: AUC: 0.878; Figure 2F ).

To better understand the underlying mechanisms of FGFRs expression in GC, we evaluated the genes associated with the expression of FGFR family members in the STAD datasets of TCGA with the Limma R package. After the selection, upregulated and downregulated genes were identified according to the FGFR high and FGFR low expression group in STAD. a correlated heatmap was drawn between the different groups, which showed the highest ten positive and ten negative genes of every FGFR member ( Figures 3A–D ) (Pearson Correlation Coefficient, p<0.05). To verify the gene-gene interaction network, STRING and GeneMANIA were also used to explore the correlated genes. The potential target gene interactions with the FGFR family are shown in Figure S3 . The interactions were generally in agreement with the genes from the association analysis. Then we performed a functional enrichment analysis to excavate the underlying biological function of FGFR based on the co-expressed genes. The KEGG analysis showed the FGFR1/3/4 were closely related to the PI3K-AKT signaling pathway and ECM receptor interaction, while FGFR2 was highly correlated with cell cycle, pathways of neurodegeneration multiple diseases, and Wnt signaling pathway ( Figures 4A–D ).

In the results above, we found that FGFR1, FGFR2, FGFR3, and FGFR4 were all associated with the prognosis of GC patients, in which FGFR4 had the most prognostic and diagnostic value among the whole family. To further elucidate the functional role of FGFR4 overexpression, NCI-N87 and MGC-803 cell lines with low FGFR4 expression were transfected with lentivirus plasmid pHAGE-FGFR4. FGFR4 overexpression stable cell line was established by puromycin selection. FGFR4 overexpression was confirmed by qRT-PCR and immunoblotting ( Figures 5A, F ). The overexpression of FGFR4 induced a significant elevation of cell proliferation according to the proliferation assay ( Figure 5B ). Similarly, plate colony formation showed considerably more colonies in the FGFR4 overexpressed group ( Figure 5C ). As to the transwell cell migration assay, high levels of FGFR4 resulted in the increase of both NCI-N87 and MGC-803 cells ( Figure 5D ). The data from scratch assay was consistent with the data above ( Figure 5E ). To ascertain the downstream regulation pathway of FGFR4 overexpression, we examined the expression levels of representative markers in the ERK/MAPK pathway, PI3K/Akt pathway, and ECM signaling. Western blot analyses showed that p-ERK, p-AKT, fibronectin, and matrix metalloproteinase 2 (MMP-2) were remarkably elevated after FGFR4 overexpression. Meanwhile, E-cadherin was significantly reduced in the N87-FGFR4 cells ( Figure 5F , Figure S4 ). The data was consistent with the result from the KEGG pathway analysis and functional experiment. Moreover, the nude mice bearing FGFR4 overexpression tumors (N87-FGFR4) showed better tumor-growth promotion than those with the control plasmid ( Figures 5G–I ). These results indicated that FGFR4 overexpression promotes the malignancy and the EMT in GC. Therefore, FGFR4 can also be expected to serve as a novel target for GC treatment.

Gene mutation was highly correlated with tumor origination and progression; therefore, we speculated that FGFR family variation in gene mutation might affect GC development. Using the cBioportal database, we initially investigated the types and frequencies of alterations in the FGFR family in GC samples. From the histogram, three types of gene alteration of FGFR family members were discovered in various GC cohorts, including mutation, deep deletion, and amplification. The mutation frequency of FGFR2 was the most significant (6%), of which the majority were amplification mutations, whereas the three types of gene alterations were more balanced in FGFR1, FGFR3, and FGFR4 ( Figures 6A, B ). However, in the survival analysis, no considerable differences were discovered between the FGFR2 altered mutation group and the unaltered one ( Figure S5 ). TIMER database was employed to further assess the association between the level of tumor infiltration and SCNA. The variations of FGFRs somatic copy number alterations were all significantly associated with immune cell infiltration levels of CD8+ T-cells, CD4+ T-cells, macrophages, and dendritic cells in GC, suggesting the involvement of FGFRs in regulating immune cells in the TME ( Figure 6C ).

To understand the correlation between FGFR family members and immune cell markers, the TCGA-STAD database was investigated using the GSVA package included in the R software. A positive correlation was found between the expressions of FGFR1 and FGFR2 in GC patients and the level of immune cell infiltration, which included NK cells, CD8+ T-cells, dendritic cells, and macrophages. However, FGFR3 and FGFR4 expressions showed a reverse negative trend ( Figure 7A ). Among all the FGFR family members, the Spearman correlation number between FGFR1 and NK cells was highest (=0.756). While FGFR4 had the most obvious negative correlation with the CD8+ T cells (=-0.240). In contrast, the correlation number of FGFR2 and FGFR3 was much lower. The correlation between enrichment of distinct immune cells and FGFRs are shown in Figures 7B, C . TIMER database analysis also revealed a similar pattern ( Figure S6 ).

To the best of our knowledge, the tumors that could attract more T-cell infiltration are considered as hot tumors, they tend to be more sensitive and effective to immunotherapy. Since we found that FGFR1 was significantly positively correlated with many immune cells and had a high correlation coefficient, we wondered whether overexpression of FGFR1 could form a hot tumor, thereby improving the immunotherapy effect. Furthermore, we evaluated whether a negative correlation with FGFR4 overexpression could diminish the immunotherapeutic consequence. To confirm this, NCI-N87 cells were transfected by pHAGE-FGFR1, pHAGE-FGFR4, and its empty vector. After selecting cells with stable FGFR expression, we confirmed their protein levels through a Western blot assay ( Figure 7D ). Then N87-FGFR1, N87-FGFR4, and N87-EV cells were subcutaneously injected into NOG mice. After tumor formation, the mice were randomly allocated to receive either an anti-PD-1 mAb or a placebo. All three kinds of xenograft NOG mouse models showed TV decrease when treated with anti-PD-1 mAb. Among them N87-FGFR1 had most significant differences, while no considerable difference was found in the N87-FGFR4 group ( Figure S7 ). In addition, N87-FGFR1 group showed significant regression of average TV after anti-PD-1 mAb treatment compared with the control group, while there was no obvious decrease in the TV of the N87-FGFR4 group compared with the control group when treating with anti-PD-1 mAb. Interestingly, a significant difference was found in mouse bodyweight between the N87-FGFR1 and N87-FGFR4 groups ( Figure 7E ). And the survival analysis revealed that the survival time was also significantly elevated in the N87-FGFR1 group after anti-PD-1 treatment compared with the N87-FGFR4 group ( Figure 7F ). As a result, overexpression of FGFR1 may improve the immunotherapeutic impact of GC, whereas overexpression of FGFR4 may impair the immunotherapeutic effect.

Overexpression of FGFR4 is not effective for immunotherapy by itself; therefore, a combination of immunotherapy and FGFR4 inhibitors potentially improve the efficacy of the treatment. We investigated the combination treatment of FGFR4 inhibitor (FGF401) and anti-PD-1 mAb in an FGFR4 overexpressed xenograft mouse model. In a mouse subcutaneous xenograft model of N87-FGFR4, the daily anti-PD-1 treatment showed only modest antitumor efficacy compared with the control mice. Following combination with oral dosing of FGF401 at 30 mg/kg, the TV of the mice decreased remarkably ( Figure 8A ). Survival curves of mice in the FGFR4 overexpression model by different treatment (placebo, FGF401, anti-PD-1 mAb, or combination) were presented in Figure 8B . The combination of FGF401 and anti-PD-1 treatment significantly prolonged the survival compared to the control group. The similar situation was seen in the FGF401 single-agent group. However, no significant difference was found between the anti-PD-1 mAb group and the control group or FGF401 treatment group. In contrast, mice treated with combination reagents survived longer than mice treated with anti-PD-1 mAb alone (HR, 0.235; 95% confidence interval (CI), 0.07–0.74; p= 0.01274). The median survival of the control, anti-PD-1 mAb, FGF401, and combination groups was 16.5 days, 19 days, 21 days, and 30 days, respectively ( Figure 7B ). Collectively, these results demonstrated that a high expression level of FGFR1 could indicate a good prognosis with the anti-PD-1 treatment, while the opposite was seen with FGFR4. Combination therapy with FGF401 and anti-PD-1 could improve the antitumor immune response in high expression FGFR4 GC.