Identification and validation of extracellular matrix-related genes in the progression of gastric cancer with intestinal metaplasia

Abstract

BACKGROUND

Gastric cancer (GC) is a highly lethal malignancy with a high incidence and mortality rate globally. Its development follows the Correa model, with intestinal metaplasia (IM) being a critical precursor to GC. However, the mechanisms underlying IM progression to GC remain unclear. This study explored extracellular matrix (ECM)-related gene changes during IM progression to GC, aiming to identify biomarkers that could improve early diagnosis and treatment strategies for GC, ultimately enhancing patient outcomes.

AIM

To analyze transcriptome sequencing data, molecular biomarkers that can predict GC risk and monitor IM progression can be identified, providing new insights and strategies for preventing IM-GC transformation.

METHODS

Weighted gene co-expression network analysis served for confirming gene modules. Upregulated ECM-related genes were further tested using univariate Cox regression and least absolute shrinkage and selection operator analysis to select hub genes and construct a survival analysis model. The intestinal cell model was established by stimulating GES-1 cells with chenodeoxycholic acid.

RESULTS

Weighted gene co-expression network analysis identified 1709 differentially expressed genes from the GSE191275 dataset, while The Cancer Genome Atlas stomach adenocarcinoma revealed 4633 differentially expressed genes. The intersection of these datasets identified 71 upregulated and 171 downregulated genes, which were enriched in ECM-related pathways. Univariate Cox regression analysis identified six genes with prognostic significance, and least absolute shrinkage and selection operator regression pinpointed secreted protein acidic and rich in cysteine (SPARC) and SERPINE1 as non-zero coefficient genes. A prognostic model integrating clinical tumor node metastasis staging, age, SERPINE1, and SPARC was constructed. Immunohistochemistry analysis confirmed an increasing expression of SPARC protein from normal gastric mucosa (-), to IM (+- to +), and to GC (+ to ++), with significant differences (P < 0.05). Western blot analysis demonstrated significantly higher SPARC expression in induced intestinal cells compared to GES-1. Furthermore, after SPARC knockdown in the human GC cell line HGC27, cell counting kit-8 and colony formation assays showed a reduction in cell proliferative ability, while the wound healing assay revealed impaired cell migration capacity.

CONCLUSION

Comprehensive analysis suggested that a model incorporating clinical tumor node metastasis staging, age, and SPARC/SERPINE1 expression served as a prognostic predictor for GC. Moreover, elevated SPARC expression in IM and GC suggests its potential as a proper biomarker to detect GC in early stage and as a novel therapeutic target, guiding clinical applications.

Key Words: Gastric cancer; Intestinal metaplasia; Extracellular matrix; SPARC; SERPINE1

Core Tip: This study aimed to identify molecular biomarkers for predicting gastric cancer (GC) risk and monitoring intestinal metaplasia progression. Using weighted gene co-expression network analysis, 71 upregulated and 171 downregulated genes related to the extracellular matrix were identified. Secreted protein acidic and rich in cysteine (SPARC) and SERPINE1 were found to be significant prognostic genes. A survival model integrating clinical tumor node metastasis staging, age, and SPARC/SERPINE1 expression was developed. Immunohistochemistry and western blot analyses confirmed the elevated expression of SPARC in intestinal metaplasia and GC. SPA0RC expression correlated with cancer progression and could serve as an early biomarker and therapeutic target for GC.

INTRODUCTION

Gastric cancer (GC) is a highly lethal malignancy. The incidence is the fifth highest, and the mortality is the fourth highest, based on epidemiological data. In 2020, over 1000000 new cases of GC were reported worldwide and 770000 died from GC, seriously threatening human health[1-3]. Advanced GC has a worse prognosis, and early diagnosing and treating of early-stage GC and precancerous lesions are crucial for reducing incidence and mortality rates.

It is widely accepted that GC development follows the Correa model: Chronic gastritis to atrophy of the gastric mucosa to intestinal metaplasia (IM) to dysplasia to GC[4,5]. The progression of GC typically starts with active chronic inflammation, often alongside non-atrophic chronic gastritis or multifocal atrophic gastritis. It advances through stages of complete metaplasia, incomplete IM, and low-grade and high-grade dysplasia, culminating in cancer[6]. Gastric IM means replacing the surface, central, and/or glandular epithelium in the gastric mucosa or glandular mucosa with intestinal epithelium. IM is a critical step in the precancerous process of GC and can greatly induce GC development[7]. However, the mechanism by which IM occurs and progresses to GC remains unclear. Understanding the specific mechanisms underlying the occurrence and progression of IM to GC can aid in identifying possible biomarkers and therapeutic targets specific to the disease, contributing to the early diagnosis and treatment to improve patient prognosis[8]. A review of existing literature indicates that there are currently no definitive biomarkers available for predicting the way IM progresses to GC. Therefore, new biomarkers shall be identified for enhancing the diagnostic process, especially for early detection, increasing the likelihood of successful treatment, and ultimately improving survival rates.

Tumors occur and develop in a complex manner and are regulated by numerous internal factors. The extracellular matrix (ECM), encompassing a network of proteins, proteoglycans, growth factors, and other secreted components, remarkably impacts tumor initiation, growth, invasion, and metastasis[9-11]. The secreted protein acidic and rich in cysteine (SPARC) gene encodes a key protein for tissue development and repair. It impacts the interaction between ECM and cells, influencing cell adhesion, proliferation, metastasis, and growth factor signaling[12]. Therefore, studying extracellular mechanisms and related proteins can help reveal the development mechanism regarding GC and assist in creating new therapeutic strategies.

This study utilized public databases and clinical data to investigate the molecular mechanisms by which IM progressed to GC. It focused on changes in ECM-related genes, a critical factor in GC development, for identifying central genes that affected the progression. The goal is to discover key biomarkers for the progression from IM to GC, providing new insights for diagnosing and treating GC in the early stage.

MATERIALS AND METHODS

Selection and preprocessing of transcriptome data

Expression data for normal gastric mucosa, IM, and GC came from Gene Expression Omnibus (https://www.ncbi.nlm.nih.gov/geo) and The Cancer Genome Atlas (TCGA) (https://www.cancer.gov/ccg/research/genome-sequencing/tcga) databases. The datasets selected were GSE191275 and TCGA-stomach adenocarcinoma (STAD), respectively. We calculated the logarithmic fold change of genes and took any greater than 1 and a P < 0.05 as the threshold for identifying differentially expressed genes (DEGs).

Construction of the co-expression network

The weighted gene co-expression network analysis (WGCNA) package served as the co-expression network analysis[13]. All genes were included based on sample grouping, with the top 50% exhibiting the highest median absolute deviation selected during preprocessing. Genes with values greater than 0.01 were further analyzed. Sample trees and trait heat maps were visualized using WGCNA. The removal of outliers was followed by the transformation of a similarity matrix into an adjacency matrix based on a β value. That was followed by the construction of a scale-free network together with the topological overlap matrix. Genes showing high absolute correlation were clustered into modules, and a cluster tree was generated using the dynamic tree cutting method. Each module had ≥ 100 genes, with a module merging threshold set at 0.3 and a sensitivity of 3.

Hub module and hub genes screening

Using WGCNA, genes were grouped into modules to explore their association with clinical phenotypes. Modular significance was defined as the relevance of module eigengenes to clinical features. The average gene significance within each module, after calculation, prioritized genes associated with IM and GC. Modules of which the correlation coefficient |R| was over 0.3 and P was less than 0.05 were identified as key modules. Within these key modules, we calculated the correlation between gene expression and module eigengene for confirming module membership. Hub genes were identified following the criteria: |Module membership| > 0.7 and |gene significance| > 0.1).

Enrichment analysis

Specific to gene set enrichment analysis, gene functions were annotated using org.Hs.eg.db (version 3.1.0) across three dimensions: Biological process; cell component; and molecular function. Kyoto Encyclopedia of Genes and Genomes (KEGG) REST API helped to obtain the KEGG pathway annotations. We had the genetic map against this background set, performed enrichment analysis under the assistance of clusterProfiler (version 3.14.3) in R, identified significant enrichments with a minimum gene set size of 5 and maximum of 5000, and used thresholds of P < 0.05 and false discovery rate < 0.1.

Prognostic analysis

Univariate and multifactorial regression analyses served for assessing the prognostic value of individual genes. Significant factors identified from the multivariate analysis, such as hub gene risk scores and disease stage, were integrated in a nomogram. The length of each line in the nomogram reflects the relative importance of the variable and the impact on the outcome. The decision curve and clinical impact curve (CIC) were plotted using the R package “rmda.” The decision curve illustrated the net benefit of using the model for decision-making at different threshold probabilities. In the CIC, the red solid line represents the number of patients predicted as “high-risk” by the model at different threshold probabilities, and the blue dashed line represents the number of patients predicted as “high-risk” who actually experienced the event (e.g., disease, death, etc.) at different threshold probabilities. The confidence interval of the curve reflects the stability of the model in identifying the event, with a narrower confidence interval indicating greater stability of the model.

Immunoinfiltration analysis

Using the single-sample gene set enrichment analysis algorithm from R-GSVA [1.46.0][14], we evaluated immune cell infiltration in TCGA-STAD based on 24 markers of immune cells as described in an article on immunity[15]. These markers include activated dendritic cells (DC), B cells, cluster of differentiation (CD) 8 T cells, cytotoxic cells, eosinophils, immature DC, macrophages, mast cells, neutrophils, natural killer (NK) CD56 bright cells, NK CD56 dim cells, NK cells, plasmacytoid DC, T cells, T helper (Th) cells, T central memory, T effector memory, T follicular helper, T gamma delta, Th1 cells, Th17 cells, Th2 cells, and regulatory T cells. We analyzed the relevance of hub genes to the immune infiltration matrix data, and lollipop charts from the ggplot2 package visualized the results.

Verification of central gene expression

For immunohistochemistry (IHC) analysis of gastric endoscopic biopsy tissue sections after fixation in 10% neutral formalin and paraffin embedding, sections were sliced to 4 μm thickness. That was followed by baking, deparaffinization, and hydration. Quenching of endogenous peroxidase activity relied on 3% hydrogen peroxide, and antigen was retrieved at room temperature. After 30 min of blocking in 5% goat serum, primary antibody targeting the SPARC protein was employed and incubated for one night at 4 °C. Then the treated sections underwent 30 min of incubation by a secondary antibody at 37 °C. After PBS wash, sections received diaminobenzidine visualization, hematoxylin counterstaining, dehydration, and microscopic examination after being mounted with neutral gum.

Cell culture and experiments

The cell lines used in this study include the human normal gastric epithelial cell line GES-1 and the human GC cell line HGC27, both purchased from Wuhan Pricella Biotechnology Co., Ltd. (China). The cell lines were cultured at 37 °C in a 5% carbon dioxide (CO₂) humidified atmosphere. GES-1 cells were cultured in RPMI 1640 medium containing 10% fetal bovine serum, while HGC27 cells were cultured in complete medium consisting of RPMI-1640 base medium, 10% fetal bovine serum, and 1% penicillin-streptomycin. According to previous literature[16], a chenodeoxycholic acid (CDCA)-induced IM model was established. GES-1 cells were cultured to 50%-70% confluence, then serum-starved in Dulbecco’s modified eagle medium for 24 h before treatment with 200 μM CDCA for 24 h. Cell protein extracts were collected to assess the expression of intestinal markers (KLF4, CDX2) in the model.

Cell transfection: HGC27 cells at log growth phase were seeded in 6-well plates at a density that ensured approximately 50% confluence the following day. For transfection, 500 μL OPTI-MEM medium, 1.5 μL Plus reagent (Lipofectamine LTX and Plus Reagent, Invitrogen, CN2481212), and 3 μL small interfering (si)-RNA were mixed and incubated for 10 min, followed by addition of 3.75 μL LTX reagent (Lipofectamine LTX and Plus Reagent, Invitrogen, CN2481212). After 30 min of incubation, the transfection complexes were prepared and added to the cells. The cells were then incubated in the transfection medium for 4-6 h before replacing it with complete medium for 1-2 days. Si-normal control (NC), si-SPARC-1, and si-SPARC-2 were synthesized by He Yuan Biotechnology Company.

Cell counting kit 8 assay: The cell counting kit 8 (CCK-8) assay was used to evaluate cell proliferation. Approximately 1000-3000 cells were seeded per well in a 96-well plate. Each group included at least five replicates, and the plates were incubated in a temperature-controlled shaker for 2 h in the dark after adding 10 μL CCK-8 and 90 μL complete medium. Absorbance was measured at 450 nm using a microplate reader (BIO-RAD, Japan) at different time points.

Colony formation assay: The colony formation assay was used to assess cell proliferation ability. Cells were seeded at a density of 400-700 cells per well in 6-well plates and incubated for 2-3 weeks in a CO₂ incubator following different treatments. The cell states and colony sizes were observed under a microscope. When colonies containing ≥ 50 cells were observed, cells were fixed with methanol for 5 min, stained with 0.1% crystal violet, and counted and analyzed using ImageJ software.

Wound healing assay: The wound healing assay was used to assess changes in cell migration ability. Two groups of cells, NC and si-RNA, were cultured for 48 h in 6-well plates. When the cell confluence reached > 90%, a straight scratch was made across the wells using a 10 μL pipette tip, perpendicular to the plate and guided by a ruler. After removing suspended cells with PBS washing, serum-free medium was added, and the cells were incubated in the CO₂ incubator. Photographs were taken at 0 h, 6 h, 12 h, and 24 h to monitor migration.

Statistical analysis

Experimental results were in the form of mean ± SD. Spearman correlation indicated the correlation between central genes and immune cells, and P < 0.05, P < 0.01, and P < 0.001 reported statistical significance. Bioinformatics analyses utilized R software (version 4.2.1), Xiantaozi (https://www.xiantaozi.com), and Sangerbox online tools (http://sangerbox.com/).

RESULTS

Hub genes selection for IM and GC

The GSE191275 dataset, which includes normal controls and IM samples, comprised 20 samples. After data cleaning and filtering, we obtained 17893 genes. The WGCNA package was utilized for mRNA co-expression network analysis, and a β value of 12 was chosen as the optimal value, leading to the identification of eight modules. The modules related to IM were evaluated, revealing that those with higher modular significance values had a stronger correlation with disease progression (Figure 1). The black membership and blue membership modules, associated with IM, were considered hub modules. Following the selection criteria, we ascertained 1709 highly connected genes in the clinically significant module as hub genes. According to the selection criteria, the TCGA-STAD dataset contained 4663 DEGs (2322 upregulated and 2341 downregulated) (Figure 2A).

Figure 1 Flowchart and weighted gene co-expression network analysis. A: Flowchart of research; B: The optimal β soft threshold was determined as 12 through scale-free topology and average connectivity analysis; C: Cluster analysis of characteristic gene modules resulted in the identification of eight modules; D-F: Heatmaps showed the correlation between phenotypes and characteristic gene modules confirmed black membership and blue membership as hub modules. GEO: Gene Expression Omnibus; TCGA: The Cancer Genome Atlas; STAD: Stomach adenocarcinoma; IM: Intestinal metaplasia; GC: Gastric cancer; ECM: Extracellular matrix; LASSO: Least absolute shrinkage and selection operator; WGCNA: Weighted gene co-expression network analysis; DEGs: Differentially expressed genes; CDCA: Chenodeoxycholic acid; WB: Western blot; IHC: Immunohistochemistry; SPARC: Secreted protein acidic and rich in cysteine.

Figure 2 The Cancer Genome Atlas-stomach adenocarcinoma volcano plot and intersection gene analysis. A: Volcano plot for differentially expressed genes (DEGs) between normal and gastric cancer in The Cancer Genome Atlas-stomach adenocarcinoma (STAD); B and C: Venn diagram of the intersection of hub genes in intestinal metaplasia with DEGs in STAD; D and E: Gene Ontology and Kyoto Encyclopedia of Genes and Genomes enrichment results of upregulated DEGs; F and G: Downregulated DEGs; H: The Venn diagram of the intersection of hub genes in intestinal metaplasia, upregulated DEGs in STAD, and genes related to the extracellular matrix; I: A heatmap of the expression of 18 intersecting genes in GSE191275; J: A heatmap of the expression of 18 intersecting genes in The Cancer Genome Atlas-STAD. DEGs: Differentially expressed genes; GO: Gene Ontology; KEGG: Kyoto Encyclopedia of Genes and Genomes; ECM: Extracellular matrix.

An intersection analysis between the hub genes from the WGCNA analysis of IM and the DEGs of GC respectively identified 71 and 171 upregulated and downregulated genes (Figure 2B and C) in TCGA-STAD. These 242 genes were considered cross genes that possibly impact IM development and IM-GC progression.

Functional annotation and enrichment analysis of cross genes

The commonly upregulated and downregulated genes in IM-GC underwent Gene Ontology (GO) and KEGG enrichment analyses for uncovering the mechanisms and processes possibly related to these genes. The GO enrichment analysis of upregulated cross genes identified 638 significant terms, including processes associated with ECM organization, collagen-containing ECM, ECM binding, and positive regulation of ECM organization (Figure 2D). These terms are associated with the ECM. The KEGG enrichment analysis of upregulated genes highlighted significant pathways such as prostate cancer, proteoglycans in cancer, and platinum drug resistance, with proteoglycans being a key component of the ECM (Figure 2E).

For the commonly downregulated genes, the GO enrichment analysis revealed 84 significant terms, including processes such as anterograde axonal protein transport, protein localization to presynapse, axo-dendritic protein transport, and protein transport along microtubules (Figure 2F). The KEGG enrichment analysis of the downregulated genes identified significant pathways like ATP-binding cassette transporters and regulation of the actin cytoskeleton (Figure 2G). Based on these enrichment analyses, it is evident that the IM-to-GC progression is closely associated with changes and regulation of the ECM. Thus, further analysis will focus on the upregulated cross genes related to the ECM.

ECM-related genes and intersecting genes

A total of 1539 ECM-related genes came from the GeneCards website (https://www.genecards.org/). The intersection of the upregulated genes in IM-GC revealed 18 genes: AGRN, BMP1, COL1A1, COL3A1, CXCL1, FERMT1, ITGA2, ITGAX, KIF23, MELTF, MMP9, PDGFRB, PLAU, SERPINE1, SOX9, SPARC, TNFRSF11B, and TOP2A. These 18 genes were identified as key ECM-related genes that are upregulated in IM-GC (Figure 2H). Heatmaps were used to visualize their expressions in the TCGA-STAD and GSE191275 datasets (Figure 2I and J).

Prognostic analysis of intersecting genes and establishment of prognostic models

We conducted single-factor Cox regression analysis on the clinical data of the TCGA-STAD dataset for the 18 intersecting genes, with a forest plot used for visualization (Figure 3A). Six genes, TNFRSF11B, COL1A1, COL3A1, SERPINE1, SPARC, and PDGFRB, showed an obvious relevance to overall survival (OS), and a Kaplan-Meier curve was plotted for these six genes (Figure 3B). Further screening of these genes was performed using prognostic least absolute shrinkage and selection operator (LASSO) regression analysis, determining the λmin value to be 0.01. At this value, the non-zero coefficients were SPARC and SERPINE1 (Figure 3C).

Figure 3 Prognostic analysis of intersecting genes and model establishment. A: Forest plot of univariate regression analysis for 18 genes in The Cancer Genome Atlas-stomach adenocarcinoma dataset; B: Kaplan-Meier curve plot of six genes exhibiting an obvious relevance to overall survival; C: Least absolute shrinkage and selection operator analysis of six genes based on prognostic information from The Cancer Genome Atlas-colon adenocarcinoma. HR: Hazard ratio; CI: Confidence interval; SPARC: Secreted protein acidic and rich in cysteine.

A multivariable regression analysis including clinical tumor node metastasis stage, age, and the two genes (SPARC and SERPINE1) was conducted to establish a prognostic nomogram. The consistency of the model (Concordance, C-index) was 0.666 (0.641-0.691) (Figure 4A). The prognostic calibration curve indicated that the model accurately predicted the OS of patients with GC. The calculation formula of the risk score that was used was (0.2 × expression of S) + (0.3 × expression of SERPINE1). This formula determined each patient’s risk score in the TCGA-STAD dataset. A risk factor plot was generated based on gene expression according to the patients who had high levels of both genes. These patients exhibited a poorer prognosis and a higher risk of death (Figure 4B). The median risk score was considered for dividing patients into two risk groups, and the two groups presented obviously different survivals. The high-risk group presented remarkably poorer OS (Figure 4C; P < 0.001). In addition, time-dependent receiver operating characteristic curves based on the risk score showed the area under the curves of 0.603, 0.636, and 0.674 at 3 years, indicating higher prognostic efficacy (Figure 4D).

Figure 4 The establishment and validation of prognostic models. A: Nomogram and prognostic calibration curve constructed based on clinical tumor node metastasis stage, age, and hub genes; B: Risk factor plot based on secreted protein acidic and rich in cysteine/SERPINE1 expression and risk score in the Cancer Genome Atlas-stomach adenocarcinoma dataset; C: Prognostic analysis and Kaplan-Meier curve for the two risk groups; D: Prognostic receiver operating characteristic curves at 1-year/2-year/3-year based on secreted protein acidic and rich in cysteine/SERPINE1 risk scores in the The Cancer Genome Atlas-stomach adenocarcinoma dataset; E: The decision curve analysis curve; F: The clinical impact curve. HR: Hazard ratio; TPR: True positive rate; FPR: False positive rate; AUC: Area under the curve; SPARC: Secreted protein acidic and rich in cysteine.

Additionally, we performed validation analysis of the model using the decision curve analysis and CIC curves. In the decision curve analysis curve, when the threshold is less than 0.5, the model (orange) demonstrates a net benefit higher than the “All” and “None” lines, indicating that the model offers higher net benefits in low-risk clinical decision-making and is suitable for moderate-risk decisions (Figure 4E). In the CIC curve, when the threshold probability is below 0.5, there is a clear separation between the predicted number of high-risk cases (red line) and the actual number of high-risk cases (blue line), indicating that the model has good clinical discriminatory value within this range. Specifically, in the 20%-40% typical clinical decision threshold range, the model effectively identifies high-risk patients who truly require intervention while avoiding overtreatment of low-risk patients. However, when the threshold probability exceeds 0.5, the two curves begin to overlap, suggesting that the model has limited ability to identify extremely high-risk patients (Figure 4F).

Immunoinfiltration and single cell analysis

This study also investigated the relevance of SPARC and SERPINE1 to immune cells in the TCGA-STAD dataset. In the lollipop chart, both genes were positively correlated with various immune cells (macrophages, NK cells, Th1 cells, DC cells, and mast cells) (Figure 5A and B; P < 0.05). Scatterplots were created for the four immune cells most strongly correlated with SPARC/SERPINE1 (Figure 5C and D), with SPARC showing the strongest correlation with macrophages (r = 0.632).

Figure 5 Secreted protein acidic and rich in cysteine/SERPINE1 immunoinfiltration analysis. A and B: Secreted protein acidic and rich in cysteine (A) and SERPINE1 (B) in the data set the Cancer Genome Atlas (TCGA)-stomach adenocarcinoma (STAD) immune infiltration lollipop plot; C and D: Scatter plots of correlation between secreted protein acidic and rich in cysteine and major immune cells in TCGA-STAD (C) and between SERPINE1 and major immune cells in TCGA-STAD (D). NK: Natural killer; DC: Dendritic cells; iDC: Immature dendritic cells; pDC: Plasmacytoid dendritic cells; TFH: T follicular helper; aDC: Activated dendritic cells; CD: Cluster of differentiation; TPM: Transcripts per million; SPARC: Secreted protein acidic and rich in cysteine; Th: T helper cells; Treg: Regulatory T cells; Tgd: T gamma delta cells; Tcm: T central memory cells.

These findings align with the immunoinfiltration results, suggesting that SPARC and SERPINE1 may influence the formation and regulation of the tumor ECM through their expression in fibroblasts. One possible explanation for these changes is that these genes interact with downstream targets to modify tumor cell cytokine secretion, playing a crucial role in regulating immune function and the ECM.

SPARC/SERPINE1 mutation in GC

Using TCGA-STAD data and analyzing it through the cBioPortal database (https://www.cbioportal.org), the study assessed the genetic variation of SPARC and SERPINE1[17]. The degree of genetic variation differed between the two genes. SPARC showed a low incidence of genetic variation (mutation %, amplification 0%) (Figure 6A). SERPINE1 exhibited a higher incidence of genetic variation (mutation, amplification, deep deletion, and multiple alteration) (Figure 6B).

Figure 6 Genetic variation analysis and immunohistochemistry results in Human Protein Atlas of secreted protein acidic and rich in cysteine/SERPINE1. A: Secreted protein acidic and rich in cysteine mutation analysis results; B: SERPINE1 mutation analysis results; C and D: Immunohistochemical images of Secreted protein acidic and rich in cysteine (C) and SERPINE1 (D) in normal gastric mucosa tissues and gastric cancer. SPARC: Secreted protein acidic and rich in cysteine.

Protein expression in GC tissues

The Human Protein Atlas (https://www.proteinatlas.org/) served for analyzing protein expression in human GC tissue samples. Comparing the IHC results of GC tissues and normal gastric mucosa tissues, SPARC was moderately or strongly expressed in GC tissues, with a significant difference from normal gastric mucosa (Figure 6C). The SERPINE1 immunohistochemical staining did not show obvious difference between them (Figure 6D).

Validation of the prognostic model in GSE15459

Dataset GSE15459 is a GC cohort from Singapore and comprises 200 patients with primary GC[18]. The prognostic efficacy of the SPARC/SERPINE1-based model was validated using GSE15459. Patient risk scores were taken into account to divide them into high-risk and low-risk groups, with the former presenting a poorer OS in the survival analysis (Figure 7A; P = 0.008). A risk factor map, combining the expression of SPARC/SERPINE1 in the GSE15459 dataset, showed that the low-risk group possessed longer survival times and lower risk (Figure 7B). The 1-year, 2-year, and 3-year receiver operating characteristic curves based on risk scores revealed area under the curves of 0.573, 0.576, and 0.609, respectively (Figure 7C).

Figure 7 Prognostic model validation and hub gene analysis. A: Prognosis model in the GSE15459 dataset risk score Kaplan-Meier curve; B: Prognostic model in the GSE15459 dataset risk factor map; C: The time-varying prognostic receiver operating characteristic curve of the prognostic model in the GSE15459 dataset; D: TNFRSF11B COL1A1, COL3A1, SERPINE1, secreted protein acidic and rich in cysteine (SPARC), PDGFRB in diagnosis of least absolute shrinkage and selection operator GSE191275 dataset regression coefficient; E: Two variables with non-zero regression coefficients generated by the least absolute shrinkage and selection operator analysis in GSE191275: SPARC/TNFRSF11B; F and G: SPARC expression in The Cancer Genome Atlas (TCGA)-stomach adenocarcinoma (STAD) unpaired samples and paired samples; H: Expression of SPARC in normal gastric mucosa and intestinal metaplasia samples of GSE191275; I and J: SPARC was expressed in TCGA-STAD dataset at different clinical stages; K: Low and high SPARC expression groups had obviously different survivals in TCGA-STAD, clinical stage IV, T3T4 and N2N3 (P < 0.05); L: Radar map of SPARC expression in pancarcinoma; M: Prognosis and overall survival of SPARC in pancarcinoma. ^aP < 0.001. HR: Hazard ratio; TPR: True positive rate; FPR: False positive rate; TPM: Transcripts per million; AUC: Area under the curve; ACC: Adrenocortical carcinoma; BLCA: Bladder urothelial carcinoma; BRCA: Breast invasive carcinoma; CESC: Cervical squamous cell carcinoma and endocervical adenocarcinoma; CHOL: Cholangiocarcinoma; COAD: Colon adenocarcinoma; DLBC: Lymphoid neoplasm diffuse large B-cell lymphoma; ESCA: Esophageal carcinoma; GBM: Glioblastoma multiforme; HNSC: Head and neck squamous cell carcinoma; KICH: Kidney chromophobe; KIRC: Kidney renal clear cell carcinoma; KIRP: Kidney renal papillary cell carcinoma; LAML: Acute myeloid leukemia; LGG: Brain lower grade glioma; LIHC: Liver hepatocellular carcinoma; LUAD: Lung adenocarcinoma; LUSC: Lung squamous cell carcinoma; MESO: Mesothelioma; OV: Ovarian serous cystadenocarcinoma; PAAD: Pancreatic adenocarcinoma; PRAD: Prostate adenocarcinoma; PCPG: Pheochromocytoma and paraganglioma; READ: Rectum adenocarcinoma; SARC: Sarcoma; SKCM: Skin cutaneous melanoma; STAD: Stomach adenocarcinoma; TGCT: Testicular germ cell tumors; THCA: Thyroid carcinoma; THYM: Thymoma; UCS: Uterine carcinosarcoma; UVM: Uveal melanoma; UCEC: Uterine corpus endometrial carcinoma; SPARC: Secreted protein acidic and rich in cysteine.

Evaluation of diagnostic efficacy of key genes in GSE191275

To further evaluate the diagnostic efficiency of six key genes (TNFRSF11B, COL1A1, COL3A1, SERPINE1, SPARC, PDGFRB) for IM prognosis, LASSO diagnostic analysis was performed on GSE191275. The λmin was 0.02 (Figure 7D). At this value, two non-zero coefficients, SPARC and TNFRSF11B, were identified (Figure 7E).

Characteristic analysis and experimental verification of SPARC

Through the previous analyses, SPARC demonstrated good prognostic efficacy in GC and showed strong diagnostic efficacy in both IM and GC. Therefore, we further analyzed the characteristics of SPARC and verified its expression in clinicopathological tissues and cells. In the unpaired and paired samples from the TCGA-STAD dataset, the GC group exhibited dramatically higher SPARC expression relative to the normal control group (Figure 7F and G, P < 0.001). In the GSE191275 dataset, the IM group presented remarkably higher SPARC expression relative to the normal gastric mucosa group (Figure 7H, P < 0.001). In the TCGA-STAD prognosis subgroup analysis, females, T3/T4, pathological stage IV, and subgroup N (N2 and N3) with high SPARC expression had significantly worse OS (Figure 7I-K, P < 0.05).

In the pan-cancer analysis, SPARC expression was significantly different across various tumor tissues, including breast invasive carcinoma, colon adenocarcinoma, kidney renal clear cell carcinoma, prostate adenocarcinoma, with most being upregulated (Figure 7L). Additionally, SPARC expression had good prognostic efficacy in several tumors, including bladder urothelial carcinoma, kidney renal clear cell carcinoma, and STAD (Figure 7M) (P < 0.05).

The validation of SPARC

We validated the expression of SPARC in normal gastric mucosal epithelium, IM, and GC tissues through IHC. SPARC protein was negative (-) in normal gastric mucosa, showed variable expression (+- to +) in IM, and was elevated (+ to ++) in GC, with statistically significant differences (P < 0.05) (Figure 8A). We established an IM model in the GES-1 cell line by CDCA treatment and used KLF4 and CDX2 as markers for IM to validate the model. The results showed increased expression of KLF4 and CDX2 in CDCA-treated GES-1 cells, indicating successful model establishment. Based on this, we observed elevated SPARC expression in the IM cell model (Figure 8B).

Figure 8 Secreted protein acidic and rich in cysteine expression in tissues and cells. A: The immunohistochemistry results of secreted protein acidic and rich in cysteine (SPARC) in normal gastric mucosa, intestinal metaplasia (IM), and gastric cancer tissue; B: The western blot results of IM markers KLF4, CDX2, and SPARC in GES-1 cell line and chenodeoxycholic acid-induced GES-1 IM cell model; C: The validation results of SPARC knockdown by small interfering (si)-RNA through WB and quantitative real-time polymerase chain reaction; D: The line graph assessing cell proliferation by cell counting kit 8 assay after SPARC knockdown with small interfering (si)-RNA; E: The results of the colony formation assay after SPARC knockdown with si-RNA; F: The results of the cell wound healing assay after SPARC knockdown with si-RNA. ^aP < 0.001, ^bP < 0.01, ^cP < 0.05. NC: Normal control; IM: Intestinal metaplasia; GC: Gastric cancer; IHC: Immunohistochemistry; DMSO: Dimethyl sulfoxide; CDCA: Chenodeoxycholic acid; OD: Optical density; si-NC: Small interfering-normal control; si-SPARC: Small interfering-secreted protein acidic and rich in cysteine; SPARC: Secreted protein acidic and rich in cysteine.

To further explore the biological function of SPARC in GC, we performed SPARC knockdown in the human GC cell line HGC27 via transfection. The knockdown efficiency of SPARC in GC cells was assessed by quantitative PCR and western blotting. Both si-SPARC sequences showed significant knockdown effects, with the si-2 sequence demonstrating the most efficient knockdown (Figure 8C). We selected si-SPARC-2 for subsequent experiments. CCK-8 assays revealed that the cell proliferation rate in the si-NC group was significantly higher than in the si-SPARC group, with a statistically significant difference (Figure 8D). Colony formation assays showed that cells in the si-SPARC group formed smaller colonies and fewer colonies compared with the control group, with a statistically significant difference (P < 0.05) (Figure 8E). Wound healing assays demonstrated that the wound healing ability (migration capacity) of cells in the si-SPARC group was significantly impaired, with a statistically significant difference compared with the control group (P < 0.01) (Figure 8F).

DISCUSSION

GC is a prevalent disease globally, posing a significant threat to human health. Often, GC is diagnosed at an advanced stage, leading to a poor survival prognosis. Recently, there has been increasing attention on the precancerous lesions of GC. Monitoring and intervention efforts are being made for these lesions, including chronic atrophic gastritis, IM, and atypical hyperplasia. Helicobacter pylori is a major pathogenic factor in GC[19]. However, according to recent studies, atrophy can be effectively reversed by Helicobacter pylori eradication therapy[20], while IM has a very low reversal rate[21]. This suggests that IM, as an intermediate stage in the Correa cascade, may involve more gene or molecular regulatory changes, significantly increasing the risk of GC. Thus, focusing research on IM as a special precancerous state can provide early warning of GC risk and benefit more patients. Exploring genetic and epigenetic changes in the development and progression of IM can provide new insights into GC and precancerous lesions. Therefore, we utilized transcriptome data to characterize the development of IM and its transformation into cancer, aiming at confirming genes with diagnostic and prognostic potential.

Our study analyzed transcriptome data from IM and GC, identifying intersecting differential genes common to both processes. Functional analysis of these shared genes revealed enrichment in cellular activities and related pathways, such as immune response and the ECM. Consequently, we propose that these related genes may contribute to the progression of IM and assist IM in developing into GC via the regulation of the ECM.

In this study, WGCNA was employed to determine the IM-related modules and hub genes in the GSE191275 dataset using WGCNA, and the common differential genes in the normal gastric mucosa-IM-GC process were obtained by intersecting with TCGA-STAD differential genes. Enrichment analysis of the intersecting genes indicated that the ECM might be a key factor in this process, particularly regarding upregulated differential genes. Based on univariate Cox regression and LASSO-Cox regression analyses, we built a risk score prognostic model comprising two marker genes (SPARC/SERPINE1). Samples with high-risk scores presented poor prognoses. The calibration curves demonstrated good accuracy for the risk score-based nomogram model, and validation across different cohorts confirmed the strong predictability of the model. The prognostic model effectively distinguished high-risk from low-risk groups based on age (> 65 years), T stage, M stage, and N stage, demonstrating good predictive power for various clinical factors.

SPARC is an ECM glycoprotein that binds to calcium ions, interacts with various proteins in the ECM, and competes with cell membrane surface receptors for growth[22]. Evidence suggests that SPARC regulates cell adhesion, cell proliferation, cell migration, and tissue remodeling during cell development and ECM renewal[23]. A 2010 study found that high SPARC expression showed a close association with GC progression and poor patient prognosis[24], indicating that SPARC can independently predict tumor risk, thereby assessing tumor invasion depth and lymph node metastasis[25]. SPARC interacts with collagen types I, III, and IV, promoting ECM stiffness and facilitating cancer cell invasion[26]. It also enhances the expression of matrix metalloproteinases, which degrade the basement membrane, thereby enabling tumor cell dissemination. However, its role in the development of gastric prelesions remains unexplored.

SERPINE1 regulates the plasminogen activation system. Highly expressed SERPINE1 reports poor prognosis across 21 different cancers in pan-cancer analyses[27]. SERPINE1 crucially impacts the immune response to malignancies, potentially contributing to immunosuppressive effects. Moreover, it has been associated with tumor mutational burden, response to immunotherapy, microsatellite instability, and sensitivity to diverse cancer treatments[28,29]. Additionally, research has shown that SERPINE1 may regulate the expression of vascular endothelial growth factor and interleukin-6 through the vascular endothelial growth factor signaling pathway and the Janus kinase-signal transducer and activator of transcription inflammatory signaling pathway, ultimately influencing the invasion and migration of GC cells[30].

These genes, SPARC and SERPINE1, directly or indirectly influence the formation of the ECM and immune cell infiltration, thereby contributing to tumor development. In our study, SPARC has demonstrated significant potential and efficacy in both prognosis and diagnosis, making it suitable for GC and gastrointestinal metaplasia. We confirmed significant step-like increases in SPARC expression from NC to IM and GC stages using IHC and different stages of cellular protein western blot. Thus, SPARC is a promising candidate target for IM-GC, effectively predicting the IM progression, early GC screening, and diagnosis and serving as an important intervention target for the IM-GC process. SPARC and SERPINE1 gene expression levels can serve as predictive markers for identifying patients who are at high-risk for GC who require more frequent monitoring. Regular oxygen and glucose deprivation screenings (e.g., every 1-3 years) should be recommended for patients with higher gene expression, with monitoring frequency adjusted based on individual risk profiles and clinical factors.

Our work characterized the transition from IM to cancer in detail, screened for genes instrumental in diagnosing IM and GC, and identified potential biomarkers for GC. These findings can inform relevant research with regard to the IM-GC transformation. Further research and analysis of SPARC will undoubtedly be beneficial. Our findings place overall emphasis on the value of the SPARC gene in GC patients and will benefit future research.

However, there were limitations. We only analyzed transcriptome data from two datasets, which could introduce selection bias. Additionally, some crucial genes may have been overlooked in the multistep selection process, potentially restricting the applicability of risk scores. Future research should focus on validating this risk model in clinical settings to assess its predictive accuracy and enhance its clinical utility. Moreover, larger sample sizes are needed for validation, along with experimental validation in cellular and animal models for more deeply exploring the action mechanism of these targets. This direction outlines our future research goals and strategies.

CONCLUSION

In this study, we identified two key markers in GC, SPARC and SERPINE1, which serve as predictors of poor progression and as prognostic biomarkers for IM and GC. Using bioinformatics methods, we constructed prognostic risk assessment models for these two gene traits using data from the TCGA and Gene Expression Omnibus databases. Whether the model could reliably and stably predict the prognosis of patients with GC was verified. This study improved our understanding of the relationship between GC and IM. These genetic markers hold promise as prognostic biomarkers and may be valuable therapeutic targets for patients with GC.