You are viewing a javascript disabled version of the site. Please enable Javascript for this site to function properly.
Go to headerGo to navigationGo to searchGo to contentsGo to footer
In content section. Select this link to jump to navigation

Construct dysregulated miRNA-mRNA interaction networks to conjecture possible pathogenesis for Stomach adenocarcinomas



Post-transcriptional regulation of mRNA induced by microRNA is known crucial in tumor occurrence, progression, and metastasis. This study aims at identifying significant miRNA-mRNA axes for stomach adenocarcinomas (STAD).


RNA expression profiles were collected from The Cancer Genome Atlas (TCGA) and GEO database for screening differently expressed RNAs and miRNAs (DE-miRNAs/DE-mRNAs). Functional enrichment analysis was conducted with Hiplot and DAVID-mirPath. Connectivity MAP was applied in compounds prediction. MiRNA-mRNA axes were forecasted by TarBase and MiRTarBase. Real-time reverse transcription polymerase chain reaction (RT-qPCR) of stomach specimen verified these miRNA-mRNA pairs. Diagnosis efficacy of miRNA-mRNA interactions was measured by Receiver operation characteristic curve and Decision Curve Analysis. Clinical and survival analysis were also carried out. CIBERSORT and ESTIMATE was employed for immune microenvironment measurement.


Totally 228 DE-mRNAs (105 upregulated and 123 downregulated) and 38 DE-miRNAs (22 upregulated and 16 downregulated) were considered significant. TarBase and MiRTarBase identified 18 miRNA-mRNA pairs, 12 of which were verified in RT-qPCR. The network of miR-301a-3p/ELL2 and miR-1-3p/ANXA2 were established and verified in external validation. The model containing all 4 signatures showed better diagnosis ability. Via interacting with M0 macrophage and resting mast cell, these miRNA-mRNA axes may influence tumor microenvironment.


This study established a miRNA-mRNA network via bioinformatic analysis and experiment validation for STAD.


Gastric cancer (GC) is the fourth common malignancy and is becoming one of the leading causes of tumor-related death, bringing about approximately 1.1 million new cases and 768 thousand deaths in 2020 globally [1]. More than half of the newly diagnosed cases were observed in developing countries [2]. Helicobacter pylori infection, excessive salt intake, tobacco smoking, alcohol exposure, familial predisposition and lack of fruits and vegetables in daily diets are the main risk factors for gastric cancer [3]. Despite the progression in early diagnosis and therapies of gastric cancer, the overall 5-year survival rate remained low, ranging from 18–58% in different races [4]. Stomach adenocarcinomas (STAD) accounted for almost 95% of gastric cancer [5]. So, it is of great value to dig in candidate therapeutic targets and novel biomarkers with prognostic efficacy for STAD.

MicroRNAs (MiRNAs) are a family of 22 nucleotide short noncoding RNAs known to play an essential role in tumor onset, development, immune escape and metastasis [6]. MiRNA deregulation in tumor cells has been verified in tumor cells [7]. By binding with the 3’-untranslated regions (3’-UTR) of mRNA, miRNAs affect the post-transcriptional regulation of messenger RNA (mRNA) [8]. Numerous studies have confirmed the significance of miRNA-mRNA axes in digestive system neoplasm. For instance, miRNA-223 targeting at EPB41L3 promotes invasion and metastasis of gastric cancer [9]; miR-221-3p and miR-15b-5p targeting at Axin2 modulate cell proliferation and invasion of liver cancer [10]. Thus, investigating the negatively regulated miRNA-mRNA pairs helps to conjecture molecular mechanism of STAD. Besides, miRNAs and their target mRNA can serve as biomarkers in diagnosis, prognosis measurement and speculation of sensitivity to treatment for STAD [11]. Pang et al. suggested that miR-15a-5p targeting at PHLPP2 promotes platinum resistance and that lower serum miR-15a-5p level usually indicates better response to oxaliplatin in gastric cancer [12]. Therefore, miRNAs functioning as potential biomarkers for STAD deserve further research.

Nowadays, databases based on high-throughput sequencing and microarrays provide opportunities for obtaining more global views of miRNAs and mRNAs in tumor. Some studies have employed some Gene Expression Omnibus (GEO) datasets to identify differently expressed miRNAs (DE-miRNAs) and differently expressed mRNAs (DE-mRNAs) [13, 14]. Researches of Zhijie Dong et al. and Yong-jun Guan revealed several key miRNA-mRNA pairs in gastric cancer by analysis of databases and PubMed publications [15, 16]. However, most of these previous researches only absorbed few GEO datasets. In present study, totally 23 mRNA and 8 miRNA GEO datasets concentrating on STAD tissue and normal stomach tissue were incorporated. TarBase and MiRTarBase, databases collecting experiment verified miRNA-mRNA networks, were applied in screening miRNA-mRNA axes. Next, we conducted reverse transcription quantitative polymerase chain reaction (RT-qPCR) in 30 paired STAD tissue and adjacent normal tissue to measure the expression level of these miRNA-mRNA pairs. The receiver operation characteristic (ROC) curve and Decision Curve Analysis (DCA) based on binary logistic regression helped to evaluate the diagnosis efficiency. Clinical subgroup analysis and survival analysis were promoted for prognosis estimation. Correlation analyzes between miRNA-mRNA axes expression and tumor immune microenvironment, immune cell infiltration, tumor mutation burden was carried out to measure their effects on phenotypic features. This study combining bioinformatic analysis and RT-qPCR may throw new a light on the pathogenesis and potential drug target for STAD.

Table 1

Information pertaining to the GEO datasets based on STAD tissues

GEO accessionNumber of tumorsNumber of non-tumorsTotal numberCountryPlatformPMID
miRNA arrayGSE2659560868South KoreaGPL817924222951
GSE54397161733South KoreaGPL1515925167801
mRNA arrayGSE10323610919RomaniaGPL413332587446
GSE33651401252South KoreaGPL289522133303
GSE51575272754South KoreaGPL1360725254613
GSE96668491160the UKGPL1055831050820

Figure 1.

The Flow chart displayed the procedures of identifying miRNA-mRNA axes and extensive analysis for gastric cancer in this study.

The Flow chart displayed the procedures of identifying miRNA-mRNA axes and extensive analysis for gastric cancer in this study.


2.1Data collection from GEO and TCGA databases

To acquire candidate miRNA and mRNA datasets, we searched in the Gene Expression Omnibus (GEO) database with ‘(gastric OR stomach) AND (tumor OR cancer OR adenocarcinoma OR carcinoma) NOT (cell line)’ as the key word. Only microarray datasets measuring miRNA or mRNA expression levels in homo sapiens between normal tissue and tumor tissue were included. Under these constrains, totally 23 mRNA datasets and 8 miRNA datasets (GSE26595; GSE23739; GSE28700; GSE63121; GSE77380; GSE2564; GSE93415; GSE54397; GSE103236; GSE118897; GSE130823; GSE13861; GSE13911; GSE19826; GSE2685; GSE26899; GSE27342; GSE29272; GSE29998; GSE33335; GSE33651; GSE49051; GSE51575; GSE55696; GSE56807; GSE63089; GSE65801; GSE66229; GSE78523;GSE79973; GSE96668) listed in Table 1 were introduced for further analysis. RNA-sequencing profile and clinical information of gastric cancer collected in TCGA database were downloaded from the GDC data portal in the National Cancer Institute (https://portal. Besides, the overlapped differently expressed miRNAs involved in public databases including Human MicroRNA Disease Database (HMDD), Database of Differentially Expressed MiRNAs in human Cancers (dbDEMC) and miRCancer were also absorbed. The flowchart of the whole research is exhibited in Fig. 1.

2.2Data processing of RNA expression profile

GEO2R, an R-based online application using the GEOquery and limma R packages, was employed in differential analysis of GEO datasets. We performed Wilcoxon rank-sun test on sequencing data obtained from TCGA database with an R-based package called psych. Fold change (FC) and P-value were calculated in the two methods mentioned above. P-value < 0.05 and |log2FC|> 0.58 were considered cut-off criterion for identifying differently expressed miRNAs (DE-miRNA) and mRNAs (DE-mRNA). Then, the intersection was taken to obtain the list of DE-miRNAs and DE-mRNAs in the GEO database and TCGA database. By summarizing GEO, TCGA, HMDD, dbDEMC, and miRCancer database, the list of DE-miRNAs and DE-mRNAs was ultimately identified.

2.3Functional enrichment analysis and speculation of agents

Gene oncology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses of DE-miRNAs were performed by DAVID-mirPath GO/KEGG analysis of DE-mRNAs were conducted on Hiplot ( and presented in bubble plot P-value < 0.05 and counts > 2 were considered statistically significant. Connectivity MAP established connection between mRNA expression signatures and approximately 5,000 small molecule compounds, which helps to offer potential agents for gastric cancer treatment. Negative log10 transformed FDR q-values (fdr_q_nlog10) > 2 and normalized connectivity score (norm_cs) <-1 were set as the screening criterion.

2.4Predicting target mRNAs of DE-miRNAs

TarBase and MiRTarBase, databases collecting experimentally verified miRNA-mRNA interactions, were applied in assessing target mRNAs. Only miRNA-mRNA axes recorded in both databases were taken into consideration. To make sure that these negatively regulatory relationship between miRNA and mRNA can be observed in gastric tissue, we searched each miRNA-mRNA pair in PubMed. MiRNA-mRNA networks were plotted with Cytoscape v3.8.2.

2.5Sample collection, RNA isolation and reverse transcription quantitative PCR

We collected 30 paired Formalin-Fixed and Paraffin-Embedded (FFPE) samples of STAD tissue and adjacent normal stomach tissue. All patients had their stomach surgery in the First Affiliated Hospital of Nanjing Medical University None of the patients underwent preoperative treatment. Patients with stage IV gastric cancer do not have indications for surgery, so the samples collected in this study did not include stage IV patients. All patients signed written consent after fully understanding this research. This study was approved by the Institutional Review Boards of the First Affiliated Hospital of Nanjing Medical University (ID: 2016-SRFA-148). The clinical features of patients were disclosed in Table 2.

Table 2

Clinicopathological and molecular features of stomach cancer patients

Stomach cancer (n= 30)Rate (%)
Age (year)
 Mean (SD)66.37 (10.14)
 Median [min, max]68 [48, 86]
Tumor location
 Fundus of stomach310
 Corpus of stomach826.67
 Antrum of stomach620
TNM stage

All stomach FFPE specimen were stored at 20C. We isolated total RNA from FFPR stomach samples with RNAprep Pure FFPE Kit (TIANGEN Biotech, Beijing, China, DP439) according to the manual. Poly (A) Polymerase Kit (Takara) was applied in adding a poly (A) tail to each RNA. Reverse transcription was performed with PrimeScript RT reagent Kit (Takara, Tokyo, Japan, RR037A) at 37C for 15 minutes, then 85C for 5 seconds. Complementary DNA (cDNA) were stored at 4C before conducting polymerase chain reaction (PCR). Quantitative PCR was carried out with SYBR Premix Ex Taq II (Takara, Tokyo, Japan, RR820A). We performed qPCR in a volume of 10 μL on the qTOWER384 (Analytik Jena) plate with the temperature set as 95C for 20 s, followed by 40 cycles composed of 10 s at 95C and 20 s at 60C. The mRNA results were normalized with GAPDH and 18S [17]. RNU6B was set as reference for miRNA [18]. We queried Primer Bank ( for the primer sequence listed in Table 3.

Table 3

The sequences of primers for candidate miRNAs and targeted mRNAs

NameForword primer sequences (5’3’)Reverse primer sequences (5’3’)
miRNAmiR-378a-3pACTGGACTTGGAGTCAGAAUniversal reverse primer
miR-26a-5pTTCAAGTAATCCAGGATAUniversal reverse primer
miR-143-3pTGAGATGAAGCACTGTAGCTUniversal reverse primer
miR-29c-3pTAGCACCATTTGAAATCGGTUniversal reverse primer
miR-1-3pTGGAATGTAAAGAAGTATGTUniversal reverse primer
miR-30a-5pTGTAAACATCCTCGACTGGAAUniversal reverse primer
miR-19a-3pTGTGCAAATCTATGCAAAACTUniversal reverse primer
miR-21-5pTAGCTTATCAGACTGATGTTUniversal reverse primer
miR-20a-5pTAAAGTGCTTATAGTGCAGGTUniversal reverse primer
miR-221-3pAGCTACATTGTCTGCTGGGTTTUniversal reverse primer
miR-301a-3pCAGTGCAATAGTATTGTCAAAUniversal reverse primer

Figure 2.

TarBase and MiRTarBase screened totally 18 negatively regulated miRNA-mRNA axes between 38 DE-miRNAs and 228 DE-mRNAs, which were plotted via Cytoscape. Orange nodes represent the upregulated miRNAs/mRNAs in STADs versus NCs, while blue nodes represent the downregulated miRNAs/mRNAs in STADs versus NCs.

TarBase and MiRTarBase screened totally 18 negatively regulated miRNA-mRNA axes between 38 DE-miRNAs and 228 DE-mRNAs, which were plotted via Cytoscape. Orange nodes represent the upregulated miRNAs/mRNAs in STADs versus NCs, while blue nodes represent the downregulated miRNAs/mRNAs in STADs versus NCs.

Table 4

Pearson’s correlation analysis of miRNA-mRNA pairs in paired FFPE STAD samples

miRNA (up)p-valueFC (2-ΔΔCT)mRNA (down)p-valueFC (2-ΔΔCT)Spearman’s correlation
hsa-miR-21-5p< 0.00015.2087LIFR0.00070.6150.10600.2107
hsa-miR-19a-3p0.06913.4478KIT< 0.00010.61520.4354-0.1026
hsa-miR-20a-5p0.01753.5172KIT< 0.00010.61520.2331-0.1563
hsa-miR-221-3p0.04672.751KIT< 0.00010.61520.2619-0.1472
miRNA (down)mRNA (up)p-valuer-value
hsa-miR-26a-5p0.00050.5768CKS20.30611.7258< 0.00010.6633
hsa-miR-1-3p0.01370.7874ANXA20.00021.6947< 0.0001-0.5775
hsa-miR-30a-5p0.93111.1117TGFBI0.76111.4706< 0.00010.6036

Figure 3.

(A) The expression levels of KIT, LIFR, and ELL2 declined in STAD tissues. ANXA2 was overexpressed in tumor samples. (B) The miRNA expression level of miR-21-5p, miR-20a-5p, miR-221-3p, and miR-301a-3p were increased in STAD specimen. MiR-26a-5p, miR-29c-3p, and miR-1-3p were down regulated in STAD tissues. Data were presented as mean ± SEM. *p < 0.05, **p < 0.01, ***p < 0.001 and ****p < 0.0001 (Wilcoxon rank sum test). (C and D) The networks of miR-1-3p (down)/ANXA2 (up) and miR-301a-3p (up)/ELL2 (down) were the negatively regulated pairs verified via RT-qPCR.

(A) The expression levels of KIT, LIFR, and ELL2 declined in STAD tissues. ANXA2 was overexpressed in tumor samples. (B) The miRNA expression level of miR-21-5p, miR-20a-5p, miR-221-3p, and miR-301a-3p were increased in STAD specimen. MiR-26a-5p, miR-29c-3p, and miR-1-3p were down regulated in STAD tissues. Data were presented as mean ± SEM. *p < 0.05, **p < 0.01, ***p < 0.001 and ****p < 0.0001 (Wilcoxon rank sum test). (C and D) The networks of miR-1-3p (down)/ANXA2 (up) and miR-301a-3p (up)/ELL2 (down) were the negatively regulated pairs verified via RT-qPCR.

Figure 4.

Receiver operating characteristic (ROC) curves of TCGA database and internal validation were separately presented in Fig. 4A and B. (A) The AUC of the complex model containing both miR-1-3p/ANXA2 and miR-301a-3p/ELL2 axes were 0.9208, the maximum among all 15 diagnosis models. (B) DCA presented great diagnostic performance of miRNA-mRNA axes in TCGA database and internal validation.

Receiver operating characteristic (ROC) curves of TCGA database and internal validation were separately presented in Fig. 4A and B. (A) The AUC of the complex model containing both miR-1-3p/ANXA2 and miR-301a-3p/ELL2 axes were 0.9208, the maximum among all 15 diagnosis models. (B) DCA presented great diagnostic performance of miRNA-mRNA axes in TCGA database and internal validation.

2.6Diagnostic efficiency assessment of MiRNA-mRNA axes and Survival analysis for STAD

We conducted binary logistic regression on miRNA and mRNA expression data with SPSS Statistics 26.0 (IBM, New York, USA). The receiver operation characteristic (ROC) curve and Decision Curve Analysis (DCA) were performed for comparing diagnostic efficiency between multiple models. ROC curves were plotted with data downloaded from TCGA database. According to survival data downloaded from TCGA database, we drew Kaplan-Meier survival curve with Hiplot. The average follow-up time was 577 days (1.58 years) with the maximal follow-up time equal to 3720 days (10.19 years). The medium divided all patients into low expression group and high expression group.

2.7Clinical subgroup analysis for STAD

Clinical information of stomach adenocarcinoma patients was downloaded from GDC data portal. According to their age (at first diagnosis), gender, tumor grade and tumor stage, cases were divided into several subgroups. Expression levels of miRNA-mRNA axes were compared between these subgroups by Wilcoxon rank-sum test. P-value < 0.05 was regarded as significant cut-off.

2.8Correlation analysis between miRNA-mRNA axes and tumor phenotypes

Wilcoxon Rank-sum test was applied in comparing the percentage of 22 immune cells between tumor samples and normal samples according to CIBERSORT ( Correlation analysis between proportion of discrepant immune cells and miRNA/mRNA expression level was exhibited in correlation heatmap plotted by Hiplot with Spearman’s correlation. ESTIMATE provided the Stromal Score, immune Score, estimate Score and tumor purity with mRNA expression signatures recorded in TCGA database. Tumor mutational burden (TMB) was measured by Masked Somatic Mutation data obtained with the ‘maftools’, an R-based package. We carried out Spearman’s correlation analysis between miRNA-mRNA axes and the tumor features mentioned above.

2.9Statistically analysis

MiRNA and mRNA expression data obtained from RT-qPCR was managed with 2 - ΔΔCT method (ΔCT = CtmiRNA or mRNA- Ctnormalizer; Ct: the threshold cycle). The correlations between miRNA and mRNA expression data were calculated by Spearman’s correlation via psych package. GraphPad prism v8.0.1 was employed in ROC curves plotting and the area under the ROC curves (AUC) calculating. DCA curves were drawn with an R-based package named ‘rmda’. P-value < 0.05 was considered statically significant. Hiplot, Cytoscape, GraphPad prism and R language were applied in figure construction.


3.1Screening DE-miRNAs and DE-mRNAs in TCGA and GEO for STAD

Differential analysis of TCGA database provided basic DE-mRNAs and DE-miRNAs for STAD. As is presented in Table 1, a total of 23 mRNA GEO datasets and 8 miRNA GEO datasets met our inclusion criterion. Table 1 displayed the GEO accession, number of samples, researchers’ nationality, platform, and PMID of these datasets. 105 upregulated mRNAs and 123 downregulated mRNAs showed statistical significance in after calculating with GEO2R. Compared with miRNA expression profile of normal tissue, 22 upregulated and 16 downregulated miRNAs were considered significant. Combining with dbDEMC, HMDD and miRCancer, 22 upregulated and 16 downregulated miRNAs were selected ultimately.

3.2GO/KEGG analysis and potential compounds prediction

Outcome of function annotation and pathway enrichment analysis based on DE-mRNAs were presented in Fig. S1A and S1B. Mitotic cell cycle phase transition, extracellular matrix structural constituent and collagen-containing extracellular matrix were the GO terms of upregulated mRNAs with the most counts. Muscle system process and oxidoreductase activity acting on CH-OH group of donors were the most notable GO term of downregulated mRNAs. KEGG enrichment analysis suggested human papillomavirus infection crucial for upregulated mRNAs and gastric acid secretion remarkable for down regulated mRNAs. Outcomes of GO/KEGG analysis based on downregulated and upregulated miRNAs were separately revealed in Fig. S1C and S1D. Cellular amino acid metabolic process and transcription coactivator activity had smallest P-value for downregulated miRNAs enrichment analysis. KEGG analysis revealed that downregulated miRNAs were enriched in ErbB signaling pathway.

Following the DE-mRNAs signature query, small molecule agents including ixazomib, gemcitabine and tipifarnib were determined as potential therapeutic. The 15 compounds with the smallest connectivity scores were displayed in Fig. S2A. Agents’ mechanism of actions and their relative gene counts were exhibited in Fig. S2B.

Figure 5.

Patients were divided into different groups according to their clinical features, such as gender, age, tumor grade, and tumor stage. Subgroup analysis based on TCGA database suggested that ANXA2 expression was differed between patients older than 65 year and patients younger than 65. ELL2 was overexpressed in tissue with higher grade and stage.

Patients were divided into different groups according to their clinical features, such as gender, age, tumor grade, and tumor stage. Subgroup analysis based on TCGA database suggested that ANXA2 expression was differed between patients older than 65 year and patients younger than 65. ELL2 was overexpressed in tissue with higher grade and stage.

Table 5

Infiltration of 22 types of immune cells in STAD versus controls

ID of immune cellsProportion in STADProportion in normalP-valueFC
B cells naive0.05440.04710.54041.1548
B cells memory0.01270.01430.76230.8881
Plasma cells0.05040.17510.00000.2879
T cells CD80.10690.10860.48160.9838
T cells CD4 naive0.00030.00010.77722.5264
T cells CD4 memory resting0.19170.22320.09620.8589
T cells CD4 memory activated0.04040.00860.00014.7134
T cells follicular helper0.02430.01150.00092.1160
T cells regulatory (Tregs)0.06260.03650.00011.7167
T cells gamma delta0.00270.00220.69041.2177
NK cells resting0.01430.00670.17612.1243
NK cells activated0.02130.02590.12710.8252
Monocytes0.00630.0222< 0.00010.2819
Macrophages M00.12730.0078< 0.000116.3585
Macrophages M10.06750.0223< 0.00013.0288
Macrophages M20.11210.12330.50830.9097
Dendritic cells resting0.01940.03380.25370.5741
Dendritic cells activated0.01270.00790.80091.5948
Mast cells resting0.03030.0978< 0.00010.3099
Mast cells activated0.02520.00760.20943.3190

Figure 6.

(A) This bar plot exhibited cell fraction of 22 immune cells in STAD versus normal tissue. (B) The relationship between differently infiltrated 8 immune cells and miRNA/mRNA axes was analyzed by Spearman’s correlation. (C) The association between expression levels of miRNA/mRNA axes and tumor microenvironment features including tumor mutation burden, tumor purity, ESTIMATE score, immune score and stromal score.

(A) This bar plot exhibited cell fraction of 22 immune cells in STAD versus normal tissue. (B) The relationship between differently infiltrated 8 immune cells and miRNA/mRNA axes was analyzed by Spearman’s correlation. (C) The association between expression levels of miRNA/mRNA axes and tumor microenvironment features including tumor mutation burden, tumor purity, ESTIMATE score, immune score and stromal score.

3.3Identification of miRNA-mRNA networks for STAD

To predict the role of DE-miRNAs in STAD, we collected their potential target mRNAs from TarBase and MiRTarBase. Figure 2 showed the miRNA-mRNA network. After searching in PubMed, totally 12 miRNA-mRNA axes were tested in gastric FFPE tissue via RT-qPCR.

3.4Corroboration of the DE-mRNAs and DE-miRNAs with RT-qPCR

Specific fold change and p-value of all miRNAs and mRNAs calculated via 2 - ΔΔCt method was listed in Table 4. MiRNA and mRNA expression data were exhibited in Fig. 3A and B separately. KIT (p< 0.0001), LIFR (p= 0.0007) and ELL2 (p= 0.0004) showed significant downregulated trend in gastric cancer tissue. Only ANXA2 (p= 0.0002) was upregulated in tumor tissue. Hsa-miR-26a-5p (p= 0.0005), hsa-miR-29c-3p (p= 0.029) and hsa-miR-1-3p (p= 0.0164) were observed downregulated in STAD compared with normal controls. Hsa-miR-21-5p (p< 0.0001), hsa-miR-20a-5p (p= 0.0175), hsa-miR-221-3p (p= 0.0467) and hsa-miR-301a-3p (p= 0.0324) expression level of gastric tumor tissue were remarkably higher than that of adjacent normal tissue. Spearman’s correlation analysis confirmed two reverse regulatory miRNA-mRNA axes. As is shown in Fig. 3C, hsa-miR-301a-3p was negatively regulated with ELL2 (p= 0.0123, r=-0.3212). Hsa-miR-1-3p expression was significantly correlated with ANXA2 (p< 0.0001, r=-0.5722) (Fig. 3D).

3.5Evaluating the diagnosis value of miRNA-mRNA axes in STAD

None of the GEO datasets contained both miRNA and mRNA expression profile of the same specimen. As a result, only expression patterns gained from TCGA data and RT-qPCR were applied in measuring predictive value. MiR-1-3p, ANXA2, MiR-301a-3p and ELL2 were combined as a panel with binary logistic regression ROC curves of 14 models and the AUC of each curve calculated with TCGA data were presented in Fig. 4A. The panel containing 4 signatures showed the best diagnosis value (AUC = 0.9208, 95% CI: 0.8804 to 0.9612). The DCA shown in Fig. 4B suggested satisfied diagnosis efficiency of the complex model in TCGA and RT-qPCR validation.

3.6Clinical subgroup analysis and survival analysis based on TCGA clinical information

GDC data portal provided clinical features and follow-up time of all clinical cases recorded in TCGA-STAD database. Patients were divided by age at first diagnosis, gender, tumor grade, and tumor stage. ANXA2 expression level was higher in patients with older first diagnosis age (p= 0.0269) (Fig. 5A). No significant difference of the miRNA-mRNA axes expression was observed between male and female patients (Fig. 5B). ELL2 expression was higher in patients with higher grade number (p= 0.0122) (Fig. 5C). ELL2 expression had a significant between early-stage (I + II) group and late-stage (III + IV) group (p= 0.0049) (Fig. 5D). Overall survival (OS) was employed to estimate tumor prognosis. As is demonstrated in Fig. S3, high ELL2 expression may suggest poor survival time (p= 0.041).

3.7Correlation analysis between miRNA-mRNA axes expression and tumor-related phenotypes

Among 22 types of immune cells, proportion of 8 immune cells showed significant difference between normal tissue and gastric tumor (Fig. 6A and Table 5). Results of Spearman’s correlation between immune cell proportion and miRNA/mRNA expression were presented in Fig. 6B. Percentages of M1 macrophages, monocytes, and activated memory CD4 T cells were correlated with miR-1-3p/ANXA2 axes. Regulatory T cells (Tregs) and follicular T helper cells were correlated with miR-301a-3p/ELL2 axes. Resting mast cells and M0 macrophages were significantly correlated with these 2 miRNA-mRNA axes mentioned above. ESTIMATE was applied to estimate the samples’ stroma and immunity levels. As is exhibited in Fig. 6C, the network of miR-1-3p and ANXA2 interacted with tumor mutational burden and tumor immune microenvironment. The p-value for Fig. 6B and C is presented in the supplementary table.


Dysregulation of miRNA in cancers have raised wide concern in the past few decades [19]. In current study, we identified 228 DE-mRNAs (105 up and 123 down) and 67 DE-miRNAs (22 up and 16 down) via bioinformatic method based on GEO datasets and TCGA database. Via binding to 3’UTR of mRNA, miRNA posttranscriptionally suppress target mRNA expression [20]. TarBase and MiRTarBase helped to define negatively regulated miRNA-mRNA axes. Previous studies on miRNA-mRNA interactions have emphasized their essential role in occurrence, invasion, immune escape, and metastasis of gastric cancer [21]. However, these researches based on bioinformatic methods did not verify the inverse regulatory relationship of DE-miRNAs and DE-mRNAs in human specimen. In present study, we measured the expression level of key miRNAs and mRNA by conducting RT-qPCR in STAD tissue and adjacent normal tissue. Finally, we considered miR-1-3p/ANXA2 and miR-301a-3p/ELL2 as specific miRNA/mRNA signatures for STAD.

Pathway enrichment analysis based on dysregulated mRNAs indicated that gastric acid secretion, mitotic cell cycle phase transition, Human papillomavirus infection, and alcohol dehydrogenase (NAD (P) +) activity may have profound effect on gastric tumor. Several well-known signaling pathways, including ErbB signaling pathway, Wnt signaling pathway, MAPK signaling pathway, and PI3K-Akt signaling pathway, were associated with STAD, which correspond with previous studies [22, 23]. Several well-known cancer suppressor drugs such as ixazomib, gemcitabine and tipifarnib were revealed in the outcome of Connectivity map. Among all these compounds, gemcitabine, dorsomorphin, glibenclamide, and tipifarnib have been reported to exert negative effects on development of gastric cancer in cell experiments, animal experiments or clinical trials [24, 25, 26, 27].

Ultimately, the networks of miR-1-3p (down)/ANXA2 (up) and miR-301a-3p (up)/ELL2 (down) were recognized as essential miRNA-mRNA signatures. Downregulations of MiR-1-3p in gastric cancer tissue and cell lines were observed in previous studies [28, 29]. MiR-1-3p was reported to function as a suppressor in cell proliferation, apoptosis, and migration of STAD [30, 31]. ANXA2, short for annexin A2, is a member of the calcium-dependent phospholipid-binding protein family annexins [32]. By interacting with several tumor-related signal pathways, including mTOR, p38MAPK and AKT pathways, ANXA2 plays a promoting role in drug resistance and tumor growth of STAD [33, 34]. It is reported that circulating ANXA2 level may serve as a candidate biomarker for diagnosis and chemotherapy sensitivity assessment in gastric cancer [35]. Overexpression of miR-301a-3p was verified in esophageal carcinoma [36], lung cancer [37] and gastric cancer [38]. In HER2-positive gastric cancer, miR-301a-3p mediate the trastuzumab resistance [39]. Elongation Factor for RNA Polymerase II 2 (ELL2) was reported to interact with the retinoblastoma pathway and thus influence cell proliferation, invasion, and migration of prostate cancer [40]. To sum up, miR-1-3p, miR-301a-3p, ANXA2 and ELL2 are of significance in pathogenesis and progress of STAD. ROC based on these two miRNA-mRNA networks indicated that the complex model containing all 4 signatures showed best diagnostic value comparing to other 14 modes, the AUC of which equals to 0.9208 (95% CI: 0.8804 to 0.9612). DCA conducted with TCGA-STAD database also affirmed the diagnostic efficacy of these 2 miRNA-mRNA axes.

Survival analysis based on overall survival (OS) suggested that low ELL2 expression tend to indicate better prognosis in gastric cancer (p= 0.041 in K-M survival curve), which has not been reported ever before, while overexpression of ELL2 indicates poor prognosis in glioblastoma multiforme [41]. According to the research of Xiaodong Xu et al., high miR-301a-3p expression levels in gastric tumor tissue may suggest worse prognosis [42]. Although K-M curve did not support the correlation between increased miR-301a-3p expression and poor prognosis in present study, it is unreasonable to totally deny the prognostic efficacy of miR-301a-3p for STAD. Likewise, the miR-1-3p/ANXA2 axes did not show significant relevancy with OS in this research. ANXA2 was included in a prognostic model containing 10 genes for gastric cancer [43]. MiR-1-3p were also absorbed in a model containing six miRNAs for bone metastasis prediction of gastric cancer [44]. Since detailed survival analyzes are often restricted by insufficient information collected in follow-up studies, the prognosis efficacy of the miRNA-mRNA axes needs to be further measured.

Correlation between miRNA/mRNA axes and immune cells proportion emphasized the crucial role of resting mast cells and M0 macrophages in gastric adenocarcinomas. Proportion of resting mast cell declined in STAD tissue, while percentage of M0 macrophage increased in tumor tissue. Mast cells activated by interleukin (IL)-33 accelerate macrophage infiltration and thus promote gastric tumor growth [45]. Enriched mast cells capable of interacting with TNF-α-PD-L1 pathway are independent risk factor of poor survival for STAD patients [46]. As for macrophages numerous studies has revealed the significant role of macrophage M2 in development, metastasis, and response to chemotherapy of gastric cancer [47, 48, 49]. Although some scholars noticed the differences in proportion of M0 macrophages between STAD and normal tissue, few studies focused on the specific role of M0-macrophage in tumor microenvironment, which deserves further research [50]. So, it is reasonable to consider that miR-301a-3p/ELL2 and miR-1-3p/ANXA2 may have effect on tumor immune microenvironment.

Current study gave an overview of miRNA-mRNA regulatory network for STAD. However, there are still some shortcomings such as insufficient sample size and lack of follow-up data in this study. Besides, patients with lymph node metastasis, distant metastasis, or vascular invasion have no surgical opportunity, so the samples collected in this study did not include patients with advanced tumors. Therefore, studies based on more clinical samples and long-term follow-up data are needed.


In summary, we carried out comprehensive bioinformatic analysis and experimental verification for miRNA-mRNA network related in STAD. ThesemiRNA-mRNA axes may function as biomarkers for early diagnosis, long-term prognosis, and drug resistance assessment of gastric cancer. Our findings may provide new ideas for STAD treatment.

Authors’ contributions

Conception: Feng Gao, Hang Yang and Wei Zhu.

Interpretation or analysis of data: Shuang Peng and Guoxin Song.

Preparation of the manuscript: Shuang Peng, Jingfeng Zhu and Hao Zhang.

Revision for important intellectual content: Shuang Peng, Shiyu Zhang and Cheng Liu.

Supervision: Wei Zhu.

Supplementary data

The supplementary files are available to download from


We thank openbiox community and Hiplot team ( for providing technical assistance and valuable tools for data analysis and visualization. This study was supported by the Science and Education Health Promotion Project of Wujiang District, Suzhou [Grant number: WWK201908].



H. Sung et al., Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries, CA Cancer J Clin 71: (3) ((2021) ), 209–249.


R. Sitarz et al., Gastric cancer: Epidemiology, prevention, classification, and treatment, Cancer Manag Res 10: ((2018) ), 239–248.


E.C. Smyth et al., Gastric cancer, The Lancet 396: (10251) ((2020) ), 635–648.


C. Allemani et al., Global surveillance of cancer survival 1995–2009: Analysis of individual data for 25 676 887 patients from 279 population-based registries in 67 countries (CONCORD-2), The Lancet 385: (9972) ((2015) ), 977–1010.


P. Rawla and A. Barsouk, Epidemiology of gastric cancer: Global trends, risk factors and prevention, Prz Gastroenterol 14: (1) ((2019) ), 26–38.


R. Rupaimoole and F.J. Slack, MicroRNA therapeutics: towards a new era for the management of cancer and other diseases, Nature Reviews Drug Discovery 16: (3) ((2017) ), 203–222.


J. Lu et al., MicroRNA expression profiles classify human cancers, Nature 435: (7043) ((2005) ), 834–838.


Z. Ali Syeda et al., Regulatory Mechanism of MicroRNA Expression in Cancer, Int J Mol Sci 21: (5) ((2020) ),


X. Li et al., miRNA-223 promotes gastric cancer invasion and metastasis by targeting tumor suppressor EPB41L3, Mol Cancer Res 9: (7) ((2011) ), 824–833.


Y. Dong et al., miR-221-3p and miR-15b-5p promote cell proliferation and invasion by targeting Axin2 in liver cancer, Oncol Lett 18: (6) ((2019) ), 6491–6500.


R.A. de Mello et al., Current and potential biomarkers in gastric cancer: A critical review of the literature, Future Oncol 17: (25) ((2021) ), 3383–3396.


K. Pang et al., miR-15a-5p targets PHLPP2 in gastric cancer cells to modulate platinum resistance and is a suitable serum biomarker for oxaliplatin resistance, Neoplasma 67: (5) ((2020) ), 1114–1121.


G. Yang et al., Identification of Potentially Functional CircRNA-miRNA-mRNA Regulatory Network in Gastric Carcinoma using Bioinformatics Analysis, Med Sci Monit 25: ((2019) ), 8777–8796.


Z. Fu et al., Construction of miRNA-mRNA-TF Regulatory Network for Diagnosis of Gastric Cancer, Biomed Res Int 2021: ((2021) ), 9121478.


Z. Dong et al., Identification of circRNA-miRNA-mRNA networks contributes to explore underlying pathogenesis and therapy strategy of gastric cancer, J Transl Med 19: (1) ((2021) ), 226.


Y.J. Guan et al., Identification of circRNA-miRNA-mRNA regulatory network in gastric cancer by analysis of microarray data, Cancer Cell Int 19: ((2019) ), 183.


P. Murthi et al., GAPDH, 18S rRNA and YWHAZ are suitable endogenous reference genes for relative gene expression studies in placental tissues from human idiopathic fetal growth restriction, Placenta 29: (9) ((2008) ), 798–801.


D. Cassol et al., Identification of reference genes for quantitative RT-PCR analysis of microRNAs and mRNAs in castor bean (Ricinus communis L.) under drought stress, Plant Physiology and Biochemistry 106: ((2016) ), 101–107.


B. He et al., miRNA-based biomarkers, therapies, and resistance in Cancer, Int J Biol Sci 16: (14) ((2020) ), 2628–2647.


M.R. Fabian et al., Regulation of mRNA Translation and Stability by microRNAs, Annual Review of Biochemistry 79: (1) ((2010) ), 351–379.


X. Zeng et al., Research progress on the circRNA/lncRNA-miRNA-mRNA axis in gastric cancer, Pathol Res Pract 238: ((2022) ), 154030.


Z.N. Lei et al., Signaling pathways and therapeutic interventions in gastric cancer, Signal Transduct Target Ther 7: (1) ((2022) ), 358.


F. Meric-Bernstam et al., Advances in HER2-Targeted Therapy: Novel agents and opportunities beyond breast and gastric cancer, Clin Cancer Res 25: (7) ((2019) ), 2033–2041.


K. Kawaguchi et al., Combination of gemcitabine and docetaxel regresses both gastric leiomyosarcoma proliferation and invasion in an imageable patient-derived orthotopic xenograft (iPDOX) model, Cell Cycle 16: (11) ((2017) ), 1063–1069.


K. Zhou et al., Bone morphogenetic protein 4 is overexpressed in and promotes migration and invasion of drug-resistant cancer cells, Int J Biol Macromol 101: ((2017) ), 427–437.


K. Magierowska et al., Oxidative gastric mucosal damage induced by ischemia/reperfusion and the mechanisms of its prevention by carbon monoxide-releasing tricarbonyldichlororuthenium (II) dimer, Free Radic Biol Med 145: ((2019) ), 198–208.


N. Egawa et al., Antitumor effects of low-dose tipifarnib on the mTOR signaling pathway and reactive oxygen species production in HIF-1α-expressing gastric cancer cells, FEBS Open Bio 11: (5) ((2021) ), 1465–1475.


N. Wang et al., HDAC6/HNF4alpha loop mediated by miR-1 promotes bile acids-induced gastric intestinal metaplasia, Gastric Cancer 24: (1) ((2021) ), 103–116.


J. Ke et al., MiR-1-3p suppresses cell proliferation and invasion and targets STC2 in gastric cancer, Eur Rev Med Pharmacol Sci 23: (20) ((2019) ), 8870–8877.


P. Deng et al., LINC00242/miR-1-3p/G6PD axis regulates Warburg effect and affects gastric cancer proliferation and apoptosis, Mol Med 27: (1) ((2021) ), 9.


X.Q. Lin et al., mir-1 inhibits migration of gastric cancer cells, Front Biosci (Landmark Ed) 25: (3) ((2020) ), 452–462.


L. Mao et al., EphA2-YES1-ANXA2 pathway promotes gastric cancer progression and metastasis, Oncogene 40: (20) ((2021) ), 3610–3623.


Y. Li et al., S100A10 Accelerates Aerobic Glycolysis and Malignant Growth by Activating mTOR-Signaling Pathway in Gastric Cancer, Front Cell Dev Biol 8: ((2020) ), 559486.


Z.D. Zhang et al., Annexin A2 is implicated in multi-drug-resistance in gastric cancer through p38MAPK and AKT pathway, Neoplasma 61: (6) ((2014) ), 627–637.


F. Tas et al., Circulating annexin A2 as a biomarker in gastric cancer patients: Correlation with clinical variables, Biomed Pharmacother 69: ((2015) ), 237–241.


N. Zhang and J.F. Liu, MicroRNA (MiR)-301a-3p regulates the proliferation of esophageal squamous cells via targeting PTEN, Bioengineered 11: (1) ((2020) ), 972–983.


L. Li et al., beta-elemene suppresses Warburg effect in NCI-H1650 non-small-cell lung cancer cells by regulating the miR-301a-3p/AMPKalpha axis, Biosci Rep 40: (6) ((2020) ).


X. Xia et al., Hypoxic gastric cancer-derived exosomes promote progression and metastasis via MiR-301a-3p/PHD3/HIF-1alpha positive feedback loop, Oncogene 39: (39) ((2020) ), 6231–6244.


J. Guo et al., miR-301a-3p induced by endoplasmic reticulum stress mediates the occurrence and transmission of trastuzumab resistance in HER2-positive gastric cancer, Cell Death Dis 12: (7) ((2021) ), 696.


X. Qiu et al., Physical and functional interactions between ELL2 and RB in the suppression of prostate cancer cell proliferation, migration, and invasion, Neoplasia 19: (3) ((2017) ), 207–215.


Q. Huang et al., Up-regulated microRNA-299 corrected with poor prognosis of glioblastoma multiforme patients by targeting ELL2, Jpn J Clin Oncol 47: (7) ((2017) ), 590–596.


X. Xu et al., Upregulation of miRNA301a3p promotes tumor progression in gastric cancer by suppressing NKRF and activating NFkappaB signaling, Int J Oncol 57: (2) ((2020) ), 522–532.


C. Liang et al., Identification and validation of a pyroptosis-related prognostic model for gastric cancer, Front Genet 12: ((2021) ), 699503.


J. Ma et al., Prognostic microRNAs associated with phosphoserine aminotransferase 1 in gastric cancer as markers of bone metastasis, Front Genet 13: ((2022) ), 959684.


M.F. Eissmann et al., IL-33-mediated mast cell activation promotes gastric cancer through macrophage mobilization, Nat Commun 10: (1) ((2019) ), 2735.


Y. Lv et al., Increased intratumoral mast cells foster immune suppression and gastric cancer progression through TNF-alpha-PD-L1 pathway, J Immunother Cancer 7: (1) ((2019) ), 54.


Y. Chen et al., Tumor-recruited M2 macrophages promote gastric and breast cancer metastasis via M2 macrophage-secreted CHI3L1 protein, J Hematol Oncol 10: (1) ((2017) ), 36.


R. Kim et al., Early tumor-immune microenvironmental remodeling and response to first-line fluoropyrimidine and platinum chemotherapy in advanced gastric cancer, Cancer Discov 12: (4) ((2022) ), 984–1001.


V. Gambardella et al., The role of tumor-associated macrophages in gastric cancer development and their potential as a therapeutic target, Cancer Treat Rev 86: ((2020) ), 102015.


H. Chen et al., Identification and validation of CYBB, CD86, and C3AR1 as the key genes related to macrophage infiltration of gastric cancer, Front Mol Biosci 8: ((2021) ), 756085.