Acquisition of LMR genes and GC mRNA expression profile
Gene Set Enrichment Analysis (GSEA) was conducted to gain LMR genes on the website (https://www.gsea-msigdb.org/gsea/index.jsp). And we used following 6 most common pathways as the key words: metabolism of lipid, sphingolipid metabolism, fatty acid metabolism, glycerophospholipid metabolism, peroxisome proliferator-activated receptor alpha and transcriptional regulation of white adipocyte differentiation. The mRNA expression profile (FPKM) of GC (tumor tissues = 375, adjacent tissues = 32) was obtained from TCGA (https://portal.gdc.cancer.gov/), then we converted ENSG_ID to Gene Symbol via GENCODE data (https://www.gencodegenes.org/human/), eliminated the genes and samples whose missing values were more than 50%, and further used the impute.knn function of R package “impute” to complete the missing values, the nearest neighbor Kappa index was 10, and finally carried out standardization processing based on log2 (X + 1).
LMR differentially expressed genes (DEGs) and Consensus clustering
We analyzed the LMR-DEGs between GC and adjacent tissues using R package “limma” (version 3.40.6) with the cut-off criteria set as P < 0.05, and fold change (FC) ≥ 1.5 [11]. Using the R package “org.HS.eg.db” (version 3.1.0), we performed the Gene ontology (GO) functional annotation and Kyoto Kyoto encyclopedia of genes and genomes (KEGG) enrichment analysis [12] to clarify the biological process, subcellular localization, molecular function and involved pathway between the DEGs and total LMR-genes. Basing on the LMR-DEGs level, we employed Consensus-ClusterPlus to identify various metabolic subtypes in 375 GC patients [13], utilizing agglomerative “Partitioning Around Medoids” (PAM) clustering with a 1-Pearson correlation distances and resampling 80% of the samples for 10 repetitions. The optimal number of clusters was determined using the empirical cumulative distribution function plot. We obtained clinical data from UCSC Xena (https://xenabrowser.net/), and drew the Kaplan-Meier curve of each subtype, indicators including overall survival (OS), PFI and disease specific survival (DSS).
Weighted correlation network analysis (WGCNA)
After analyzing the pathway and prognosis differences in metabolic subtypes, we used R passage WGCNA to screening LMR-module [14]. Firstly, the Pearson’s correlation matrices and average linkage method were both performed for all pair-wise genes, then a weighted adjacency matrix was constructed. A soft-thresholding parameter, β, was set to 5 with R2 = 0.86 in our study to highlight strong correlations between Genes while penalizing weaker correlations. The adjacency was transformed into a topological overlap matrix (TOM), which could measure the network connectivity of a gene defined as the sum of its adjacency with all other genes for network gene ration, and the corresponding dissimilarity (1-TOM) was calculated. Genes with comparable expression profiles were grouped into gene modules using average linkage hierarchical clustering based on the TOM-derived dissimilarity measure, with a minimum size of 60 for the genes dendrogram. With setting the sensitivity to 3 and module merge threshold to 0.25, we finally obtained several modules and identified the best module correlated to metabolic subtypes. To establish a predictive model rather than explore mechanism network, we further screened predictive genes of PFI (P < 0.01) via univariate cox regression rather than analyse the hub genes.
Construction and optimization of predictive model
According to the proportion of 3:2, we divided GC patients into the training set (N = 214) and verification set (N = 141) by SPSS 23.0 's randomized generator, and used t or Chi-square test to verify the balance of clinical characteristics. Then we performed the LASSO regression using the “glmnet” package [15] to construct a original model to predict PFI in the training set, and the efficacy was evaluated by time-ROC curve even in the verification set. If the overall efficiency or the stability in internal verification was poor, multivariate cox regression (step back wald method) and clinical characteristics were used to increase efficiency, portability and internal stability. In different data sets, time-ROC, KM and decision curve analysis (DCA) curve were used to evaluate the efficiency, stability and clinical benefit of the final model.
Single sample gene set enrichment analysis (ssGSEA) to identify the changeable pathways in GC subtypes
We divided GC patients into high and low groups based on the final model score and NPR3 expression, and performed ssGSEA (version 3.0) to find enrichment pathways in different subgroups. Then we set the minimum gene set to 5, the maximum gene set to 5000, and resample a thousand times. Pathways with P < 0.05, false discovery rate (FDR) < 0.25, and normalized enrichment score (NES) > 1.5 or < -1.5 were considered to be related. Before Vitro experiment, we retrieved the online data of NPR3 from The Human Protein Atlas (HPA, https://www.proteinatlas.org/) and Gene Set Cancer analysis (GSCA, http://bioinfo.life.hust.edu.cn/GSCA/#/) to determine the subsquent verification method.
Cell culture and transfection
Cell line GES-1 were obtained from the cell bank of SAIBAIKANG biotechnology company (Shanghai, China), while other 5 human GC cell lines were obtained from the Chinese Academy of Sciences (SNU-1, AGS, NCI-N87, Beijing; HGC27, Kunming; MKN45, Shanghai). All cells were cultured in RPMI-1640 medium (Biological Industries, Cromwell, CT, United States) supplemented with 10% fetal bovine serum (FBS, Biological Industries, Cromwell, CT, United States), and maintained in a 37°C incubator with a 5% CO2. We used the cloning vector (pcDNA3.1-NPR3, insert size: 1638bp, Tsingke, Nanjing) to overexpress NPR3, while empty vector (pcDNA3.1-vector, Tsingke, Nanjing) was used as the negative control. The AGS cells with 60-80% confluence in 6 wells dishes were transfected with 2ug recombinant plasmid and 5ul Hieff Trans® Liposomal Transfection Reagent (Yeasen, Shanghai, China) without FBS. After 6 hours, Changing the medium containing 10% FBS and continue to culture to 48h.
RNA extraction and RT-qPCR assay
Total RNA was extracted from human cell lines using TRIzol reagent (Invitrogen, Carlsbad, CA, United States) according to the manufacturer’s instructions. The reverse transcription reaction (RT) was performed with Reverse Transcription kit (Vazyme, Nanjing, China). The RT-qPCR reactions were performed with a SYBR Green PCR Kit (Vazyme, Nanjing, China), measured in Quadruplicate and performed on an Applied Biosystems 7900HT sequence detection system (Applied Biosystems). β-actin was used as an internal control for mRNA. The relative expression levels of the target genes were calculated using the comparative 2−ΔΔCt method. All primers used in this study were listed in Supplementary Table S1, and repeated three times.
Western blotting
Samples were homogenized in RIPA buffer, sonicated, and prepared for Western blotting. Protein were separated by electrophoresis on 10% SDS-PAGE gels before semi-drytransfer to the PVDF membrane. NPR3 and β- actin were detected using antibodies purchased from Abcam (EPR12716, Berlin, Germany) and Proteintech. All results were repeated three times.
Wound healing assay and transwell assay
AGS cells were seeded into six-well plates and grown to 90-100% confluence. The cell layer was scrapping with a 1000μL sterile pipette tip and washed with three times with PBS to remove the floating and detached cells, then cultured without FBS medium and images were acquired under a microscope at two time points (0, 24 h). Cell migration abilities were evaluated by transwell chambers (Corning Life Sciences, Bedford, MA, United States). Briefly, a total of 3 × 104 AGS cells suspended in 200ul media without FBS were seeded in the upper chambers. Then, 500 μL of culture medium containing 20% FBS was added to the lower chambers. After incubation at 37°C for 24 h, the cells in the lower chambers were fixed with 4% paraformaldehyde for 30 min, stained with crystal violet (Beyotime, Shanghai, China) for 30 min. Finally, two random fields were microscopically examined. All results were repeated three times.
Immunohistochemistry (IHC), H&E and sirius red staining
A total of 42 pathologically confirmed GC tissues were obtained from Nanjing Drum Tower Hospital, and the expression of NPR3 (Rabbit, 26706-1-AP, proteintech) was detected by IHC staining. Staining of NPR3 was scored by the two pathologists blind to the clinical data by applying a semi-quantitative immunoreactivity score (IRS), which determined by the product of staining intensity and staining area under imageJ. In our study, IRS of 0–1 and IRS of 2–4 were classified as low and high expression of NPR3, respectively. In addition, we stained these tissues by the Sirius Red Kit (PHYGENE,PH1098) to determine the degree of fibrosis.
Statistical analysis methods
Generate graphics using GraphPad Prsim 8 software. In addition to R analysis, the data were analyzed using SPSS (version 22.0, Chicago, IL, USA), with statistical tests including chi square test, log rank test, and student t-test. All tests were conducted bilaterally, with P<0.05 indicating a statistically significant difference.