An integrated analysis of key microRNAs, regulatory pathways and clinical relevance in bladder cancer.

OBJECTIVE: The aim of this study was to identify the key microRNAs (miRNAs) and their regulatory networks in bladder cancer (BC).
MATERIALS AND METHODS: Three miRNA and three gene expression microarray datasets were downloaded for analysis from Gene Expression Omnibus database. The differentially expressed miRNAs (DEMs) and differentially expressed genes (DEGs) were accessed by the use of GEO2R. Gene ontology process and Kyoto Encyclopedia of Genes and Genomes pathway enrichment analyses were performed by using the Database for Annotation, Visualization and Integrated Discovery program. Protein-protein interaction (PPI) and miRNA-mRNA regulatory networks were established by using the Search Tool for the Retrieval of Interacting Genes/Proteins and Cytoscape tool. Besides, the results and clinical significance were validated in The Cancer Genome Atlas (TCGA) dataset.
RESULTS: A total of 18 significant DEMs, 121 upregulated DEGs and 199 downregulated DEGs were identified. Functional enrichment analysis showed that significant DEGs were related to cell cycle and MAPK pathway in BC. Key DEGs such as CDK1, CCNB1, VGL and PRKCA were found as the hub genes in PPI networks. TCGA analysis supported our results, and the miRNAs were correlated with the pathological stages and survival of BC patients.
CONCLUSION: In this study, we found 18 DEMs that may play key roles in the regulatory networks of BC. The higher expression of miR-99a, miR-100, miR-125b, miR-145, miR-214 and miR-487b or the lower expression of miR-138 and miR-200a can indicate poor survival in the prognosis of BC. Further experimental studies are required to test our results.

Introduction

Bladder cancer (BC) is the leading cause of death among urological cancers, with an estimated 81,190 new cancer cases and 17,240 deaths in the US in 2018.1

Approximately 70% of BC is non-muscle-invasive bladder cancer (NMIBC) at first diagnosis; however, roughly 20% of NMIBC will progress to muscle-invasive BC (MIBC).2 Owing to the aggressive metastasis and delayed diagnosis, the prognosis of MIBC is still poor.3 Although some conventional methods such as surgery and cytotoxic chemotherapy can be applied, the recurrence rate can be up to 70% within 2 years, which makes the treatment of BC unsatisfactory.4 The thorough understanding of molecular mechanism and biological markers can help in the early detection, risk stratification, appropriate therapeutic intervention and prognosis in patients with BC.

MicroRNA (miRNA) is a group of single-stranded short RNA molecules (generally 21 nucleotides) that take part in the posttranscriptional gene regulation.5 Emerging evidence reveals that miRNAs are involved in the various stages such as tumor initiation, proliferation, migration and invasion.6,7 MiRNAs can also functionally act as oncogene or tumor suppressors.8 The detection technique of miRNA in serum, urine and biopsy specimen is available; therefore, miRNA is an attractive source of noninvasive biomarkers.9 The quantitative reverse transcription polymerase chain reaction (qRT-PCR), small RNA sequencing and microarray are widely applied in BC researches in the past decade. In addition, the utility of miRNA as diagnosis or prognosis markers has generated a great interest. However, the results of different studies are inconclusive due to the heterogeneous ethnicity, sample size, different sample types and less clinical parameters.

To the best of our knowledge, few studies have reported the miRNAs and miRNA–mRNA regulatory networks by integrating multiple microarray datasets. Therefore, a comprehensive analysis was administrated in this study to further elucidate the key miRNAs, regulatory pathways and clinical significance in BC.

Materials and methods

Microarray data

Three miRNA expression profiles (GSE36121, GSE39093 and GSE68594) and three gene expression profiles (GSE27448, GSE40355 and GSE52519) were obtained by the search of the Gene Expression Omnibus (GEO, http://www.ncbi.nlm.nih.gov/geo). The profiles based on BC cell lines or xenograft models were excluded; thus, only BC patients and healthy controls were included in this study. Table S1 lists the detailed characteristics of the profiles. Every included dataset from the abovementioned profiles contains at least 10 BC samples and normal bladder tissues (Table 1).

Table 1

Characteristics of the microarray datasets from GEO database

Accession/ID	PMID	Platform	Number of BC tissues	Number of normal bladder tissues (control)	Gene/miRNA
GSE36121	22745731	GPL8227	n=16 (Ta–T1=13, T2–T4=3)	n=8	miRNA
GSE39093	23270926	GPL8786	n=10 (Ta–T1=5, T2–T4=5)	n=10	miRNA
GSE68594	26943042	GPL20120	n=33 (Ta–T1=18, T2–T4=15)	n=17	miRNA
GSE27448	21483740	GPL2895	n=10 (Ta–T1=6, T2–T4=4)	n=5	Gene
GSE40355	–	GPL13497	n=16 (Ta–T1=13, T2–T4=3)	n=8	Gene
GSE52519	–	GPL6884	n=9 (Ta–T1=3, T2–T4=6)	n=3	Gene

Abbreviations: BC, bladder cancer; GEO, Gene Expression Omnibus; miRNA, microRNA; PMID, PubMed ID.

Data processing

We performed the comparison between the two groups of samples (BC versus normal prostate tissues) in each GEO dataset to identify differentially expressed genes (DEGs) and differentially expressed miRNAs (DEMs). The comparison was launched by a Limma R package based on online program, GEO2R (http://www.ncbi.nlm.nih.gov/geo/geo2r/). The adjusted P-value from the Benjamini–Hochberg method can correct the false-positive results, and the “adjusted P-value<0.05” and “|log FC|>1” parameters were chosen as primary cutoff criteria to test the results. Among the DEMs from the four datasets, only those involved in two or more datasets could be considered the statistically significant DEMs.

Identification of miRNA targets

The MiRWalk 2.0 database provides a large collection of predicted and experimentally verified information about the binding targets of miRNAs.10 We downloaded the data of the significant DEMs and their target mRNAs from the MiRwalk 2.0 retrieval system, which was based on experimental literature reports. Despite the fact that validated target genes were from various disease models, the expression of those genes in BC might not be in consistent with other diseases. Hence, it is important to select the genes from the intersection between the significant DEM–mRNA intersection data and the DEGs from the three GEO datasets. These genes were finally defined as the significant DEGs.

Functional enrichment analysis

The Database for Annotation, Visualization and Integrated Discovery (DAVID) is an online platform that offers the functional annotation and analysis of enormous quantity of genes derived from various genomic resources.11 We used the DAVID database to conduct gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis on significant DEGs. The species was limited to homo sapiens, and the “adjusted P-value (from the Benjamini–Hochberg method) <0.05” was considered statistically significant.

Interaction and regulatory network establishment

The protein–protein interaction (PPI) was visualized, and its network was constructed by the use of Search Tool for the Retrieval of Interacting Genes/Proteins (STRING), whose database can collect and predict the interaction information from genomic context predictions, high-throughput laboratory experiments, co-expression, automated text mining and previous knowledge in databases.12 The data annotated from high-throughput laboratory experiments and previous knowledge in databases were imposed to access the more objective and reliable results, and the minimum interaction value was set at medium confidence (0.400).

The Targetscan (http://www.targetscan.org/) is a famous miRNA prediction database by searching for the presence of conserved 8mer, 7mer and 6mer sites that match the seed region of each miRNA. The miRanda algorithm from the “microRNA.org-Targets and Expression” is freely available and can be applied to the whole genome sequences using identified miRNA sequences.13 The two programs mentioned earlier were used to predict target mRNAs of the significant DEMs. We note that two miRNAs shared a common target mRNA; they might exist in a similar regulatory pathway. The miRNA–mRNA regulatory network was visualized by using Cytoscape 3.6.0.14

Validation and survival analysis

The Cancer Genome Atlas (TCGA) database provides abundant clinical information from huge sample size. Therefore, the dataset “TCGA-BLCA (https://portal.gdc.cancer.gov/)” was used to find if the significant DEMs were in correlation with survival outcome. The dataset contains 412 patients who were diagnosed as having BC in the period of 1999–2013. Among them, 304 patients are male, while 108 are female. As for the race, 327 BC patients are white, 44 are Asian and 23 are black or African-American, and 23 were not reported. Among the 404 clinical outcome events available, 227 patients were living and 177 were dead when the follow-up ended. We examined the significant DEMs’ expression, relationship with BC stages and survival. The “independent two-sample Student’s t-test” was applied to check the expression level, while the chi-square test was used on different stages. The Kaplan–Meier plots were constructed using Oncolnc.15 The “P-value<0.05” was considered statistically significant.

Results

DEG/DEM identification

After compared with normal bladder tissues, the results of GEO2R analysis showed that there were totally 159 DEMs (Figure 1A), 4,013 upregulated DEGs (Figure 1B) and 4236 downregulated DEGs (Figure 1C). Among them, 18 DEMs appeared in the results of two datasets at least, so they could be identified as significant DEMs. Interestingly, miR-145 was the only DEM in all three datasets simultaneously. The following DEMs were found in every two different datasets: miR-429, miR-487b, miR-494, miR-222, miR-335, miR-202, miR-212, miR-200a, miR-125b, miR-210, miR-214, miR-182, miR-99a, miR-100 and miR-138. Subsequently, we further compared clinical stages T2–T4 versus Ta–T1 and high histological grade versus low grade to identify DEMs which might be related to aggressive BC. Interestingly, only miR-601 and hsa-let-7c appeared in two dataset intersections (Figure S1).

Figure 1

Venn diagram of DEM/DEG selection in different datasets.

Notes: (A) DEMs; (B) upregulated DEGs; (C) downregulated DEGs.

Abbreviations: DEGs, differentially expressed genes; DEM, differentially expressed miRNA; miRNAs, microRNAs.

Using the MiRWalk 2.0-validated-target miRNA-gene retrieval system, we got 9,053 candidate genes of the 18 significant DEMs. The intersection number of these candidate genes and the upregulated DEGs which appeared in two or more datasets was 121, while the intersection number of these candidate genes and the downregulated DEGs which appeared in two or more datasets was 199. Therefore, the 121 upregulated and the 199 downregulated genes were identified as the final sets of significant DEGs.

GO and pathway enrichment

The GO ontology contains three terms: cellular component, molecular function and biological process. The three terms mentioned earlier were used to identify the key biological functions. The results demonstrated that the most significant GO terms for upregulated significant DEGs were “nucleoplasm,” “protein binding” and “cell division.” Whereas for the downregulated significant DEGs, the terms “costamere,” “heparin binding” and “cell–matrix adhesion” turned to be important (Table 2).

Table 2

Results of significant GO biological processes and the KEGG pathway analysis

Category	Description	Gene counts	Adjusted P-value
Upregulated gene
GO:0005654, CC	Nucleoplasm	49	4.40E-09
GO:0005634, CC	Nucleus	65	2.20E-06
GO:0000777, CC	Condensed chromosome kinetochore	7	8.70E-04
GO:0005819, CC	Spindle	8	9.30E-04
GO:0030496, CC	Midbody	8	1.10E-03
GO:0005515, MF	Protein binding	89	1.40E-06
GO:0003682, MF	Chromatin binding	16	6.00E-06
GO:0005524, MF	ATP binding	26	1.10E-03
GO:0051301, BP	Cell division	19	1.40E-08
GO:0007067, BP	Mitotic nuclear division	14	4.70E-06
GO:0007062, BP	Sister chromatid cohesion	9	1.30E-04
GO:0016925, BP	Protein sumoylation	8	3.00E-03
GO:0007052, BP	Mitotic spindle organization	5	8.40E-03
KEGG: hsa04110	Cell cycle	12	6.20E-08
Downregulated gene
GO:0043034, CC	Costamere	7	7.10E-06
GO:0043034, CC	Focal adhesion	18	1.00E-04
GO:0015629, CC	Actin cytoskeleton	12	1.60E-03
GO:0005576, CC	Extracellular region	35	3.90E-03
GO:0070062, CC	Extracellular exosome	51	4.30E-03
GO:0008201, MF	Heparin binding	11	3.20E-03
GO:0003779, MF	Actin binding	13	8.50E-03
GO:0007160, BP	Cell–matrix adhesion	8	4.30E-02
KEGG: hsa04110	MAPK signaling pathway	14	1.80E-03
KEGG: hsa05200	Pathways in cancer	15	2.10E-02
KEGG: hsa05205	Proteoglycans in cancer	10	3.80E-02

Abbreviations: BP, biological process; CC, cellular component; GO, gene ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; MF, molecular function.

Furthermore, the KEGG pathway analysis indicated that “cell cycle,” “MAPK signaling pathway,” “pathways in cancer” and “proteoglycans in cancer” pathways played an essential role in BC pathogenesis (Table 2).

PPI network

The PPI networks were established separately by the significant upregulated DEGs (Figure 2A) and significant downregulated DEGs (Figure 2B). As for the former, 64 nodes and 172 edges in total constituted the PPI network. As for the latter, the network was composed of 73 nodes and 80 edges.

Figure 2

PPI network of DEGs.

Notes: (A) Upregulated DEGs; (B) downregulated DEGs.

Abbreviations: DEGs, differentially expressed genes; PPI, protein–protein interaction.

In a single PPI network, the more edges a gene has, the more important role it plays (like a seed). We utilized the parameter “degree” to calculate edge counts of every single gene in a PPI network. The top 5% degree genes are listed in Table 3, which were assessed as the hub genes.

Table 3

The top 5% hub genes in the PPI networks

Ensembl gene ID	Gene symbol	Full gene name	Degree
Upregulated gene
ENSG00000170312	CDK1	Cyclin-dependent kinase 1	25
ENSG00000134057	CCNB1	Cyclin B1	21
ENSG00000072501	SMC1A	Structural maintenance of chromosomes 1A	16
ENSG00000117399	CDC20	Cell division cycle 20	14
Downregulated gene
ENSG00000035403	VCL	Vinculin	8
ENSG00000154229	PRKCA	Protein kinase C alpha	8
ENSG00000170345	FOS	Fos proto-oncogene, AP-1 transcription factor subunit	7
ENSG00000168610	STAT3	Signal transducer and activator of transcription 3	5

Abbreviation: PPI, protein–protein interaction.

MiRNA–mRNA regulatory network

As illustrated in Figure 3, the target genes were predicted using miRanda and miRTarbase, among which two or more miRNAs might target the same mRNA. For instance, both miR-92a and miR-145 had the interactions with cyclin-dependent kinase 4 (CDK4). There were 12 mRNA nodes that might interact with at least two miRNAs.

Figure 3

The miRNA–mRNA regulatory network.

Note: Round nodes, miRNA; square nodes, mRNA; green line, from Targetscan database; blue line, from miRanda prediction.

Abbreviation: miRNA, microRNA.

Validation and Kaplan–Meier curves

After comparing the expression level of miRNAs in the TCGA dataset, we found that 13 among the 18 significant DEMs were according to the results in GEO databases. In contrast, the expression of miR-92, miR-494, miR-429, miR-202 and miR-212 was not significantly different between BC samples and normal control bladder tissues. The relationship between miRNAs and the clinical TNM stages was also tested in the TCGA dataset by using chi-square analysis. The results of the significant DEMs are summarized in Table 4.

Table 4

Results of the significant DEM expression validation and the relationship with bladder cancer clinical stages in TCGA dataset

miRNA ID	Higher expression	P-value (t-test)*	Pathological stages	P-value (chi-square)
hsa-miR-145-3p	Normal	1.73e-06	Pathological T/N/M status	4.76e-04/8.90e-05/9.62e-04
hsa-miR-21-3p	Tumor	1.79e-04	Pathological M status	7.27e-03
hsa-miR-487b-3p	Normal	1.92e-11	Pathological T status	3.87e-03
hsa-miR-494-3p	–	–	Pathological T status	2.65e-02
hsa-miR-222-3p	Tumor	4.10e-02	Pathological M status	5.54e-03
hsa-miR-335-3p	Tumor	1.64e-06	–	–
hsa-miR-212-3p	–	–	Pathological T/M status	2.21e-03/1.60e-02
hsa-miR-200a-3p	Tumor	8.14e-05	Pathological T/M status	1.41e-03/7.23e-06
hsa-miR-210-3p	Tumor	3.43e-08	Pathological T/N status	3.01e-02/4.59e-03
hsa-miR-214-3p	–	–	Pathological T/M status	1.84e-07/1.74e-03
hsa-miR-182-3p	Tumor	1.59e-02	Pathological M status	2.34e-06
hsa-miR-100-3p	–	–	Pathological T status	4.15e-02
hsa-miR-145-5p	Normal	1.33e-04	Pathological T/N/M status	7.14e-03/1.79e-03/7.26e-05
hsa-miR-21-5p	Tumor	2.42e-06	Pathological M status	5.83e-06
hsa-miR-222-5p	Tumor	1.43e-02	–	–
hsa-miR-335-5p	Tumor	5.22e-05	Pathological M status	1.85e-02
hsa-miR-200a-5p	Tumor	7.68e-05	–	–
hsa-miR-125b-5p	Normal	4.25e-02	Pathological T/N/M status	3.26e-05/8.07e-04/5.08e-06
hsa-miR-214-5p	Tumor	9.80e-03	Pathological T/N/M status	5.29e-09/2.15e-03/5.31e-07
hsa-miR-182-5p	Tumor	1.03e-06	Pathological M status	6.18e-05
hsa-miR-99a-5p	Normal	1.94e-02	Pathological T/N/M status	1.63e-07/1.10e-03/4.17e-04
hsa-miR-100-5p	Normal	7.01e-04	Pathological T/N/M status	6.86e-04/6.83e-04/5.99e-07
hsa-miR-138-5p	Tumor	9.07e-03	Pathological T/N status	1.24e-02/3.73e-02

Note:

Independent two-sample Student’s t-test.

Abbreviations: DEMs, differentially expressed miRNAs; miRNAs, microRNAs; TCGA, The Cancer Genome Atlas.

We subsequently drew Kaplan–Meier plots (Figure 4) based on the TCGA survival data. The plots indicated that eight significant DEMs, such as miR-99a, miR-100, miR-125b, miR-138, miR-145, miR-200, miR-214 and miR-487b, were found to be related to overall survival (OS) in BC patients. The higher expression of the abovementioned miRNAs meant the poor OS in patients with BC except miR-138 and miR-200a.

Figure 4

The Kaplan–Meier plots of patients from TCGA dataset.

Notes: Survival plots of the bladder cancer patients with high and low expression of miR-145-3p (A), miR-126-5p (B), miR-99a-5p (C), miR-100-5p (D), miR-214-5p (E), miR-487b-3p (F), miR-200a-3p (G) and miR-138-5p (H).

Abbreviation: TCGA, The Cancer Genome Atlas.

Discussion

In the current study, 18 significant DEMs, 121 upregulated DEGs and 199 downregulated DEGs were identified after integrating and screening multiple datasets. The results of functional enrichment analysis indicated that the significant DEGs in BC were related to classical cancer pathways such as cell cycle and MAPK signaling pathways. The key DEGs such as CDK1, CCNB1, VCL and PRKCA were defined as the hub genes in the PPI networks. Significantly, most of the DEMs were correlated with tumor stages and OS in the TCGA datasets, which meant that the miRNAs have high clinical application value in the prognosis of BC.

BC shares many acknowledged pathways with other human cancers such as p53 and MAPK pathways. With the growing in-depth studies of BC, people have found a great deal of key molecules at gene, transcription regulation and protein interaction level. In these complicated processes, miRNAs nearly take part in all key cellular pathways because one miRNA can interact with many mRNAs and one mRNA can also interact with a lot of miRNAs. Intriguingly, miR-145 was the only one significant DEM that was identified in all three different datasets. We assumed that the source of heterogeneity might be ethnic and sample size. The expression level of miR-145 was higher in BC tissue and cell line than normal control.16,17 Minami et al18 reported that miR-145 could negatively regulate Warburg effect by silencing KLF4 and PTBP1. Long noncoding RNA (lncRNA) can also interact with miRNA in cancer pathways. LncRNA UCA-1 can promote BC cell migration and invasion by hsa-miR-145-ZEB1/2-FSCN1 pathway.19 miR-21 regulates cell proliferation and migration by targeting PTEN and p53 in BC.20 miR-202,21 miR-212,22 miR-222,23 miR-429,24 miR-200,25 miR-210,26 miR-214,27 miR-182,28 miR-99a,29 miR-10030 and miR-13831 were reported in BC studies, which were mostly dependent on cell lines. The studies on miRNA relationship with BC survival were insufficient. Urine is an ideal source of biomarkers due to the easily acquired and noninvasive property.32 Recently, it was reported that inflammation could drive tumor growth and progression.33,34 Mearini et al34 reported that although the expression of miR-145 in urine was higher in BC than in inflammation, and even higher in the high histological grade or high clinical stages, the results were not statistically significant. The expression of miR-145 was significantly higher in BC tissues in our results. It is assumed that the difference may be owing to low secretion amount or degradation. Therefore, we need the microarray studies concerning miRNAs in serum or urine as biomarkers in BC based on large sample sizes in the future.

Our studies may also provide hints for new target genes and pathways in the basic experiments of BC. Except for the previously reported miRNAs, there are also some miRNAs not reported in BC such as miR-335, miR-125b, miR-494 and miR-487b. We found that the higher miR-125b and miR487b expression suggested the poor survival based on 404 patients in TCGA. The mechanisms and in-depth pathways of these miRNAs should be investigated by experiments in the future.

By constructing the PPI network, we can recognize the key genes which miRNAs may interact with. Since the genes were filtered by the 18 significant DEMs’ potential targets, there were still 121 upregulated and 199 downregulated genes. Therefore, when it comes to all the 159 DEMs, the enormous and complicated miRNA–mRNA regulatory network can be imaginable. The hub genes of a network are always important like “seeds,” which connect the different signal pathways. In the current study, we identified several hub genes. Among the downregulated DEGs, the vinculin and protein kinase C alpha were identified. The vinculin was reported as a potential prognostic biomarker and the target for prostate cancer cell migration.35 Whereas in the upregulated DEGs, the hub gene CDK1 and cyclin B1 are associated with cell cycle and proliferation.36 Within the miRNA–mRNA network, we constructed that the key mRNAs such as CDK6 and CDK4 were cell-cycle-related modulators, while BCL2L11 was related to autophagy and apoptosis.37 In conclusion, the miRNA regulatory network may provide potential targets for new drug development and treatment.

Although this study first investigated miRNAs’ regulatory role in BC by integrating multiple microarray profiles, some limitations should be clarified. First, the gender, races of BC patients, tumor stages or grades were different in the included datasets, which may cause heterogeneity. Second, the source of microarray is only tissue. Body fluids such as serum, urine, may contain circulating miRNAs, which are noninvasive to be applied in the clinical practice. Third, although our results were validated in TCGA dataset, the integrating analysis results were not validated by the experiments such as qRT-PCR and functional experiments in vitro. Further studies are required to validate our results by the experiments on larger sample sizes and clinical follow-up in different ethnics.

Conclusion

In this study, we found a total of 18 DEMs which may play key roles in the pathways such as cell cycle and MAPK pathway in the regulatory networks of BC. The higher expression of miR-99a, miR-100, miR-125b, miR-145, miR-214 and miR-487b or the lower expression of miR-138 and miR-200a can indicate poor survival in the prognosis of BC. Further experimental studies are required to test our results.

Venn diagram of DEMs (A: Ta–T1 vs T2–T4; B: G3 vs G1) in different datasets.

Notes: (A) GSE36121: hsa-miR-135b, hsa-miR-508, hsa-miR-7, hsa-miR-21, hsa-miR-615, hsa-miR-124, hsa-miR-106b, hsa-miR-98, hsa-miR-132, hsa-miR-30e, hsa-miR-130a, hsa-miR-193a, hsa-miR-582, hsa-miR-151, hsa-let-7c. GSE39093: hsa-mir-1308, hsa-mir-192, hsa-mir-574, hsa-mir-720, hsa-mir-92b. GSE68594: hsa-let-7b, hsa-let-7c, hsa-let-7d, hsa-let-7f, hsa-let-7i, hsa-mir-155, hsa-mir-196, hsa-mir-429, hsa-mir-579. (B) GSE36121: hsa-miR-1224-5p, hsa-miR-1225-5p, hsa-miR-1228, hsa-miR-125a-3p, hsa-miR-134, hsa-miR-135a, hsa-miR-135b, hsa-miR-136, hsa-miR-146a, hsa-miR-150, hsa-miR-188-5p, hsa-miR-192, hsa-miR-193a-3p, hsa-miR-194, hsa-miR-197, hsa-miR-21, hsa-miR-224, hsa-miR-31, hsa-miR-345, hsa-miR-362-5p, hsa-miR-363, hsa-miR-371-5p, hsa-miR-423-5p, hsa-miR-508-5p, hsa-miR-518c, hsa-miR-532-5p, hsa-miR-574-5p, hsa-miR-575, hsa-miR-582-5p, hsa-miR-598, hsa-miR-601, hsa-miR-622, hsa-miR-629, hsa-miR-663, hsa-miR-671-5p, hsa-miR-7, hsa-miR-801, hsa-miR-940. GSE68594: hsa-let-7d, hsa-let-7f, hsa-let-7i, hsa-mir-122a, hsa-mir-184, hsa-mir-196a1, hsa-mir-196a2, hsa-mir-203, hsa-mir-206, hsa-mir-223, hsa-mir-34a, hsa-mir-595, hsa-mir-601, hsa-mir-611-P.

Abbreviations: DEMs, differentially expressed miRNAs; miRNAs, microRNAs.

Table S1

Detailed information of patients in each included GEO datasets

Accession/ID	Ethnic group	Number of males/females	Range of age (years)	Number of BC tissues with different histological grades	Number of NMIBC/MIBC tissues
GSE36121	Caucasian	22/2	47–80	G1: 3, G2: 6, G3: 7	NMIBC: 13, MIBC: 3
GSE39093	Chinese Han	10/0	55–73	NR	NMIBC: 5, MIBC: 5
GSE68594	Caucasian	39/5	36–86	G1: 4, G2: 14, G3: 15	NMIBC: 18, MIBC: 15
GSE27448	Caucasian	15/0	44–86	G1: 0, G2: 3, G3: 7	NMIBC: 6, MIBC: 4
GSE40355	NR	20/2	47–80	G1: 3, G2: 6, G3: 7	NMIBC: 13, MIBC: 3
GSE52519	NR	11/1	43–77	G1: 0, G2: 4, G3: 5	NMIBC: 3, MIBC: 6

Abbreviations: BC, bladder cancer; GEO, Gene Expression Omnibus; MIBC, muscle-invasive bladder cancer; NMIBC, non-muscle-invasive bladder cancer; NR, not reported.