Exploration of the clinicopathological and prognostic significance of BRCA1 in gastric cancer.

A total of 2536 DEGs, including 1432 upregulated genes and 1104 downregulated genes, were screened in GSE13911 (Fig.  3 A). In GSE29272 , a total of 6280 DEGs, 3011 of which were upregulated and 3269 of which were downregulated, were identified (Fig.  3 B). Additionally, 11,258 DEGs, including 10,704 upregulated genes and 554 downregulated genes, were identified in the TCGA dataset. To identify the common DEGs between the TCGA and GEO datasets, Venn diagram analysis was performed; the results revealed 792 common DEGs (Fig.  3 C). Using these common DEGs, a PPI network was subsequently constructed; the network contained 748 nodes and 14,009 edges (Fig.  3 D). According to the degree values of the nodes, the target genes were ranked and the top nine core DEGs, (GADPH, CDK1, TOP2A, CCNB1, CDC20, AURKB, BRCA1, RAD51 and AURKA) were further screened (Fig.  3 E). All our results support that these nine core genes are indeed differentially expressed in GC. Fig. 3 Screening of the DEGs. A Volcano plot of the DEGs in GSE13911 . B Volcano plot of the DEGs in GSE29272 . C Common DEGs in the GSE13911 , GSE29272 and TCGA datasets. D Protein‒protein interaction network of the DEGs. The nodes represent the target genes. E The PPI network of the top nine core DEGs Screening of the DEGs. A Volcano plot of the DEGs in GSE13911 . B Volcano plot of the DEGs in GSE29272 . C Common DEGs in the GSE13911 , GSE29272 and TCGA datasets. D Protein‒protein interaction network of the DEGs. The nodes represent the target genes. E The PPI network of the top nine core DEGs To test the hypothesis described above, the mRNA expression levels of the nine core DEGs in GC and normal tissues were compared. All nine core DEGs, including GADPH, CDK1, TOP2A, CCNB1, CDC20, AURKB, BRCA1, RAD51 and AURKA, presented significantly greater expression in GC tissues than in normal tissues (p value < 0.001, F i g.  4 A–I), validating our assumption. Fig. 4 The mRNA expression levels and the results of the GO and KEGG analyses of the core DEGs. A – I mRNA expression levels of AURKA, AURKB, BRCA1, CCNB1, CDC20, CDK1, GADPH, RAD51, and TOP2A in normal and GC tissues. J GO enrichment analysis of the core DEGs. K KEGG pathway analysis of the core DEGs. ***p < 0.001 The mRNA expression levels and the results of the GO and KEGG analyses of the core DEGs. A – I mRNA expression levels of AURKA, AURKB, BRCA1, CCNB1, CDC20, CDK1, GADPH, RAD51, and TOP2A in normal and GC tissues. J GO enrichment analysis of the core DEGs. K KEGG pathway analysis of the core DEGs. ***p < 0.001 Furthermore, to explore the biological functions and pathways of these nine core DEGs in the development of GC, GO and KEGG enrichment analyses were performed. Three categories, i.e., biological process (BP), cellular component (CC) and molecular function (MF), were included in the GO analysis. BP analysis revealed that the core DEGs are associated mainly with the regulation of the cell cycle and the G2/M phase of mitosis. Additionally, CC analysis revealed that these core DEGs are significantly related to condensed chromosomes, the germ cell nucleus and the spindle. Moreover, the MF analysis revealed that the nine core DEGs are mainly involved in the abnormal expression of histone kinase activity, ubiquitin-like protein ligase binding, protein serine kinase activity, protein serine/threonine kinase activity and ATP-dependent activity (Fig.  4 J). In addition, the KEGG results demonstrated that cellular senescence, platinum drug resistance and the p53 signalling pathway were the major potential mechanisms associated with the nine core DEGs (Fig.  4 K). These results demonstrate the crucial roles of these core DEGs in the development of GC. To investigate the impact of the identified core genes on the survival of patients with GC, Kaplan‒Meier survival analysis was carried out. Among these core DEGs, the expression level of BRCA1 was closely related to the survival time of GC patients, with statistically significant differences (p < 0.05), and the other eight genes did not significantly affect survival (Supplementary Fig. 1A, B and D–I). Furthermore, through the observation of the Kaplan‒Meier survival curve of BRCA1 depicted in Supplementary Fig. 1C, we found that the OS of GC patients in the low-risk group was significantly worse than that of those in the high-risk group. These findings suggest that BRCA1 might play an essential role in GC development and progression and may merit further study. To investigate the biological pathways and functions of BRCA1 in the development of GC, single-gene GSEA was performed. The results revealed that the five top functional enrichment terms significantly associated with BRCA1 were chromosome segregation, DNA-dependent DNA replication, DNA replication, chromosomal region, and helicase activity (Fig.  5 A). In addition, pathway enrichment analysis of genes related to BRCA1 was carried out, and the results revealed that the five pathways strongly associated with BRCA1 included the cell cycle, DNA replication, olfactory transduction, the spliceosome and ubiquitin-mediated proteolysis (Fig.  5 B). Notably, high expression of BRCA1 promotes the cell cycle, DNA replication, olfactory transduction, and the spliceosome, whereas low expression of the BRCA1 gene can promote ubiquitin-mediated proteolysis. These functions and pathways are strongly correlated with GC progression. For example, DNA and RNA sequencing analysis revealed that BRCA1, which is involved in olfactory transduction, participates in tumour cell proliferation, migration, and invasion, demonstrating that olfactory transduction may provide insights into new anti-GC treatments. Fig. 5 Biological function and relationship of BRCA1 with the clinical characteristics of patients with GC. A GO and B KEGG pathway enrichment analysis of BRCA1. C Immunohistochemistry analysis of BRCA1 expression in GC and normal tissues via HPA. D Relationship between the BRCA1 expression level and PFS in GC patients. E The ROC curve and AUCs of BRCA1 at 1, 3 and 5 years, in which the higher the AUC value is, the stronger the prediction ability of the gene is. p < 0.05 was considered to indicate a statistically significant difference Biological function and relationship of BRCA1 with the clinical characteristics of patients with GC. A GO and B KEGG pathway enrichment analysis of BRCA1. C Immunohistochemistry analysis of BRCA1 expression in GC and normal tissues via HPA. D Relationship between the BRCA1 expression level and PFS in GC patients. E The ROC curve and AUCs of BRCA1 at 1, 3 and 5 years, in which the higher the AUC value is, the stronger the prediction ability of the gene is. p < 0.05 was considered to indicate a statistically significant difference To further assess the protein distribution of BRCA1 within human tissues and cells, immunohistochemistry (IHC) images of BRCA1 in GC and normal tissues were retrieved from the Human Protein Atlas (HPA), and the results are depicted in Fig.  5 C. BRCA1 was clearly overexpressed at the protein level in GC tissue compared with normal tissue, which is consistent with the results from the TCGA database. Additionally, BRCA1 expression and its correlation with clinical characteristics and the PFS in the TCGA cohort were also analysed. The results of PFS analysis revealed that BRCA1 expression had no significant effect on prognosis (Fig.  5 D). Moreover, ROC analysis revealed that BRCA1 has disappointing predictive ability in GC patients according to the AUCs at 1, 3, and 5 years (Fig.  5 E). To further explore the associations of BRCA1 with the clinical features of patients with GC, a nomogram analysis was carried out for disease stage, age, sex, grade and the TNM staging system. The results demonstrated that BRCA1 was significantly differentially expressed between different age groups (p  65 years was greater than that in patients aged ≤ 65 years (Fig.  6 A). However, there was no difference in the mRNA expression level of BRCA1 between the sexes (p > 0.05, Fig.  6 B). In addition, the expression level of BRCA1 was not significantly different among the grade groups. Conversely, the expression level of BRCA1 in Grade 2 tumours was slightly greater than that in Grade 1 and Grade 3 tumours, but the difference was not statistically significant (p > 0.05, Fig.  6 C). Moreover, among the different pathological stages, BRCA1 was not differentially expressed in stages I, II, III or IV (p > 0.05, Fig.  6 D). There was no significant difference in the expression of the BRCA1 gene among the TMN stages (Fig.  6 E–G). In addition, the Nomogram predictive model, which integrates multiple clinical and pathological features, predicted a 1-year survival rate of 0.884, a 3-year survival rate of 0.687, and a 5-year survival rate of 0.594 for GC patients (Fig.  6 H). Subsequently, the model was calibrated using the calibration curve (F i g.  6 I). The correlations between the expression of BRCA1 and clinical characteristics, including age, sex, stage, T stage, M stage, and N stage, in patients with GC were also visualized via a heatmap (Fig.  6 J). These results indicate that the expression of BRCA1 was not related to clinical features except for age, pathological stage or T, N, or M stage. Fig. 6 Relationships between BRCA1 expression and the clinicopathological factors of GC patients. A – G Correlations between BRCA1 mRNA expression and clinical characteristics in patients with GC. Age, sex, grade, stage, T stage, M stage, and N stage. H Construction of a BRCA1-based nomogram for predicting the OS of GC patients. I Calibration curve and Hosmer–Lemeshow test of BRCA1-based nomograms in the TCGA–GC cohort for OS. J Heatmap of the relationships between the expression of BRCA1 and the clinicopathological factors of GC patients in the high- and low-risk groups. The factors from top to bottom are age, sex, stage, T stage, M stage, and N stage (*p < 0.05) Relationships between BRCA1 expression and the clinicopathological factors of GC patients. A – G Correlations between BRCA1 mRNA expression and clinical characteristics in patients with GC. Age, sex, grade, stage, T stage, M stage, and N stage. H Construction of a BRCA1-based nomogram for predicting the OS of GC patients. I Calibration curve and Hosmer–Lemeshow test of BRCA1-based nomograms in the TCGA–GC cohort for OS. J Heatmap of the relationships between the expression of BRCA1 and the clinicopathological factors of GC patients in the high- and low-risk groups. The factors from top to bottom are age, sex, stage, T stage, M stage, and N stage (*p < 0.05) To investigate the causal relationship between BRCA1 and GC, a bidirectional MR approach was used, and the results were estimated using odds ratios (ORs) and 95% confidence intervals (CIs). As shown in Supplementary Fig. 2A and B, forward MR analysis suggested that BRCA1 is not significantly associated with GC (OR 0.982 and 95% CI 0.936–1.031; p = 0.470). At the same time, similar estimates were obtained via four other methods, i.e., MR Egger, weighted median, simple mode and weighted mode, which are summarized in Supplementary Table 1. Additionally, the results of the Cochran’s Q test in IVW (Table  3 , Q = 46.149, p = 0.233) and MR-Egger (Supplementary Table 2, Q = 44.288, p = 0.258) methods suggest that there is no obvious heterogeneity in these SNPs (Supplementary Fig. 2C). Additionally, we found that when any single SNP was removed, no significant directional horizontal pleiotropy was observed via the leave-one-out approach (Supplementary Fig. 2D), which was confirmed by the results of MR-Egger regression (Supplementary Table 2, intercept =  − 0.013, p = 0.208). Similarly, the reverse MR analysis based on IVW also confirmed that there was no causal relationship between GC and BRCA1 (OR 0.998 and 95% CI 0.946–1.054; p = 0.958), which is depicted in Supplementary Table 1 and Supplementary Fig. 3A, B. In addition, no significant associations were observed when the other four methods were employed. Additionally, the results of the Cochran’s Q test in IVW (Supplementary Table 2, Q = 51.998, p = 0.191) and MR-Egger (Supplementary Table 2, Q = 51.993, p = 0.164) methods demonstrated that there was no heterogeneity in these SNPs (Supplementary Fig. 3C). Moreover, there was no obvious directional horizontal pleiotropy in the observed association, as indicated by Supplementary Fig. 3D and Supplementary Table 2 (intercept = 5.2 * 10 −4 , p = 0.949). To further validate the association of NADPH with GCs, MR analysis was also carried out. As expected, the IVW results revealed that NADPH significantly increased the risk of GC (OR 1.113, 95% CI 1.013–1.223; p = 0.026) in the forward MR analysis (Supplementary Fig. 4A, B; Table  2 ). In addition, a similar estimate was obtained via the weighted median approach (OR 1.124 and 95% CI 1.001–1.262; p = 0.048), indicating that NADPH is strongly correlated with GC. According to the MR-Egger regression analysis, there was no incidence of potential pleiotropy (Supplementary Table 3, intercept = 0.034, p = 0.195). Additionally, the results of the Cochran’s Q test in IVW (Q = 6.840, p = 0.336) and MR-Egger (Q = 4.604, p = 0.466) methods combined with the funnel plot collectively suggested that there was no obvious heterogeneity in the observed associations (Supplementary Table 3; Supplementary Fig. 4C). Furthermore, even when any single SNP was removed, there was no significant change in the observed association in the leave-one-out analysis (Supplementary Fig. 4D), suggesting the stability of the result. Table 2 NADPH and its association with Gastric Cancer in the MR analyses Exposure Outcome No. of SNPs Methods OR (95% CI) β (SE) p Forward MR analysis  NADPH GC 7 MR Egger 0.946 (0.750–1.192)  − 0.056 (0.118) 0.656 Weighted median 1.124 (1.001–1.262) 0.117 (0.059) 0.048 IVW 1.113 (1.013–1.223)  − 0.107 (0.048) 0.026 Simple mode 1.189 (0.964–1.467)  − 0.173 (0.107) 0.157 Weighted mode 1.138 (1.006–1.286) 0.129 (0.063) 0.086 Reverse MR analysis  GC NADPH 15 MR Egger 0.777 (0.617–0.979)  − 0.252 (0.118) 0.052 Weighted median 0.899 (0.781–1.033)  − 0.107 (0.071) 0.134 IVW 0.906 (0.812–1.011)  − 0.099 (0.056) 0.079 Simple mode 0.812 (0.646–1.021)  − 0.208 (0.117) 0.096 Weighted mode 0.897 (0.776–1.036)  − 0.109 (0.074) 0.162 NADPH and its association with Gastric Cancer in the MR analyses Conversely, the IVW results from the reverse MR analysis (Table  2 , OR 0.906 and 95% CI 0.812–1.011; p = 0.079) indicate that the GC is not exposed to NADPH (Supplementary Fig. 5A, B). Moreover, similar results were verified in other models (Supplementary Fig. 5A). There was no significant heterogeneity (Supplementary Fig. 5C; Supplementary Table 3) or pleiotropy (Supplementary Fig. 5D; Supplementary Table 3) identified in these models. To further investigate the association between BRCA1 and NADPH, MR analysis was also performed. After the relationship between NADPH and GC was determined, the MR method was also utilized to analyse the causal relationship between BRCA1 and NADPH. As described in Table  3 , the scatter plot (Supplementary Fig. 6A) and forest plot (Supplementary Fig. 6B) revealed that the expression of BRCA1 was strongly associated with NADPH (OR 1.682 and 95% CI 1.479–1.913; p = 2.590 * 10 −15 ), which is consistent with the results of the MR Egger (OR 1.711 and 95% CI 1.195–2.449; p = 0.022), weighted median (OR 1.621 and 95% CI 1.347–1.951; p = 3.310 * 10 −7 ), simple mode (OR 1.593 and 95% CI 1.143–2.219; p = 0.025) and weighted mode (OR and 95% CI 1.587, 1.235–2.037; p = 0.007) analyses. Table 3 BRCA1 and its association with NADPH in the MR analyses Exposure Outcome No. of SNPs Methods OR (95% CI) β (SE) p Forward MR analysis  BRCA1 NADPH 9 MR Egger 1.711 (1.195–2.449) 0.537 (0.183) 0.022 Weighted median 1.621 (1.347–1.951) 0.483 (0.095) 3.31 * 10 −7 IVW 1.682 (1.479–1.913) 0.520 (0.066) 2.59 * 10 −15 Simple mode 1.593 (1.143–2.219) 0.465 (0.169) 0.025 Weighted mode 1.587 (1.235–2.037) 0.462 (0.128) 0.007 Reverse MR analysis  NADPH BRCA1 4 MR Egger 0.838 (0.345–2.036)  − 0.177 (0.453) 0.734 Weighted median 0.983 (0.870–1.111)  − 0.017 (0.063) 0.786 IVW 1.254 (0.702–2.241) 0.227 (0.296) 0.444 Simple mode 1.005 (0.865–1.167) 0.005 (0.077) 0.956 Weighted mode 0.983 (0.872–1.109)  − 0.017 (0.062) 0.803 BRCA1 and its association with NADPH in the MR analyses Additionally, MR-Egger regression analyses revealed that there was no incidence of potential directional horizontal pleiotropy between BRCA1 and NADPH (Supplementary Table 4, intercept =  − 0.0047, p = 0.9226). Furthermore, the results of the Cochran’s Q test with the MR-Egger method (Q = 9.8255, p = 0.1987) and the IVW method (Q = 9.8397, p = 0.2764) suggested that there was no absence of heterogeneity in the association between BRCA1 and NADPH (Supplementary Table 4; Fig.  6 B, C ). In addition, the results of leave-one-out sensitivity analyses (Supplementary Fig. 6D) indicated that removing any single SNP had no obvious effect on the results. In contrast, the results of the five methods of reverse MR analysis indicated that NADPH was not causally related to BRCA1, which was confirmed by the results of the main IVW method (OR 1.254 and 95% CI 0.702–2.241; p = 0.444, Table  3 ; Supplementary Fig. 7A, B). Additionally, the results of heterogeneity analysis (Supplementary Fig. 7C) combined with the Cochran’s Q test in IVW (Supplementary Table 4, Q = 89.471, p = 2.85 * 10 −19 ) and MR-Egger (Supplementary Table 4, Q = 54.279, p = 1.63 * 10 −12 ) methods revealed that NADPH and BRCA1 have significant heterogeneity and horizontal pleiotropy (Supplementary Fig. 7C). In addition, the leave-one-out sensitivity analyses and MR-Egger regression analyses (intercept = 0.181, p = 0.373) suggested that there was no obvious pleiotropy of NADPH or BRCA1 (Supplementary Table 4; Supplementary Fig. 7D).

For this study, a total of seven datasets [i.e., data for cohorts A, B, C, D, E, F and G derived from The Cancer Genome Atlas (TCGA), the Gene Expression Omnibus (GEO), the NHGRI-EBI Catalogue of Human Genome-wide Association Studies (GWAS Catalogue), and the IEU OpenGWAS database] were used. Among them, datasets A and B were obtained from the TCGA; dataset A (for cohort A) contains information on 36 normal and 412 cancer samples, while a total of 443 clinical data records (with several attributes, such as survival time, survival status, age, sex, and cancer stage) are included in dataset B. To retrieve all available GC transcriptome data for cohort A and clinical data for cohort B, the ‘tcgabiollinks’ R package was employed. In addition, dataset C ( GSE13911 ) and dataset D ( GSE29272 ) were obtained from the GEO; dataset C comprises information on 31 normal samples and 38 tumour samples, and dataset D includes information on 134 normal samples and 62 GC tumour samples. Moreover, a total of 167,122 records, including 7921 from the control group and 159,201 from the experimental group, are contained in dataset E obtained from the GWAS Catalogue. In addition, datasets F and G were generated from the IEU OpenGWAS database, and the corresponding cohorts included 3301 different individuals. To extract the differentially expressed genes (DEGs) between the GC and control groups, the limma package in R software (version 4.3.0) was used for Cohort A, with a filtering threshold of p < 0.001 and a false discovery rate (FDR) < 0.001. With respect to cohorts C and D in the GEO database, the “limma” package was also used, and then log2(x + 1) transformation and the “normalizeBetweenArrays” function were employed. By setting the filtering criteria of p value < 0.001, adj-p value  0, the DEGs in GSE13911 and GSE29272 were obtained. Volcano plots of the DEGs were generated via the ggplot2 software package. Using the “VennDiagram” package in R software, the common DEGs between the TCGA and GEO databases were obtained. By setting the minimum required interaction score to 0.4, the protein‒protein interaction (PPI) network for these common DEGs was subsequently built via the STRING database and further visualized via Cytoscape software (version 3.82). In this network, the nodes represent the common DEGs, and the edges represent the associations between common DEGs. To evaluate network centrality and identify the hub genes, an important topological parameter, i.e., degree, was further analysed via the CytoHubba plugin within Cytoscape. The greater the degree is, the more important the node. To further explore the associations between the obtained hub genes and GC, a biological function analysis of these key genes was carried out. First, log2(x + 1) transformation was performed on the expression values of these hub genes. The clusterProfiler package within R software (version: 4.3.0) was subsequently used to perform enrichment analysis of the Gene Ontology (GO) functions and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways of these genes. Finally, the data were subsequently visualized via the ggplot2 package. In an attempt to assess the impact of the identified hub genes on the survival outcome and survival time of GC patients, a clinical risk prognostic model was established with cohort B based on the following formula. 1 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text{risk }\text{score}={\sum }_{i=1}^{n}coefi\times {\text{id}},$$\end{document} risk score = ∑ i = 1 n c o e f i × id, where coefi represents the coefficient and id represents the normalized count for the gene. Then, on the basis of the calculated median risk score, high- and low-risk groups were created for the GC patients. Using Kaplan–Meier survival analysis with the “limma” and “survival” packages within R software, the expression of the identified hub genes was investigated to evaluate the overall survival rates of the high- and low-risk groups. The key genes were defined as those included in the risk score model when the difference in survival rate between the high- and low-risk groups was significant (p value < 0.05). To further clarify the mechanisms of action of the key genes in GC, GSEA, an in silico analysis approach was performed on the basis of the expression levels of the genes in the high- and low-risk groups of patients with GC via the “clusterProfiler” package. GSEA reveals differences in the enrichment of predefined gene sets between two biological states. Additionally, Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were carried out via the “org. Hs.Eg.Db” R package and visualized with the “ggplot2” package. The significantly enriched GO or pathway terms were screened according to a threshold p value < 0.05. To assess the reliability and accuracy of the obtained risk prognostic model, progression-free survival (PFS) analysis was carried out, and receiver operating characteristic (ROC) curves were generated. PFS analysis can be used as an indicator of the extent of cancer progression to some degree, while ROC curves enable the determination of sensitivity and specificity of key genes in prognosis prediction. By employing the “timeROC” package, the prognostic effects of the key genes on patient survival at 1, 3, and 5 years were obtained. The larger the area under the curve (AUC) is, the greater the accuracy of the obtained prognostic risk model. To reveal the correlations of the key genes with clinical characteristics, including disease stage, age, sex, grade and TNM stage, a clinical nomogram model for patients with GC was built using the “survival” and “RMS” packages in R software. Through nomogram analysis, the 1-, 3-, and 5-year OS probability of patients with GC was evaluated. The nomogram and calibration curves were subsequently depicted by using the R package “regplot”. Furthermore, to explore whether there is a causal relationship between key genes and GC, MR analysis was performed. Owing to the strong correlation between nicotinamide adenine nucleotide phosphate oxidase (NADPH) and GC identified in previous studies [ 16 , 17 ], MR analysis of the key genes and NADPH was also carried out. Subsequently, the pairwise random combinations of the key genes, NADPH, and GC were performed with these factors as either exposure or outcome variables (Fig.  2 ). Fig. 2 Workflow of bidirectional MR in this study Workflow of bidirectional MR in this study Additionally, to retain appropriate single nucleotide polymorphisms (SNPs) in the exposure and outcome data, a genome-wide significance threshold and the linkage disequilibrium (LD) were set. The filtering criteria for the SNPs are summarized in Table  1 . After palindromic SNPs with intermediate allele frequencies were removed, the pertinent genetic variants were integrated into our MR analysis and used as instrumental tools for further research. Table 1 The filter criteria for SNPs Exposure Outcome p value r 2 Number of retained SNPs Key gene GC 5 * 10 −5 0.001 41 GC Key Gene 5 * 10 −5 0.001 45 NADPH GC 5 * 10 −6 0.001 7 GC NADPH 5 * 10 −6 0.001 15 Key gene NADPH 5 * 10 −6 0.001 9 NADPH BRCA1 5 * 10 −7 0.001 4 The filter criteria for SNPs Presently, a total of five different methodologies are employed, including MR Egger, weighted median, random-effects inverse variance-weighted (IVW), simple mode and weighted mode approaches. Among them, the random-effects IVW method is the principal statistical approach for estimating potential bidirectional causality between the key gene, NADPH and GC, which relies on all SNPs being valid instrumental variables and completely independent of each other. The remaining four methods served as supplementary approaches to assess causality, thereby enabling us to obtain precise estimations of causal relationships. To evaluate the heterogeneity between the key genes, NADPH and GC, Cochran’s Q test was performed via the IVW and MR-Egger methods [ 18 , 19 ]. If the p value of Cochran’s Q test was > 0.05, heterogeneity was not detected. Additionally, using the intercept of MR-Egger, the degree of bias stemming from directional pleiotropy was assessed. If the p value was < 0.05, possible pleiotropy in IVs was considered. Moreover, to assess whether individual SNPs influence the primary causal relationship, sensitivity analysis was also performed via leave-one-out and single-SNP analyses. A significance threshold of p < 0.05 was used for all analyses in the study.

In conclusion, the expression of BRCA1 in GC significantly differs between normal gastric tissues and GC tissues. Correlation analysis of gene expression revealed that BRCA1 is notably expressed in GC tissues and expressed at low levels in normal tissues. According to survival analysis and clinical sample analysis, the survival rate of patients with GC in the high-BRCA1 expression group was significantly higher than that in the low-BRCA1 expression group, and proportion of patients with high BRCA1 expression was higher among patients younger than 65 years than among those older than 65 years. However, most other clinical factors were not found to be correlated with BRCA1. These findings led us to speculate that BRCA1 affects the occurrence of GC through other pathways. Through MR analysis, we found that there was no direct causal relationship between BRCA1 and GC, but BRCA1 had a s notable causal relationship with NADPH, and NADPH also had a notable causal relationship with GC. Our conjecture was also supported by the results of reverse MR analysis. The results support that BRCA1 indirectly regulates the expression of NADPH and promotes the occurrence of GC, so targeting these factors may represent a new strategy for GC treatment. Overall, we propose that BRCA1 may serve as a new therapeutic target, but the function of this gene has not been explored in depth. Further additional research and experimental verification are necessary.

Gastric cancer (GC) is a significant driver of cancer-related death worldwide. With the ongoing enhancement of our comprehension of the molecular facets of GC, oncologists are now realizing that these malignancies constitute a heterogeneous group characterized by remarkably distinctive pathogenic mechanisms and active oncogenic pathways [ 20 , 21 ]. Hence, a more comprehensive understanding of GC pathogenic mechanisms and oncogenic pathways is important for prevention and management. To discover prognostically valuable tumour biomarkers, we employed bioinformatics to identify significant differentially expressed genes (DEGs) in GC from TCGA and GEO datasets [ 9 ]. Utilizing a set of rigorous filtering criteria, we successfully pinpointed several genes demonstrating notable distinctions between healthy individuals and GC patients (F i g.  4 A–I). These genes included GADPH, CDK1, TOP2A, CCNB1, CDC20, AURKB, BRCA1, AURKB and RAD51. An analysis of the biological processes revealed that these core genes are mainly related to the regulation of the cell cycle; notably, dysregulation of the cell cycle can lead to the occurrence of cancer. An analysis of enriched cell components revealed that these core genes act mainly on the spindle, further highlighting the role of these genes in cell cycle regulation and the formation of cancer and thus to the risk of carcinogenesis. The molecular functions of these core genes were analysed, and they were found to regulate histone protease kinase activity, protein serine kinase activity and ATP-dependent activity. These activities affect gene transcription and multiple processes of cell cycle regulation and cell metabolism, highlighting the functions of the core genes at different levels [ 22 ]. In addition, pathway enrichment analysis of the core genes revealed that the core genes were significantly correlated with cellular senescence, platinum drug resistance and the p53 signalling pathway. For example, the BET inhibitor JQ1 suppresses the proliferation and invasion of GC cells by inducing cellular senescence [ 23 ]. The effectiveness of platinum-based chemotherapeutic drugs, such as cisplatin and oxaliplatin, for treating platinum-resistant GC is limited [ 24 ]. The p53 signalling pathway is an important pathway in the treatment of GC [ 25 ]. By regulating the p53 signalling pathway, the growth of GC can be effectively inhibited. Next, we conducted additional survival analysis on the top nine core DEGs. To our surprise, we found that the expression of these genes had no significant effect on the overall survival (OS) of GC patients. These results deviated from our initial expectations and were disappointing. However, BRCA1 stands out as a special case. Our analysis of the OS of patients with GC revealed significant differences between BRCA1 expression groups. Breast cancer gene 1 (BRCA1) is a crucial tumour suppressor gene consisting of an 1863-residue nuclear polymorphism featuring a prominent N-terminal RING domain and an evolutionarily conserved C-terminal BRCT domain [ 26 ]. Notably, 20–30% of familial breast and ovarian cancers are linked to BRCA1 mutations [ 27 ]. Evidence indicates that the BRCA1 gene plays roles and exhibits differential expression in cancers other than breast and ovarian cancers [ 28 , 29 ]. Single-sample gene enrichment analysis of BRCA1 revealed that multiple biological functions and pathways were correlated (Fig.  5 A, B ). For example, olfactory transduction was found to be a potential anti-GC treatment target. Additionally, we performed immunohistochemical analysis of BRCA1 gene expression in GC tissues and paired normal tissues via the Human Protein Atlas (HPA) as a reference (Fig.  5 C). These investigations revealed differences in cell morphology and BRCA1 staining between normal and GC tissues. However, the PFS (Fig.  5 D) and the AUC (Fig.  5 E) analyses did not lead to conclusive results. Although BRCA1 was not correlated with most clinicopathological factors (Fig.  6 A–J), interestingly, we found that BRCA1 expression was significantly different across age groups. The BRCA1 gene is expressed at higher levels after 65 years of age and is a tumour suppressor gene that can potentially inhibit the occurrence and metastasis of GC. This finding is consistent with the results of previous survival analyses. The multiple population components of GC are reportedly younger, and the incidence of GC in individuals younger than 50 years is increasing annually. This leads us to speculate that there are other pathways through which the BRCA1 gene may be indirectly regulated, thereby influencing the onset or progression of GC. Maintaining balanced redox homeostasis is of utmost importance for human health, as its disruption serves as a pivotal factor in the aetiology of several major diseases, including type 2 diabetes, atherosclerosis, chronic obstructive pulmonary disease, Alzheimer's disease, and more [ 30 ]. Additionally, the occurrence and progression of cancer are intricately linked to oxidative stress [ 31 ]. When cells are exposed to elevated levels of oxidative stress, the risk of cancer development increases significantly [ 32 ]. Once oxidative stress reaches a certain threshold, it stimulates the unrestricted proliferation of cancer cells. NADPH has traditionally been regarded as a pivotal cofactor in antioxidant defence and reductive biosynthesis, supplying high-energy electrons for these biological processes. During cellular oxidative stress, the levels of reactive oxygen species (ROS) increase. In cases where the production of NADPH fails to effectively mitigate the accumulation of excessively high ROS levels, a sustained increase in oxidative stress occurs, promoting carcinogenesis [ 33 ]. Notably, several studies facilitated notable advancements in the investigation of cellular oxidative stress. Yan et al . revealed that the overexpression of SLC7A11 combined with glucose starvation depletes intracellular NADPH and can trigger disulfide stress and potentially cell death [ 34 ]. Additionally, Sang-Min Jeon reported that, under stress conditions, the activation of AMPK can prolong cell survival through redox regulation [ 35 ]. These findings collectively suggest that oxidative stress might offer a pathway for uncontrolled cell proliferation. On the basis of existing research, Lu YX and coauthors explored the modulation of NADPH through the ME1 gene to facilitate the development and metastasis of GC [ 36 ]. Additionally, Lovren F. and colleagues discovered that BRCA1 is a novel and previously unexplored target within vascular smooth muscle cells, VSMCs [ 37 ]. It functions by inhibiting the production of NADPH Nox1-dependent ROS, effectively dampening intracellular oxidative stress. In alcoholic fatty liver disease, the regulatory role of BRCA1 in NADPH has also been discovered. Danielle M. et al. revealed that BRCA1 can play a therapeutic role in nonalcoholic fatty liver disease by inhibiting NADPH generation [ 38 ]. Consequently, we propose the following hypothesis: NADPH, as a key contributor of energy for oxidative–reductive defence and reductive biosynthesis, might require rapid replenished during the rapid proliferation of tumour cells. The expression of NADPH might impact the proliferation and migration of cancer cells. Additionally, BRCA1 might indirectly influence cancer occurrence by modulating NADPH expression. Therefore, BRCA1 may be a potential therapeutic target in GC patients, and targeting agents may suppress the proliferation and spread of cancer cells. To test the hypotheses outlined earlier, we conducted bidirectional Mendelian randomization experiments to empirically assess pairwise interactions between BRCA1, NADPH, and GC. Five methods were employed in this study, including MR Egger regression, the weighted median approach, the IVW method, the simple mode approach and the weighted mode approach. Specifically, MR Egger regression serves as a robust method for detecting and correcting directional pleiotropy in MR analyses. By allowing for a non-zero intercept in the regression model, MR Egger regression provides robust estimates of pleiotropy; however, violations of the assumption of instrument independence may introduce bias. The weighted median approach enables the estimation of causal effects even when over 50% of genetic instruments are invalid. This method proves effective in the presence of significant heterogeneity among genetic variants but, due to its limited statistical power, serves primarily as a supplementary reference. The IVW method is widely utilized in MR analyses. Assuming all genetic variants function as valid instruments, IVW method provides more precise estimates with higher weights. However, the IVW method is highly sensitive to violations of MR assumptions, particularly pleiotropy, which may bias causal estimates when pleiotropy exists. The simple mode approach is relatively straightforward, relying on the assumption that the distribution of causal effect estimates represents the true causal effect, thereby obviating the need for complex statistical modeling. However, this method may fail to account for the precision of individual genetic variants in multimodal distributions, potentially leading to biased results. The principles and limitations of the weighted mode approach are similar to those of the simple mode approach, though it estimates causal relationships with greater weighting. Like other methods, the weighted mode approach does not account for potential pleiotropy or other sources of bias [ 39 – 41 ]. Considering the advantages and limitations of these methods, the IVW method was employed as the primary criterion in our study, with other approaches served as supplementary references. Our findings suggested that there does not appear to be a direct causal link between BRCA1 and the development of GC (Supplementary Figs. 1A, 2A; Supplementary Table 1; Table  2 ). In MR analysis, unmeasured or uncontrolled confounding factors in publicly available datasets, or weak associations of exposure variables with IVs, may yield unexplained negative results despite the large sample size of this study [ 42 , 43 ]. Therefore, although clinical studies have reported differences in BRCA1 expression levels among various GC cell lines, MR analysis cannot establish a direct causal relationship between BRCA1and GC [ 44 ]. However, BRCA1 expression is significantly upregulated in GC patients, so we speculated that there may be intermediates between BRCA1 and GC that indirectly affect the occurrence of disease, thereby reinforcing the conclusions drawn from our earlier survival analysis [ 45 , 46 ]. Furthermore, we discovered a strong correlation between BRCA1, an exposure factor, and NADPH, an outcome (Supplementary Fig. 6A). With increasing risk linked to the exposure factor of the BRCA1 gene, the probability of elevated NADPH content as an outcome steadily increases. A thorough assessment of pleiotropy and heterogeneity was subsequently carried out. Pleiotropy suggests that instrumental variables do not affect the outcome through exposure, whereas heterogeneity indicates variations among different instrumental variables. This analysis indicated the absence of pleiotropy and heterogeneity between BRCA1 and NADPH (Supplementary Figs. 6C, D, 7C, D; Supplementary Table 4). Similarly, when NADPH was considered the exposure factor and GC was considered the outcome of the Mendelian analysis, a significant correlation was also identified (Supplementary Fig. 4A). Furthermore, horizontal pleiotropy and heterogeneity analyses were performed concurrently, revealing that there was no connection (Supplementary Fig. 4C, D; Supplementary Table 3). In contrast, when reverse inference was employed, neither the influence of NADPH on BRCA1 nor the effect of GC on NADPH exhibited significant associations (Supplementary Figs. 5A, 7A). This additional evidence supports our initial hypotheses. In summary, our research supports the interplay between BRCA1 and NADPH, as well as the relationship between NADPH and GC. These findings suggest that BRCA1 might influence the onset and progression of GC through its regulation of NADPH. On the basis of our earlier research, increased activity in the antioxidant pathway was observed in cells transfected with BRCA1. This phenomenon may be a result of oxidative phosphorylation activation, leading to the generation of ROS [ 47 ]. Hence, we concluded that BRCA1 may serve as a promising therapeutic target in GC and is worthy of further research. The remaining eight genes examined in this study may also be involved in GC development through indirect pathways. Further investigations can be conducted in subsequent experiments. However, this study has several limitations. The primary limitation stems from the absence of experimental validation. While our bioinformatics analyses yielded useful insights, experimental validation is helpful to provide evidence for confirming the biological significance of BRCA1 and NADPH in GC. Second, despite the relatively large sample size and the application of statistical methods to minimize selection bias and batch effects, the data may still fall short of capturing the full diversity of the disease. Third, clinical validation is needed. Conducting clinical studies to evaluate the correlation between BRCA1, NADPH, and GC could contribute to improving the reliability of our findings. Therefore, future work will focus on integrating experimental validation and collecting additional clinical GC cases to corroborate these findings.

In 2020, gastric cancer (GC) accounted for 1.1 million cancer cases and 770,000 fatalities [ 1 ]; GC is the fifth most common malignancy and the fourth leading cause of cancer-related mortality [ 2 ]. The incidence of GC exhibits geographical heterogeneity on a global scale, with the highest rates observed in East Asia, particularly in Japan and Mongolia, and Eastern Europe. In contrast, Northern Europe, North America and various African regions consistently have lower incidence rates [ 3 ]. Since Billroth performed the inaugural GC surgery in the nineteenth century, the landscape of GC treatment has undergone a profound transformation [ 4 ]. The cornerstone of GC management is curative surgery. To mitigate the risk of recurrence and increase long-term survival, an array of therapeutic modalities, including perioperative chemotherapy, adjuvant chemotherapy, and adjuvant chemoradiotherapy [ 5 ], has also been developed. These approaches have become common and are recommended in current guidelines for GC treatment [ 6 , 7 ]. Additionally, ongoing investigations are being performed to assess the application of targeted therapy and immune checkpoint inhibitors (ICIs) in the adjuvant setting. Despite the continual expansion of treatment options for GC, the absence of well-defined clinical indications frequently leads to late-stage diagnosis in the majority of patients, thereby leading to an unfavourable prognosis. The incidence of GC in young people (under 50 years of age) has been steadily increasing in both low-risk and high-risk countries [ 8 ]. It is imperative to adopt a fresh perspective on the management of GC. In recent years, there have been many reports on machine learning in medicine. A mounting body of evidence underscores the crucial importance of cancer biomarkers in clinical practice. Early cancer diagnosis can be facilitated by biomarkers, effectively enabling the detection of potential cancer transformation through the assessment of biomarker presence, levels, and alterations. For example, EIF4G1 expression levels may indicate breast cancer risk, and PMEPA1 isoforms have emerged as promising biomarkers and therapeutic targets in prostate cancer [ 9 , 10 ]. Furthermore, biomarkers have notable research value in cancer prognosis assessment. The expression status or levels of specific biomarkers can indicate disease prognosis, patient survival duration or the propensity for disease onset, potentially providing opportunities to mitigate the risk of postoperative recurrence. Chen et al. reported that GPC2 expression levels differ across various tumour types and normal tissues, indicating its increased sensitivity and specificity in cancer diagnosis [ 11 ]. Consequently, GPC2 has emerged as a prospective prognostic marker for a spectrum of malignancies. Mendelian randomization (MR) is a statistical method in biostatistics used to infer causal relationships between two genetically complex traits from observational studies. Owing to the random allocation of alleles in genetics and their free recombination and the stability of genotypes after birth, MR can reduce confounding factors and eliminate reverse causality, making it akin to the “most natural” randomized controlled trial, RCT [ 12 ]. MR utilizes instrumental variable (IV) analysis to infer causality between exposure and the use of genetic variations, particularly single-nucleotide polymorphisms, SNPs [ 13 ]. Hence, MR studies are more compelling than traditional epidemiological studies. There is evidence to suggest that an elevation in tissue protease H levels causes an increase the risk of lung cancer in patients, indicating its potential value in lung cancer assessment. Yan et al . also identified a causal relationship between endometriosis (EMS) and a reduced overall risk of breast cancer [ 14 ]. In subgroup analysis on the basis of immunohistochemical type, EMS was correlated with a lower risk of oestrogen receptor-positive breast cancer, and no causal relationship was observed between oestrogen receptor-negative breast cancer and survival rates. Numerous studies suggest that MR is a widely used research approach [ 15 ]. However, to our knowledge there is a noticeable lack of MR research exploring the prospective causal association between gene expression and GC risk exists. The principal aim of this investigation was to perform an exhaustive bidirectional MR analysis, with the intention of revealing the concealed mechanisms by which differentially expressed biomarkers in GC promote disease onset. This study leveraged GC data obtained from public databases to identify nine potential biomarkers, which might serve as valuable therapeutic targets for precise clinical interventions in patients with GC. A bidirectional MR analysis of the BRCA1 gene, characterized by substantial differential expression within the GC context, was subsequently conducted. The results indicate that BRCA1 may exert an indirect effect on GC development by regulating NADPH levels, which provides potential drug targets for GC. Figure  1 depicts the workflow of this study. Fig. 1 Workflow of this study Workflow of this study

Additional file1 (DOCX 1150 kb) Additional file1 (DOCX 1150 kb)