Comprehensive Analysis of Peritoneal Metastasis Sequencing Data to Identify LINC00924 as a Prognostic Biomarker in Gastric Cancer

Introduction

Gastric cancer peritoneal metastasis was less sensitive to chemotherapy, and the patients diagnosed with peritoneal metastasis had a poor clinical prognosis and high mortality.¹ Cytoreductive surgery combined with hyperthermic intraperitoneal chemotherapy could delay tumor progression; however, the current medical treatment was still limited for prolonging patients’ survival.² Therefore, it is necessary to study the molecular mechanism of gastric cancer peritoneal metastasis which could improve early diagnosis and provide basis for targeted therapy.

Long non-coding RNAs (lncRNAs) were a type of RNAs longer than 200 nucleotides and could not be translated into protein.³ lncRNA-involved in gene regulation were associated with many biological processes including tumor growth, cancer metastasis and immune escape.^4–6 Huang reported that LncRNA AK023391 promoted gastric cancer invasion through PI3K/Akt signaling pathway.⁷ LncRNA SNHG11 activated Wnt/beta-Catenin axis to facilitate cell migration and proliferation in gastric cancer.⁸ However, the function of lncRNAs in gastric cancer peritoneal metastasis was still unclear.

Immune checkpoint therapy – as a promising anti-cancer therapy method – had benefited many patients and extended overall survival.⁹ LncRNA was associated with immunity closely.¹⁰ Li’s research showed that the positive feedback loop of LncRNA SNHG14/miR-5990-3P/ZEB1 promoted diffuse large B cell lymphoma progression through regulating PD-1/PD-L1 checkpoint.¹¹ However, due to individual differences in gene expression, most cancer patients did not respond to immunotherapy.¹² Therefore, the study of gastric cancer peritoneal metastasis gene expression profile and the relationship between different genes and immune infiltrations was beneficial to improve the curative effect of immunotherapy.

Our study results showed that compared with primary foci, there were 3433 lncRNAs and 1210 mRNAs up-regulated, 1341 lncRNAs and 3192 mRNAs down-regulated in peritoneal metastasis lesion. Thirty-eight up-regulated lncRNAs and 9 down-regulated lncRNAs with known function and fold change >2 were selected to construct CNC network. Adipocytes were a major cell type of peritoneal microenvironment¹³ and we used adipocytes to co-culture with tumor cells to mimic the tumor’s internal environment and exclude the influence of cell diversity on sequencing results. Our results showed that LINC01151, FAM27B and LINC00924 were up-regulated in co-culture group and LINC00924 was associated with poor survival. Furthermore, immune checkpoint blockers (ICBs) had a worse effect on the LINC00924 high expression group. We also found that LINC00924-related PEX5L in CNC network was an independent prognostic factor by univariate and multivariate Cox regression analysis in gastric cancer progression.

Materials and Methods

Patient Tissue Samples

In this study, we collected three paired gastric cancer with peritoneal metastasis tissue samples at the Zhongnan Hospital of Wuhan University (Additional File, Table 1). The patient had not experienced any chemotherapy or radiotherapy. All samples were obtained with informed patients’ consent before collection, and approved by the Zhongnan Hospital of Wuhan University Ethics Committee. Moreover, the collected samples need to be fixed in formalin and embedded in paraffin for further study.

Microarray of Primary Lesion and Peritoneal Metastasis Lesion

Biopulverizer was used for cryopulverization of three paired tissue samples, then Mini-Bead-Beater-16 was used for homogenization. Quick Amp Labeling Kit, One-Color (Agilent p/n 5190–0442) was used for labeling reaction, RNeasy Mini Kit (Qiagen p/n 74104) for RNA purifying, NanoDrop ND-1000 for labeled cDNA QC. After hybridization and microarray wash, Agilent Microarray Scanner (Agilent p/n G2565BA) was used for scanning. All data were extracted using Agilent Feature Extraction Software.

Coding-Non-Coding (CNC) Gene Co-Expression Network

The normalized data of mRNA was taken to calculate the medium value of different transcripts of the same gene to represent expression value of gene. Then, the Pearson’s correlation coefficient (PCC) was calculated between selected lncRNAs and all mRNAs according to normalized data. The mRNAs with abs (PCC) ≥0.95, P-value ≤0.05 and FDR ≤0.1 were recorded. Cytoscape was used to draw the CNC gene co-expression network of lncRNAs and mRNAs.

Cell Culture and Adipocyte Differentiation Induction

The human gastric cancer cell lines (AGS, BGC823) culturing in RPMI 1640 medium (Gibco, USA) and mouse preadipocyte cell line (3T3-L1) culturing in DMEM (Gibco, USA) were purchased from Chinese Academy of Sciences in Shanghai. For 3T3-L1 induction, the medium of 3T3-L1 cells was changed with DMEM containing 10% FBS, 10µg/mL insulin, 1µmol/L dexamethasone, 0.5mmol/L 3-isobutyl-1-methylxanthine (IBMX, HY-12318; MCE) for 2 days cultivation after 2 days of contact inhibition, and subsequently, the medium was replaced with high glucose DMEM containing 10ug/mL insulin and 10% FBS for 2 days, then, the medium was changed every day with DMEM supplemented with 10% FBS until 3T3-L1 was induced into mature adipocyte.¹⁴ The gastric cancer lines were cocultured with the supernatant of adipocytes to mimic peritoneal microenvironment.¹⁵

RNA Isolation and Quantitative Real-Time PCR (RT-qPCR)

According to the manufacturer’s instructions, the Trizol Reagent was used to isolate the total RNA from gastric cancer cells. RNA concentration was measured and 1µg RNA was quantified for reverse transcription using the PrimeScript^TM RT reagent kit (Toyobo, Osaka), then RT-qPCR was run on the Bio-Rad IQ5 Real PCR machine (Bio-Rad, USA) in a 96-well plate containing reaction mixture of pre-primer, post-primer, cDNA, SYBR-Green PCR Master Mix (Takara, Osaka). The 2^−ΔΔCt method was used to calculate the relative expression. The primer sequences used in the study are presented in Additional File, Table 2.

Prognostic Analysis

Gastric cancer sequencing data and corresponding clinical information were obtained from the Cancer Genome Atlas (TCGA) dataset (https://portal.gdc.com). Log rank was used to verify KM survival analysis to compare the differences in survival and time-ROC analysis was performed to compare the predictive accuracy and risk score. Univariate and multivariate Cox regression analyses through survival R package were performed to identify the proper terms to build the nomogram. The forest was used to show the P value, HR and 95% CI of each variable through forestplot R package.

Immunity Analysis

According to gastric cancer sequencing data of TCGA, we divided the overall sample into high expression group and low expression group of LINC00924. Immunedeconv R package was used to assess the percentage abundance of tumor-infiltrating immune cells in each sample. TIDE algorithm was applied to predict potential ICB response.

lnCAR

lnCAR (https://lncar.renlab.org/) was a cancer-related microarray data interactive tool and used for lncRNAs differential expression, co-expression network, KEGG pathway and survival analysis.¹⁶

cBioPortal

cBioPortal (https://www.cbioportal.org/), an open-source platform, could be used to analyze and visualize cancer genomics and clinical data.¹⁷ In this study, cBioPortal was used for mRNAs survival analysis.

Statistics Analysis

Results expressed as mean ± standard error (m ± SE) were analyzed by two-tailed Student’s t-test. All experiments were performed in triplicate and statistical significance was defined as p < 0.05 using the GraphPad Prism software (version 6.0, GraphPad Software, USA). Wilcox test was used for immune checkpoint analysis.

Results

Microarray Data of lncRNAs and mRNAs

To distinguish the differential expression of lncRNAs and mRNAs, we performed microarray assay between peritoneal metastasis lesion and primary gastric cancer foci from three patients without any treatment who were diagnosed with gastric cancer peritoneal metastasis. Compared with primary foci, there was a total of 26,793 lncRNAs and 18,985 mRNAs, of which 3433 lncRNAs and 1210 mRNAs up-regulated significantly, 1341 lncRNAs and 3192 mRNAs down-regulated significantly (Figure 1A–D). Among up-regulated lncRNAs, 38 up-regulated and 9 down-regulated lncRNAs with fold change >2 had known structural and functional studies (Figure 1E). According to the sequencing results of mRNAs and lncRNAs, we further used 38 up-regulated and 9 down-regulated lncRNAs to construct CNC network between lncRNAs and mRNAs (Figure 2A–G). To analyze the function of lncRNAs and mRNAs, GO analysis were used. The GO analysis results of up-regulated lncRNA-related mRNAs showed that most mRNAs were located at cytosol and associated with positive regulation of transcription from RNA polymerase II promoter, positive regulation of gene expression, regulation of apoptotic process, response to drug and negative regulation of gene expression during the biological process. Furthermore, the molecular function of most mRNAs was protein binding, protein homodimerization activity, ubiquitin protein ligase binding, protease binding and methyltransferase activity (Figure 3A–C). The GO analysis results of down-regulated lncRNA-related mRNAs demonstrated that most mRNAs were in cytoplasm and had protease binding and protein homodimerization activity function. The top six enriched GO terms of the biological process in down-regulated lncRNA-related transcripts were positive regulation of transcription from RNA polymerase II promoter, positive regulation of gene expression, response to drug, positive regulation of sequence-specific DNA binding transcription factor activity, positive regulation of I-kappaB kinase/NF-kappaB signaling and extracellular matrix organization (Figure 3D–F).

Figure 1 The expression profile comparison of lncRNAs and mRNAs between peritoneal metastasis lesion and primary foci. (A) The lncRNAs expression variation between peritoneal metastasis lesion and primary foci. (B) The mRNAs expression variation between peritoneal metastasis lesion and primary foci. The red and green dots in (A) and (B) indicated the expression levels of lncRNAs and mRNAs with fold change >2. (C) The heat map of lncRNAs. (D) The heat map of mRNAs. (E) The heat map of lncRNAs with fold change >2. The lncRNAs were selected according to whether their biological function was known.

Figure 2 The coding-non-coding gene co-expression network between lncRNAs and mRNAs. (A–D) The up-regulated lncRNAs and their related mRNAs. (E) The heat map of up-regulated lncRNA-related mRNAs. (F) The down-regulated lncRNAs and their related mRNAs. (G) The heat map of down-regulated lncRNA-related mRNAs. The circle represented mRNA and the triangle represented lncRNA. Red meant up-regulated gene and green mean down-regulated gene.

Figure 3 Gene ontology (GO) enrichment analysis. (A–C) The biological process, cellular component and molecular function of up-regulated lncRNA-related mRNA. (D–F) The biological process, cellular component and molecular function of down-regulated lncRNA-related mRNA.

LINC00924 Was Upregulated Significantly in Simulated Peritoneal Microenvironment

Due to the diversity of cell types in the microenvironment, microarray assay of peritoneal metastasis lesion and primary gastric cancer foci did not represent the expression level of tumor genes accurately. Adipocytes, the most abundant cell type of peritoneum, were used to coculture with gastric cancer cell lines to perform in vivo environment simulation. The results of RT-qPCR showed that NR_126348 (LINC01151), NR_027422 (FAM27B) and NR_027133 (LINC00924) had a significant increase in AGS, BGC823 co-culture group (Figure 4A–F). Furthermore, we used lnCAR website to analyze and visualize the expression, co-expression and pathway analysis. The results showed that LINC00924 was up-regulated in gastric cancer and the related genes were associated with metabolic pathways according to GSE65801 dataset (Figure 5A–C). GSE30727 and GSE13911 analysis showed that LINC01151 and FAM27B were also up-regulated and related-mRNAs involved in metabolic pathways (Figure 5D–I); however, LINC01151 and FAM27B had no effect on the survival of patients with gastric cancer surprisingly (data was not shown). Above all, the LINC00924 gene expression profile was peritoneal metastases > primary foci > normal tissue and high LINC00924 expression was closely associated with low overall survival.

Figure 4 LncRNAs expression of gastric cancer cells after coculture with adipocytes. (A–C) The RT-qPCR was used to analyze the lncRNA expression differences of AGS with or without coculture. (D–F) The RT-qPCR was used to analyze the lncRNA expression differences of BGC823 with or without coculture.

Figure 5 The function analysis of LINC00924, LINC01151 and FAM27B by cBioPortal. (A) The transcription levels of LINC00924. (B) The co-expression network of LINC00924. (C) The KEGG pathway analysis of LINC00924-related mRNAs. (D) The transcription levels of FAM27B. (E) The co-expression network of FAM27B. (F) The KEGG pathway analysis of FAM27B-related mRNAs. (G) The transcription levels of LINC01151. (H) The co-expression network of LINC01151. (I) The KEGG pathway analysis of LINC01151-related mRNAs.

LINC00924-Related PEX5L in CNC Network Was an Independent Prognostic Factor and Associated with Poor Survival

According to the CNC network based on our sequencing data, we took the intersection of LINC00924 and FAM27B-related mRNAs. CD28, NHLRC1, CYREN (C7orf49), COLQ, PEX5L, PRMT8, SNX32, HAX1, E2F1, AADACL4, PITX1, ATG10, ELOF1, RGR, SYP, ISLR2, CDS1, TAZ and RP1 were in the intersection. Furthermore, we studied the effects of lncRNA-related mRNAs, age, gender and p-TNM stage on prognosis of gastric cancer. In univariate Cox regression analysis, PEX5L, AADACL-4 and TNM stage were high-risk factors. The hazard ratio (95% CI) of PEX5L, AADACL-4 and TNM stage was 1.91265 (1.30176,2.81020), 2.60071 (1.34318, 5.03559) and 1.27432 (1.06347, 1.52698), respectively (Figure 6A–D). For multivariate Cox regression analysis, the hazard ratio (95% CI) of PEX5L, AADACL-4 and TNM stage was 1.91455 (1.24433, 2.94577), 2.56854 (1.32327, 4.98567) and 1.40575 (1.16321, 1.69886), respectively (Figure 6E–H). Furthermore, we used stomach adenocarcinoma (TCGA, Firehose Legacy) clinical and genomic data in cBioPortal website and found that the patients with high PEX5L expression had a poor survival (Figure 6I–L). The above results demonstrated that LINC00924-related PEX5L was an independent prognostic factor in gastric cancer progression.

Figure 6 The univariate, multivariate analysis and survival curve of lncRNA-related mRNA derived from sequencing data of gastric cancer peritoneal metastasis in gastric cancer. (A–D) The univariate analysis. (E–H) The multivariate analysis. (I–L) Overall survival curves and different colors represented different genes.

LINC00924 Was Associated with a Worse Response Rate to ICB Treatment

According to the median of expression, we divided 375 gastric cancer samples into two groups: high LINC00924 expression group and low LINC00924 expression group. To explore the relationship between LINC00924 and immunity, we used EPIC algorithms to make immune infiltration estimations. In high LINC00924 expression group, the percentage abundance of tumor-infiltrating B cells, CD4+ T cells, endothelial cells, macrophages and NK cells were up-regulated (Figure 7A). CD274, CTLA4, HAVCR2, LAG3, PDCD1, PDCD1LG2, TIGIT and SIGLEC15 were immune checkpoints. In high LINC00924 expression group, CD274, CTLA4, HAVCR2, PDCD1, PDCD1LG2 and TIGIT were up-regulated significantly (Figure 7B), so we further used TIDE algorithm to analyze ICB response. In our study, high LINC00924 expression was associated with high TIDE score; moreover, there were 48 responders in high expression group and 92 responders in low expression group (Figure 7C). To further explored the effect of LINC00924 on gastric cancer prognosis, we plotted a scatter plot of gene expression and the corresponding survival time in different samples (Figure 8A). The prognosis of LINC00924 gene was significant (Log rank P =0.0989), and the higher expression was related to the worse prognosis (Figure 8B). Time-ROC analysis showed that 1-Year, AUC=0.584, 3-Year, AUC=0.594 and 5-Year, AUC=0.676 (Figure 8C). All results showed that LINC00924 was associated with a worse response rate to ICB treatment and LINC00924 could be used to predict the prognosis of gastric cancer.

Figure 7 The relationship between LINC00924 expression and immune escape. (A) The percentage abundance of tumor-infiltrating immune cells in LINC00924 high and low expression group. (B) Immune checkpoint-related gene expression in LINC00924 high and low expression group. (C) Sample statistical table and distribution of immune response scores in LINC00924 high and low expression group in the prediction results. *P < 0.05; **P < 0.01; ***P < 0.001.

Figure 8 The prognostic analysis of LIN00924 in the TCGA set. (A) The curve of risk score, survival status of the patients and heat map of LINC00924 expression were shown. (B) Kaplan–Meier survival analysis of LINC00924. (C) Time-dependent ROC analysis of LINC00924.

Discussion

In this study, we found the differential lncRNAs and mRNAs of peritoneal metastases foci versus primary foci through microarray sequencing. LINC00924 played an important role in gastric cancer progression. LINC00924 was associated with poor survival and immune escape. LINC00924-related PEX5L was an independent prognostic factor. Our study was helpful to provide the basis for molecular targeted therapy of gastric cancer patients with peritoneal metastasis.

Adipocytes was major cell components of peritoneum and adipocyte-derived cytokines and fatty acids promoted cancer progression in varieties of molecular mechanisms. On the one hand, FABP4 (fatty acid binding protein 4)¹⁸ and PITPNC1 (Phosphatidylinositol transfer protein, cytoplasmic 1)¹⁹ were up-regulated after induction of adipocytes which accelerated fatty acids transport to satisfy energy demands of rapid tumor growth; On the other hand, adipokines derived from adipocytes activated the intracellular signaling pathway to promote oncogenes expression to increase the malignant potential.²⁰ In our results, lncRNAs were up-regulated after coculture with adipocytes, but restricted by many factors, we could not find the upstream molecular of these lncRNAs secreted by adipocytes.

Gastric cancer peritoneal metastasis was a complex biological process concerned with abnormal expression of many coding and non-coding RNAs. In our study, 38 lncRNAs were upregulated. Xiaorong Lin’s study showed that lncRNA BDNF-AS (NR_033312) promoted RNH1 ubiquitination degradation in TRIM21-dependent manner to induce poor outcomes and endocrine therapy resistance of breast cancer, thereby leading to malignant progression.²¹ FGD5-AS1 (NR_046252) enhanced non-small cell lung cancer cells’ proliferation through regulating FGD5-AS1/has-miR-107/FGFRL1 axis,²² and Fu et al’s study also demonstrated that FGD5-AS1 improved cisplatin resistance of non-small cell lung cancer cells through miR-140-5p/WEE1 signaling pathway.²³ Moreover, other malignant tumors have also been reported to correlate with high FGD5-AS1 expression. FGD5-AS1/miR-383/SP1 axis promoted esophageal squamous cell carcinoma progression and FGD5-AS1/ miR-153-3p-CITED2 enhanced gastric cancer cell malignant characteristics both in a competing endogenous RNA mechanism.^24,25 Our results showed that many lncRNA-related mRNAs participated in positive regulation of gene expression. PITX1 (paired like homeodomain 1), IKZF1 (IKAROS family zinc finger 1), CREB3L4 (cAMP-responsive element-binding protein 3 like 4) and E2F1 (E2F transcription factor 1) encoded a transcription factor to function as a regulator of gene. PITX1 as a known tumor suppressor gene cis-regulated PDCD5 expression to promote tumor cell apoptosis in gastric cancer.²⁶ IKZF1 was a hub gene of pediatric acute lymphoblastic leukemia and its effect on cancer was still obscure.²⁷ High CREB3L4 expression inhibited preadipocytes 3T3-L1 differentiation into mature adipocytes by reducing expression of C/EBPα and PPARγ2,²⁸ so we speculated that low CREB3L4 expression in peritoneal metastasis lesion was associated with omentum adipogenesis and energy supply of cancer cells. E2F1-a tumor promoter gene-derived breast cancer progression by regulating FGF13.²⁹ Summarily, lncRNAs and lncRNA-related mRNAs promoted tumor development through metabolic programming change induced by microenvironmental remodeling and proto-oncogene activation.

Tumor immune escape referred to the phenomenon that tumor cells survived and proliferated in vivo by escaping from the recognition and attack of the body’s immune system through a variety of mechanisms. The immune system of the body had the function of immune monitoring. When the malignant cells appeared in the body, the immune system could recognize and eliminate these “non-self” cells specifically through the immune mechanism to resist the occurrence and development of tumors. However, malignant cells could evade the immune surveillance of the body through a variety of mechanisms in some cases, and proliferate rapidly in the body to form tumors.³⁰ Cytotoxic T-lymphocyte-associated antigen 4 (CTLA4) was constitutively expressed on T cells. When bound to CD80 or CD86 on the surface of antigen-presenting cells, it acted as a “off” switch to downregulate the immune response.³¹ In LINC00924 high expression group, the immune checkpoints were also up-regulated which meant the immunosuppressive status of gastric cancer. CD4+T cells not only killed tumor cells directly through IFN-gamma-dependent mechanisms but also played an important role in the generation of cytotoxic T lymphocytes (CTL).³² However, the ratio of CD4+ T cells was upregulated in LINC00924 high expression group in our results. The abnormal results might be due to the fact that tumor was a heterogeneous mass and cancer progression was the result of competition between antitumor cells and suppressor cells.³³

Conclusion

In conclusion, LINC00924 was associated with poor survival, immune escape and the related-mRNA, PEX5L was a high-risk factor for gastric cancer progression. LINC00924 might be the therapeutic targets for gastric cancer.