Screening and validation of platelet activation-related lncRNAs as potential biomarkers for prognosis and immunotherapy in gastric cancer patients.

Background: Platelets (PLT) have a significant effect in promoting cancer progression and hematogenous metastasis. However, the effect of platelet activation-related lncRNAs (PLT-related lncRNAs) in gastric cancer (GC) is still poorly understood. In this study, we screened and validated PLT-related lncRNAs as potential biomarkers for prognosis and immunotherapy in GC patients.
Methods: We obtained relevant datasets from the Cancer Genome Atlas (TCGA) and Gene Ontology (GO) Resource Database. Pearson correlation analysis was used to identify PLT-related lncRNAs. By using the univariate, least absolute shrinkage and selection operator (LASSO) Cox regression analyses, we constructed the PLT-related lncRNAs model. Kaplan-Meier survival analysis, univariate, multivariate Cox regression analysis, and nomogram were used to verify the model. The Gene Set Enrichment Analysis (GSEA), drug screening, tumor immune microenvironment analysis, epithelial-mesenchymal transition (EMT), and DNA methylation regulators correlation analysis were performed in the high- and low-risk groups. Patients were regrouped based on the risk model, and candidate compounds and immunotherapeutic responses aimed at GC subgroups were also identified. The expression of seven PLT-related lncRNAs was validated in clinical medical samples using quantitative reverse transcription-polymerase chain reaction (qRT-PCR).
Results: In this study, a risk prediction model was established using seven PLT-related lncRNAs -(AL355574.1, LINC01697, AC002401.4, AC129507.1, AL513123.1, LINC01094, and AL356417.2), whose expression were validated in GC patients. Kaplan-Meier survival analysis, the receiver operating characteristic (ROC) curve analysis, univariate, multivariate Cox regression analysis verified the accuracy of the model. We screened multiple targeted drugs for the high-risk patients. Patients in the high-risk group had a poorer prognosis since low infiltration of immune killer cells, activation of immunosuppressive pathways, and poor response to immunotherapy. In addition, we revealed a close relationship between risk scores and EMT and DNA methylation regulators. The nomogram based on risk score suggested a good ability to predict prognosis and high clinical benefits.
Conclusion: Our findings provide new insights into how PLT-related lncRNAs biomarkers affect prognosis and immunotherapy. Also, these lncRNAs may become potential biomarkers and therapeutic targets for GC patients.

Introduction

According to the most recent statistics from the American Cancer Society, the quantity of new cases and deaths cases of gastric cancer (GC) still remain a high level, and GC is the most common malignant tumor of the digestive system (Cao et al., 2020; Siegel et al., 2021). Although the 5-year survival rate of patients with early GC can reach more than 90%, due to the lack of effective biomarkers and specific clinical appearances, GC patients often present with an advanced-stage tumor at the time of diagnosis (Song et al., 2017), losing their chance to undergo surgery (Tan, 2019). Therefore, the search for new cancer-related prognostic molecular biomarkers and new targets is still needed to enhance the individual evaluation and survival rate of GC.

Studies have indicated that platelet (PLT) regulates tumorigenesis and tumor progression, such as GC, prostate cancer, lung cancer, breast cancer, and colorectal cancer etc (Oh et al., 2019; Garmi et al., 2020; Meikle et al., 2020; Plantureux et al., 2020; Rudzinski et al., 2020; Singla et al., 2020; Wang et al., 2020). For examples, PLT directly promote epithelial-mesenchymal transformation (EMT) of malignant tumors by producing TGF-β, leading to poor prognosis (Labelle et al., 2011; Heldin et al., 2012; Guo et al., 2019; Chong et al., 2020). In tumor angiogenesis, PLT can also produce vascular endothelial growth factor (VEGFR) and promote angiogenesis, providing oxygen and nutrients to malignant cells (Sabrkhany et al., 2011). PLT-produced particles contain a large number of bioactive substances that promote tumor progression. On the contrary, tumor cells can also produce a variety of bioactive substances, such as ADP, TL-6, and TGF-α, to promote PLT activation, and the activated PLT gather around tumor cells and promote tumor progression and metastasis, resulting in a vicious cycle (Schlesinger, 2018). In addition, Zaslavsky et al. found that PLT-generated PD-L1 can induce tumor cells that do not express PD-L1 to avoid being cleared by T cells and evade immune surveillance, thus leading to the progression of malignant tumors (Zaslavsky et al., 2020). All this evidence indicates that PLT have prognostic and immunotherapeutic values.

Long non-coding RNAs (lncRNAs) are about 200 nt or more non-coding protein RNAs, which significantly affect tumor immunity (Chandra Gupta and Nandan Tripathi, 2017). Recently, Ye et al. suggested using circulating lncRNAs between tumor-educated platelets (TEPs) and serum can be used as a potential diagnostic and discriminative biomarkers for colorectal cancer (Ye et al., 2022). Bioinformatics research indicated that the dysregulation of PLT-related genes is involved in cancer (Xie et al., 2021). Yet, the specific effect of platelet activation-related lncRNAs (PLT-related lncRNAs) is still unclear. Exploring the effect and mechanism of PLT-related lncRNAs in the development and progression of GC may help predict prognosis and therapy targets.

In this study, we first extracted 14,087 lncRNAs expression matrix of GC patients from the Cancer Genome Atlas (TCGA) database, and ninety-four genes related to PLT activation were extracted from the Gene Ontology (GO) Resource Database. Then, bioinformatics analysis was performed to identify PLT-related lncRNAs using Pearson’s correlation analysis, after which prognostic risk models were established, and related signaling pathways were screened. Then, we screened for candidate drugs through the Connectivity Map (CMap) database. In addition, we explored the relationship between EMT markers, DNA methylation regulators, and immunotherapy responses and the risk model. Finally, we constructed a nomogram that can predict the overall survival (OS) of GC patients. The study workflow showed in Figure 1.

FIGURE 1

The main process of this study.

Materials and methods

Data and samples collection

A total of 417 cases (375 cases of gastric cancer group and 32 cases of normal tissue) with clinical data and RNA sequencing dataset were downloaded from the Cancer Genome Atlas (TCGA). lncRNAs and mRNAs were classified by the Ensemble Human Genome browser GRCh38.p13. Ten gastric cancer and adjacent tissue specimens (specimen collection time: June 2021 to December 2021) were additionally collected from Jinan Central Hospital affiliated to Shandong First Medical University, which had been approved by the Ethics Committee of Jinan Central Hospital affiliated to Shandong First Medical University, and all patients had signed informed consent. Supplementary Table S1 shows the clinicopathological characteristics of the included patients. Upon collection, fresh tumor and adjacent normal tissues were snap frozen and stored at −80°C until they were taken out. The data from TCGA is public and therefore does not require ethical approval from the relevant authorities.

Identification of PLT activation-related lncRNAs

A total of 94 genes related to PLT activation were collected from the Gene Ontology (GO) Resource database (http://geneontology.org/). For the purpose of evaluating the relationship between PLT activation genes and lncRNAs, Pearson correlation analysis was conducted with R software (R 4.2.1), and the intersection of lncRNA expression in GC patients with a correlation coefficient of 0.4 and p value <0.001 was obtained. A total of 848 lncRNAs associated with PLT activation and their co-expression networks were obtained using “limma” package (Supplementary Figure S1) (Ritchie et al., 2015).

Construction and validation of PLT-related risk model

We combined lncRNAs expression with survival data using “limma” packages to obtain a prognostic lncRNAs expression matrix associated with PLT activation (p < 0.05). Using “survival” package and pfilter = 0.05, univariate Cox analysis showed that 18 PLT-related lncRNAs were significantly correlated with OS (Simon et al., 2011). The “ggpubr” package was then used for differential analysis to obtain the related heatmaps of lncRNAs expression levels in normal and tumor tissues (Whitehead et al., 2019). Lasso regression was performed on these prognostic lncRNAs, and seven lncRNAs associated with PLT in GC were extracted to construct a prognostic risk model (Simon et al., 2011). After excluding normal patients and patients with incomplete clinical data, 372 patients with GC were randomly divided into a testing group and a training group. We used the following algorithm to calculate the risk score for each patient: where coef (lncRNA) represents the prognostic lncRNAs coefficient, and expression (lncRNA) indicates the expression level of lncRNAs (Huang et al., 2021a). GC patients were divided into high- and low-risk groups based on the median risk score. Kaplan-Meier survival analysis used “survival” and “survminer” R package (The R Foundation for Statistical Computing, Vienna, Austria) to estimate the survival difference between the two groups. Then we used the receiver operating characteristic (ROC) curves to evaluate the accuracy of the model (Kim and Hwang, 2020).

Drug screening in risk model

Based on risk scores, effective medicine was screened using CMap (http://portals.broadinstitute.org/camp/) to obtain drugs that reduce risk in high-risk patients (Subramanian et al., 2017). Enrichment score >0 indicated that drugs couldn promote the expression of high-risk genes; a score <0, showed that drugs could suppress the expression of high-risk genes, and p < 0.05 showed that drugs could be significantly enriched (Gns et al., 2019). PubChem website (https://pubchem.ncbi.nim.nih.gov/) was used to obtain the molecular structure of the effective drugs (Kim et al., 2021).

Gene set enrichment analysis

To reveal Gene expression data by sequencing the degree of difference between genes in two groups of samples by using the Gene Set Enrichment Analysis (GSEA) (Subramanian et al., 2005). GC patients were divided into high- and low-risk groups based on the median risk score. For studying the differences in biological functions between risk groups, the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis was conducted with GSEA software, and the pathways enriched by high- and low-risk genes were obtained, respectively. FDR < 0.25 or p < 0.05 were considered statistically significant.

Estimation of the tumor microenvironment using the PLT-related lncRNAs model

Since GSEA results are mostly immune-related, we planned to analyze the tumor microenvironment (TME) in risk model we constructed. CIBERSORT was used to count the immune infiltration statuses of GC patients (Chen et al., 2018). Differences in the content of 22 types of immune cells of high- and low-risk groups were analyzed by the “vioplot” R package (Hu, 2020). Then “ggpubr” package was used to analyze the differences in the TME scores (Estimate-Scores, Immunity-Scores, and Stromal-Scores) of patients in different risk groups, and patients with high TME scores have poorer prognosis (Yoshihara et al., 2013). Exploring the immunotherapy for the model’s applicability can promoting more effective immunotherapy strategies. Then, we analyzed the microsatellite instability (MSI) status (MSS, MSI-H, and MSI-L) of the high- and low-risk groups. MSI status files are from TCIA (http://tcia.at). Also, TIDE (http://tide.dfc-i.harvard.edu/) algorithm was used to assess the different responses to immune checkpoint inhibitors in high- and low-risk groups. When the TIDE score increased, tumors were more likely to have immune escape (Jiang et al., 2018).

Acquisition of DNA methylated regulators and EMT markers

A 5-methylcytosine (5mC) methylated regulator was used to assess the correlation between risk models and DNA methylation. Eleven methylated tuning nodes were obtained from the literature (Chen et al., 2020). EMT-associated genes were used to evaluate the relationship between EMT and the risk model. EMT-related genes were obtained from the EMTome website (Vasaikar et al., 2021). We selected the top 10 markers from the EMTome website for which we could find the expression level for correlation analysis.

Quantitative reverse transcription-polymerase chain reaction analysis

Total RNA was extracted from 10 gastric cancer patients by TRIzon method. cDNA synthesis was then performed using reverse transcription reagents. Quantitative reverse transcription-polymerase chain reaction (qRT-PCR) was performed using 2× SYBR Green HS Premix (AG) on Roche 480 instrument with β-actin as an internal reference. Gene expression levels were calculated using the 2-ΔΔCT method (Livak and Schmittgen, 2001). Supplementary Table S2 shows the primer sequences used to amplify the seven lncRNAs.

Statistical analysis

R software 4.1.2 and GraphPad Prism 8 were used to analyze the data. The R software package “survival” and “survminer” were used for univariate Cox proportional risk regression analysis, multivariate Cox proportional risk regression analysis, and nomogram analysis. The Wilcoxon rank sum test or Kruskal–Wallis rank sum test was used to analyze the differences between the two groups; logarithmic rank testing was used to calculate the statistical difference of OS between high-and low-risk group. The R software package “glmnet” was used for lasso Cox proportional regression, and the R package “survival ROC” was used as the ROC curve (Simon et al., 2011). The p value <0.05 was considered to be statistically significant.

Results

Identification of PLT-related lncRNA of GC patients

Figure 1 shows the detailed workflow of the study. Firstly, 94 PLT activation-related genes were extracted from the Gene Ontology (GO) Resource Database (Supplementary Table S3), and 14,086 lncRNAs expression matrices were extracted from GC from the TCGA database. PLT-related lncRNAs were defined as those that were significantly correlated with one of the 94 PLT-related genes (|PearsonR| > 0.4 and p < 0.001). PLT-related genes and lncRNAs co-expression network is shown in Supplementary Figure S1, and 848 PLT-related lncRNAs were identified (Supplementary Table S3). Eighteen PLT-related lncRNAs were significantly correlated with OS by using univariate Cox regression analysis (Figure 2A). Then, we analyzed the expression levels of these lncRNAs in GC and corresponding normal tissues (Figure 2B). The results showed that among the 18 PLT-related lncRNAs, most lncRNAs (AL355574.1, AC037198.1, LINC01094, LINC02773, LINC00592, AL139147.1, AC002401.4, AL356417.2, AC245041.1, LINC02657, AL139147.1, AC002401.4, AL356417.2, LINC02657, LINC02657, AL355574.1, AC037198.1, LINC01094, LINC02773, LINC00592, AL139147.1, AC002401.4, LINC01711, LINC01614, AL513123.1) were up-regulated and four lncRNAs (AL161785.1, LINC01697, AC129507.1, AP001528.1, AC005165.1) were down-regulated in GC compared to the normal tissues (Figure 2B, p < 0.05). LASSO-penalized Cox analysis was then performed on the 18 lncRNAs, and vertical dashed lines were drawn at the optimal value when the order of Log(λ) was the least likely deviation for OS-related adjustment parameters, and seven lncRNAs related to the prognosis of PLT activation in GC were extracted (Figures 2C,D). These seven PLT-related lncRNAs (AL355574.1, LINC01697, AC002401.4, AC129507.1, AL513123.1, LINC01094, AL356417.2) were used to build a risk model to evaluate the prognostic risk of GC patients.

FIGURE 2

Identification of PLT-related lncRNAs in patients with GC. (A) Univariate Cox regression analysis was used to extract the prognostic lncRNAs. (B) Heatmaps of 18 prognostic lncRNAs expression of patients (***p < 0.001 * *p < 0.01 *p < 0.05). (C) The LASSO coefficient profile of 18 PLT-related lncRNAs. (D) The LASSO coefficient distributions of OS-related lncRNAs and vertical dashed lines were plotted with the values selected for 10x cross-validation. (E) The results of qRT-PCR of PLT-related lncRNAs of 10 pairs GC patients.

In addition, we further verified model-related seven lncRNAs (AL355574.1, LINC01697, AC002401.4, AC129507.1, AL513123.1, LINC01094, AL356417.2) in GC patient tissues and corresponding adjacent tissues using qRT-PCR. We observed that the expression levels of AL355574.1, AC002401.4, LINC01094, and AL356417.2 were up-regulated in most GC tissues, while LINC01697 and AC129507.1 were down-regulated (Figure 2E), which is consistent to the results of TCGA data.

Construction and validation of risk model based on PLT-related lncRNAs

For further testing, the predictive value of the model and the risk scores for each patient were calculated by using a unified formula. Patients were divided into testing set and training set for analysis and validation. Then, based on the median risk score, patients were divided into high- and low-risk groups. The distribution of clinical characteristics of patients in each group is shown in Table 1. The distribution of PLT-related lncRNAs risk scores in the training set and testing set are shown in Figures 3A,B. There were significant differences in the living conditions in survival status among different risk groups. Red dots indicate death and green dots indicate survival. Many cases died in the high-risk group, while most patients in the low-risk group survived (Figures 3D,E). Heatmaps showed seven prognostic lncRNAs expressions for each patient (Figures 3G,H). Figure 3C depicts the distribution of risk levels across all samples for the entire set. The survival status and duration of patients in the entire set are shown in Figure 3F. The prognostic value expression criteria for seven PLT-related lncRNAs risk patterns per patient in the risk model are shown in Figure 3I. Survival analysis of the training set and testing set showed that the high-risk group had a significantly lower survival rate than the low-risk group (Figures 3J,K, p < 0.05). However, the survival analysis of the entire set showed the same results (Figure 3L). ROC curve analysis was used to assess the accuracy of the prognostic model. The results showed that the area under the ROC curve (AUC) of the training set was 0.716 (Figure 3M), the AUC of the testing set was 0.655 (Figure 3N), and the entire set was 0.665 (Figure 3O). The ROC analysis results suggested that the risk model we constructed has high reliability (AUC>0.5). Collectively, these results suggested the good performance of the risk model for survival prediction.

TABLE 1

Distribution of patients’ characteristics.

	Entire set		Train set		Test set
Characteristics	Number	Percentage	Number	Percentage	Number	Percentage
Age
<60 years	111	29.92	57	30.48	54	29.35
≥60 years	257	69.27	128	68.45	129	70.11
Not available	3	0.81	2	1.07	1	0.54
Gender
Female	133	35.85	66	35.29	67	36.41
Male	238	64.15	121	64.71	117	63.59
Grade
Grade 1–2	144	38.81	71	37.96	73	39.67
Grade 3	218	58.76	115	61.5	103	55.98
Not available	9	2.43	1	0.54	8	4.35
Stage
Stage I–II	161	43.4	87	46.52	74	40.22
Stage III–IV	187	50.4	91	48.66	96	52.17
Not available	23	6.2	9	4.82	14	7.61
T
T1–T2	96	25.88	53	28.34	43	23.37
T3–T4	267	72	132	70.59	135	73.37
Not available	8	2.12	2	1.07	6	3.26
M
M0	328	88.41	166	88.77	162	88.04
M1	25	6.74	10	5.35	15	8.15
Not available	18	4.85	11	5.88	7	3.81
N
N0	108	29.11	58	31.02	50	27.17
N1–3	245	66.03	120	64.17	125	67.94
Not available	18	4.86	9	4.81	9	4.89

FIGURE 3

Prognostic value of the risk model of the seven PLT-related lncRNAs in the testing set and training set. (A–C) Distribution of PLT-related lncRNAs model presented based on a training set, testing set, and entire set risk scores. (D–F) Survival time and survival status of low- and high-risk groups for the training set, testing set and entire set. (G–I) Heat-maps of seven LncRNA expressions in the training set, testing set, and entire set. (J–L) Kaplan-Meier survival curves of the OS of patients in the training set, testing set, and entire set. (M–O) ROC curve of the training set, testing set, and entire set.

Then, we conducted univariate and multivariate Cox regression analyses to study whether the prognostic characteristics were independent risk factors. The univariate Cox regression hazard ratio (HR) and 95% confidence interval (CI) of the training set were 2.034 and 1.357–3.049 (p < 0.001); in the testing set, HR was 1.152,95% CI was 1.087–1.221 (Figures 4A,C). HR and 95% CI of multivariate Cox regression in the training set were 2.734 and 1.707–4.374 (p < 0.001) respectively; HR was 1.161, 95% CI was 1.092–1.234 (p < 0.001) in the testing set (Figures 4B,D). For the entire set, we acquired similar results (Figures 4E,F). This result indicated that the risk model was an independent prognostic factor that was not correlated with clinicopathological parameters such as gender, age, tumor grade, and tumor stage.

FIGURE 4

Correlation Analysis between risk score and Clinicopathological Features. (A,B) Univariate- and multivariate-Cox analyses of clinical characteristics and risk score with OS of the training set. (C,D) Univariate- and multivariate-Cox analyses of clinical characteristics and risk score with OS of the testing set. (E,F) Univariate- and multivariate-Cox analyses of clinical characteristics and risk score with OS of the entire set. (G) Kaplan-Meier survival curves of the OS of patients between the risk model and clinical characteristics (age, sex, TNM stage, grade, and survival status). (H) Heat-map of correlation between high- and low-risk and patient clinical characteristics (***p < 0.001 * *p < 0.01 *p < 0.05).

Correlation analysis between risk score and clinicopathological features

Based on the TCGA clinical data, differences in OS stratified by common clinicopathological features were analyzed between the low-risk and high-risk groups. In subgroups divided by gender, age, stage, or tumor stage, the OS of the low-risk group was significantly better than that of the high-risk group (Figure 4G and Supplementary Figure S2A). In addition, OS difference curves were stratified between high-risk and low-risk groups by age, gender, tumor grade, or TNM stage. Risk and clinical correlation heatmap showed that risk score is related to Grade, N (p = 0.0087) and immune score (p < 0.001), but not to age, gender, and TM stage, Stages (p < 0.05) (Figure 4H and Supplementary Figure S2B).

Identification of drugs targeting PLT-related lncRNAs model

In order to determine the effective drug for the PLT-related lncRNAs model, we used the CMap drug screening website (https://portals.broadinstitute.org/cmap/). For enrichment scores, negative values indicate that the drugs can inhibit the expression of high-risk genes and improve the survival rate of patients. Positive values represent that it can promote the expression of high-risk genes (Gns et al., 2019). Seventy compounds were screened out (p < 0.05). All screened compounds could reduce the death risk in high-risk patients, and thus deserve further analysis in GC patients (Table 2). The secondary structure and tertiary structure of some drugs are shown in Supplementary Figure S3.

TABLE 2

The compounds screened that can reduce GC patients’ risk.

Rank	Cmap name	Mean	n	Enrichment	p
1	heptaminol	−0.331	5	−0.842	0.00026
2	etiocholanolone	−0.36	6	−0.738	0.00068
3	trimethobenzamide	−0.37	5	−0.791	0.00078
4	thapsigargin	−0.61	3	−0.893	0.00236
5	sulfamonomethoxine	−0.273	4	−0.8	0.0031
6	pheneticillin	−0.296	4	−0.777	0.00511
7	amantadine	−0.263	4	−0.767	0.00599
8	colistin	−0.448	4	−0.764	0.00635
9	3-acetamidocoumarin	−0.35	4	−0.762	0.00656
10	alprostadil	−0.231	7	−0.578	0.00976
11	N-acetylmuramic acid	−0.533	4	−0.715	0.01333
12	Prestwick-1103	−0.45	4	−0.701	0.01671
13	pyrithyldione	−0.462	4	−0.677	0.02397
14	vincamine	−0.294	6	−0.564	0.02543
15	Prestwick-857	−0.281	4	−0.661	0.02988
16	terazosin	−0.332	4	−0.655	0.03284
17	aconitine	−0.281	4	−0.648	0.03638
18	indoprofen	−0.431	4	−0.639	0.04048
19	acebutolol	−0.385	5	−0.577	0.04075
20	2-aminobenzenesulfonamide	−0.326	4	−0.631	0.04538
21	nifuroxazide	−0.247	4	−0.63	0.04611
22	nimodipine	−0.187	4	−0.625	0.0487

Pathway enrichment analysis

To further explore the potential molecular mechanism of PLT-related lncRNAs and study the differences in biological functions between risk groups, each clinical sample was divided into high-risk (C2) and low-risk (C1) groups. Then, KEGG pathway enrichment analysis was performed with GSEA software. The pathways enriched by high- and low-risk genes were obtained, respectively. Pathways enriched in the C2 group mainly included complement and coagulation cascades, hematopoietic cell lineage, neuroactive ligand-receptor interaction, ECM receptor interaction, and other signaling pathways (Supplementary Figure S4). The pathways enriched in the C1 group mainly included spliceosome, RNA degradation, RNA polymerase, spliceosome, neuroactive tRNA biosynthesis, base excised repair, nucleotide excised repair, homologous recombination, P53 signaling pathway, et al. (Supplementary Figure S4). Figure 5 shows the top five pathways with the highest correlation in the high- and low-risk group. Details of the GSEA results are listed in Supplementary Tables S4, S5. We found that the high-risk group had more pathways related to immunosuppression, such as extracellular matrix (ECM) receptor interaction, which is a complex network of ECM molecules (Zeltz et al., 2020).

FIGURE 5

Pathway Enrichment Analysis. (A) GSEA analysis of the high- and low-risk groups.

TME and immunotherapy response evaluation using risk model

CIBERSOPT was used to analyze the correlation between TME and tumor immunotherapy in the PLT-related lncRNA model. Next, we analyzed the differences of 22 immune cell subtypes in the high- and low-risk groups (Figure 6A). Lower-risk patients had higher enrichment levels of immune killer cells. For example, B cells naive, Plasma cells, T cells follicular helper T cells regulatory, and Macrophages M0 cells were significantly increased in the low-risk group (Figure 6A). In addition, we also validated the correlation of the risk model with immune cells using other algorithms (Supplementary Figure S5). The results of the TME scores assessment showed that indicated that the immune, stromal, and estimate scores of the high-risk group were higher than the low-risk group (p < 0.05) (Figure 6B). Besides, more and more studies show that microsatellite instability (MSI) status affects the TME and patients with microsatellite instability-high (MSI-H) are more sensitive to immunotherapy (Lin et al., 2020). Our studies showed that low-risk scores were associated with MSI-H, which predicted that low-risk patients are more likely to benefit from immunotherapy (Figures 6C,D). Furthermore, differences in TIDE scores between high- and low-risk groups were obtained, and the results showed high-risk patients had higher TIDE scores, predicting poorer immunotherapy outcomes (Figure 6E).

FIGURE 6

Potential role of the risk model in the TME and immunotherapy. (A) The content of 22 immune cells between the high- and low-risk groups. (B) TME Estimate-Scores, Immunity-Scores, and Stromal-Scores measured between high- and low-risk groups. (C,D) Differences in microsatellite instability (MSI) between patients in the high- and low-risk groups. (E) Differences in immunotherapy scores between high- and low-risk groups.

Correlation analysis of PLT-related lncRNAs with DNA methylation and EMT

DNA methylation and lncRNA regulation are generally considered to be important factors in cell differentiation and development (Tang, 2018). Some studies indicated that DNA methylation at the same locus is associated with PLT activation variability in well-defined populations (Izzi et al., 2019). DNA methylation involved in general research mainly refers to the methylation process that occurs at the 5th-carbon atom of cytosine in CpG dinucleotides, a product also called 5-methylcytosine (5mC), which is the earliest methylation type excavated in eukaryotes (Ye and Li, 2014). As one of the important epigenetic markers, 5mC has a significant effect on various physiological and pathological processes (Ye and Li, 2014). Next, we explored whether there is a link between the risk model and DNA methylation. We analyzed the relationship between these 5mC regulators (DNMT1, DNMT3A, DNMT3B, MBD1, MBD2, MBD3, MBD4, MECP2, NEIL1, NTHL1, SMUG1, TDG, UHRF1, UHRF2, UNG, ZBTB33, ZBTB38, Z BTB4, TET1, TET2, TET3) in high- and low-risk groups (Chen et al., 2020). We discovered that most of the regulators were different between high- and low-risk groups (Figure 7A, p < 0.05), indicating that the risk model we constructed is correlated with DNA methylation. This result suggests that DNA methylation is one of the major biological characteristics of the high-risk group.

FIGURE 7

Correlation analysis of PLT-related lncRNAs with methylation and EMT. (A) Differences of 5mC Regulator expression between patients in the high- and low-risk groups. (B) Differences of EMT-related gene expression levels among high- and low-risk groups.

PLT could produce TGFβ, a cytokine highly related to EMT, which has an extremely important role in EMT(13). Therefore, we tried to explore the correlation between the risk model and EMT markers. EMT-related genes came from EMTome. Ten genes (VIM, CDH1, FN1, ZEB1, CDH2, SNAI2, SPARC, SNAI1, CCN2 and TWIST1) were selected for correlation analysis. We found that all the EMT-related genes we picked were significantly correlated with the risk model and the high-risk group patients had higher EMT gene expression (Figure 7B). The results indicated a strong correlation between the risk model and EMT, which may explain the poor prognosis of high-risk groups.

Construction and assessment of the novel nomogram

We also used 1-year, 2-year, and 3-year calibration charts to prove that the nomogram was in good agreement with the prediction of 1-, 2-, and 3-year OS (Figure 8A). Nomogram including risk grade and clinical risk characteristics were used to predict the incidence of OS at 1-, 2-, and 3-year. The risk level of the prognostic model showed outstanding predictive power in the nomogram compared to clinical factors (Figure 8A). The observed ratios of 1-year, 2-year, and 3-year OS showed definitive agreement with the predicted ratios (Figures 8B–D).

FIGURE 8

Construction and Assessment of the Novel Nomogram. (A) The nomogram that predicted 1 -, 2 -, and 3-year survival probabilities. (B–D) The calibration curve for 1 -, 2 -, and 3-year OS.

Discussion

More and more studies have been conducted on lncRNAs in recent years, which have an important role in cancer progression (Chandra Gupta and Nandan Tripathi, 2017). However, there are few studies on the role of lncRNAs in GC. The study of PLT-related lncRNAs can provide a new direction for exploring GC pathogenesis and targeted therapy. It is particularly important to study the prognostic significance of PLT-related lncRNAs in GC. Over the years, research on the effect of lncRNAs and PLT in tumors has gradually become a hot topic in medical research.

Currently, growing evidence suggests that PLT play an important part in the occurrence and development of GC. Activated PLT can promote thrombus formation, thereby accelerating tumor progression (Suzuki-Inoue, 2019). Therefore, exploring the mechanism of PLT activation in the progression of GC has an important meaning in improving the survival rate of GC patients and improving the effect of immunotherapy. Molecular markers associated with PLT activation may also have an important role in predicting the clinical risk and prognosis of GC patients. Xie et al. proposed a novel PLT-related gene signature as a practical tool for patients with triple-negative breast cancer (TNBC) with independent value in assessing clinical prognosis (Xie et al., 2021). In addition, there are more and more research on the effect of lncRNAs in tumor progression, and they have attracted more and more attention (Li et al., 2016; Wang et al., 2022). A previous study has shown that TEPs derived lncRNAs occupy an important position in the diagnosis and treatment of colorectal cancer and may elucidate the underlying molecular mechanism of PLT-tumor cell interaction, which may be related to circulating lncRNAs in the blood (Ye et al., 2021). Therefore, both PLT and lncRNAs are closely related to the occurrence and development of GC. Still, there is little research on the role of PLT-related lncRNAs in the prognosis of GC patients. Xu et al. established an m6A-related lncRNAs model and confirmed the model’s important effect in predicting the prognosis of patients with lung adenocarcinoma (LUAD), providing guidance for the immunotherapy of patients with LUAD (Xu et al., 2021). Furthermore, a recent study on the effect of autophagy on GC constructed a prognostic model containing five autophagy-related lncRNAs, indicating the key role of autophagy-related lncRNAs in GC and suggesting that these lncRNAs may be effective targets for immunotherapy point (Chen et al., 2021a). In addition, models of ferroptosis-related lncRNAs, necroptosis-related lncRNAs, and pyroptosis-related lncRNAs have also been established, providing new targets for the study of the molecular mechanism, and the immunotherapy of malignant tumors (Song et al., 2021; Xiao et al., 2021; Zhao et al., 2021). In this study, we firstly constructed an independent prognostic model based on PLT-related lncRNAs.

We first extracted and identified 848 PLT -related lncRNAs from the TCGA database and performed a series of analytical validations to explore the value of PLT-related lncRNAs in GC prognosis. We verified the prognostic value of 18 PLT-related lncRNAs in GC by univariate COX regression analysis. By univariate COX regression analysis, we verified the prognostic value of 18 PLT-related lncRNAs in GC. Seven PLT-related lncRNAs (AL355574.1, LINC01697, AC002401.4, AC129507.1, AL513123.1, LINC01094, AL356417.2) were identified by LASSO regression analysis and used to construct the prognostic model for predicting OS in GC patients. AL355574.1 was identified as a protective lncRNA associated with autophagy in GC, which can be used as a promising therapeutic target for immunotherapy in GC patients (Chen et al., 2021a). Zhang et al. pointed out that LINC01697, as a ceRNA, could be used as a biomarker for the prognosis of GC patients and was up-regulated in GC cells, while its knockdown can inhibit the proliferation of GC cells (Zhang et al., 2021). Moreover, Li et al. suggested that LINC01697 as a prognostic biomarker for oral squamous cell carcinoma (Li et al., 2020a; Zha et al., 2021). Zha and others showed AC129507.1 as a DElncRNA was upregulated in GC and significantly associated with the prognosis of GC patients (Zha et al., 2021). Sun et al. found that AL513123.1 was upregulated in a high-risk group and could be used as a DElncRNA closely related to the prognosis of breast cancer (BRCA) (Sun et al., 2019). In addition, Tuersong et al. obtained similar results, pointing out that AL513123.1 in BRCA may be involved in the regulation of the complex ceRNA network and identified as a potential prognostic biomarker and therapeutic target for BRCA diagnosis and treatment (Tuersong et al., 2019). LINC01094 is associated with the prognosis of ovarian cancer, pancreatic cancer, glioma, and renal clear cell carcinoma (Xu et al., 2020a; Jiang et al., 2020; Chen et al., 2021b; Luo et al., 2021; Liu et al., 2022). Li et al. found that AL356417.2, as an immune-related lncRNA in BRCA, is closely associated with the prognosis of BRCA and can be used as a prognostic molecular marker and immunotherapy target for BRCA patients (Li et al., 2020b). In this study we discovered and verified AC002401 for the first time. Therefore, the seven PLT-related lncRNAs obtained in our study may also become the important biomarkers and therapeutic targets for GC and even other cancer types.

We divided patients into high- and low-risk groups according to the risk scores of the model and evaluated the mechanism of regulating GC progression by GSEA. The GSEA results indicated that the complement and coagulation pathway was the most upregulated gene-enriched signaling pathway in a high-risk group. In addition, immune-related pathways were also enriched, such as ECM receptor interactions. The complement system participates in multiple pathological processes such as thrombotic diseases, immune responses, autoimmune diseases, and cancer (Afshar-Kharghan, 2017). Firstly, it is involved in various tumorigenesis and cancer progression stages by mediating inflammatory responses. Secondly, complement activation may have a role by modulating T cell response to a tumor. Markiewski et al. showed that activation of the classical complement pathway promotes in situ tumor growth in mice (Markiewski et al., 2008). The immunomodulatory effect of the classical complement pathway activated in the tumor can promote tumor growth. In addition, we also noted that the ECM receptor interaction pathway was significantly enriched in a high-risk group. The role of ECM has been demonstrated in several cancers. Bao et al. showed that ECM-related proteins or genes might be potential biomarkers for breast cancer diagnosis and treatment (Bao et al., 2019). Studies have also shown that ECM participates in the invasion and metastasis of GC and promotes EMT in colorectal cancer (Rahbari et al., 2016; Yan et al., 2018).

Because the signal pathways enriched in this model are concentrated on immune-related signal pathways, the correlation between high- and low-risk groups and TME was analyzed. It was found that the model was closely related to immune cell infiltration. CIBERSPOT algorithm was then used to calculate the correlation of different immune cell infiltration. We noted that M2 macrophages, monocytes, and dendritic cells resting had a significantly higher expression in a high-risk group, which indicates that high-risk patients have higher immune cell infiltration. As for the low-risk group, we observed more infiltration of B cell naive and T cell follicular helper. In addition, we found higher immune scores, stromal scores, and estimated scores in the high-risk group. This is consistent with the results of previous studies that high immune score, stromal score, and macrophage infiltration are associated with poor prognoses (Deng et al., 2020). NK cell consumption significantly promotes cancer metastasis in mice (Shimaoka et al., 2017). PLT have been found to protect tumor cells from NK cells, and this effect is mainly due to the transfer of PLT-derived MHC CLASS I molecules to tumor cells after the interaction between PLT and tumor cells, which reduces the anti-tumor reactivity of NK cells and thus avoids immune surveillance (Placke et al., 2012). MSI analysis showed that MSI-H patients would have a better immunotherapy prognosis. In conclusion, the immunotherapy response-related prediction marker showed that patients in the high-risk group had a better response to immunotherapy. Based on this analysis, we concluded that the risk model could contribute to identifying reliable molecular biomarkers for the immunotherapy of GC.

DNA methylation modification is the most common covalent modification method. Many recent studies have confirmed the correlation between methylation modification and malignancy (Xu et al., 2020b). 5mC is the only form of DNA methylation found in mammals, and 5mC methylation regulators are associated with tumor proliferation and metastasis (Huang et al., 2021b). Benedetta et al. reported that PLT-endothelial aggregation receptor 1 (PEAR1), driven by DNA methylation, is a marker of PLT activation variability (Izzi et al., 2019). However, our study confirmed a correlation between the risk model and DNA methylation. Our study found that most DNA methylation regulators were differentially expressed in high-risk and low-risk groups. These results showed that there might be an association between our findings and DNA methylation, which reflects an important biological feature of the model. EMT has been shown to play an important role in tumorigenesis, invasion, and metastasis (Nieto et al., 2016). This study found that the EMT markers we selected were differentially expressed in both high- and low-risk groups. The results indicate a correlation between the risk model and EMT, reflecting another important biological feature of the model. These results suggest that DNA methylation and EMT are responsible for the poor prognosis of high-risk patients.

However, further experiments are needed to prove the effect of PLT-related lncRNAs on the prognosis of GC and related molecular mechanisms. The related signaling pathways screened out in this study and the effectiveness of immunotherapy drugs should be further investigated. In this study, we only analyzed and validated the data in the TCGA database. Although we have carried out some experimental verification using the collected specimens, there may still be deviations and deficiencies. Therefore, the risk model we constructed needs more external data for verification. We plan to collect more clinical samples to further validate the value of these lncRNAs in a future study.

Conclusion

In this study, we constructed a model containing seven PLT-related lncRNAs. This study provides new clues for predicting the prognosis of GC patients and may help to elucidate the process and mechanism of PLT-related lncRNAs. In addition, small-molecule drugs were found to target PLT-related lncRNAs, and risk models showed sensitivity in distinguishing GC patients who benefited from immunotherapy. Our study further explored the role of PLT-related lncRNAs in TME, drug screening, and immunotherapy prediction in GC, providing new directions and therapeutic targets for further research and clinical practice.