3.1 Identification of DEMs, DELs and DEMis
After standardization of the GEO datasets, 56, 70 and 34 DEMis between benign nevus tissues and primary melanoma tissues, were identified respectively in GSE24996, GSE35579 and GSE62372 (Tab. 1, Fig. 1A-1F). The candidates 18 miRNAs were shared in at least two datasets (Fig. 2A), including hsa-miRNA-378a-3p, hsa-miRNA-23b-3p, hsa-miRNA-140-3p, hsa-miRNA-99a-5p, hsa-miRNA-100-5p, hsa-miRNA-204-5p, hsa-miRNA-211-5p, hsa-miRNA-205-5p, hsa-miRNA-224-5p, hsa-miRNA-200b-3p, hsa-miRNA-200c-3p, hsa-miRNA-125b-5p, hsa-miRNA-149-5p, hsa-miRNA-21-5p, hsa-miRNA-20b-5p, hsa-miRNA-424-5p, hsa-miRNA-203a-3p and hsa-miRNA-1826. According to method 2.3, 2,361 mRNAs and 277 lncRNAs were predicted using these miRNAs. And, we ruled out two of these 18 DEMis, hsa-miRNA-203a-3p and hsa-miRNA-1826, because no predicted gene was found in starbase according to method 2.3. In addition, 5,953 DEMs and 665 DELs between benign nevus tissues and primary melanoma tissues, were identified in GSE112509 (Fig. 1G and 1H). As a result, a total of 898 DEMs and 53 DELs were selected for further analysis according to method 2.3 (Fig. 2B and 2C). Finally, 898 DEMs,53 DELs and 16 DEMis were selected for further reconstruct the lncRNA-miRNA-mRNA (ceRNA) network.
3.2 Reconstruction of lncRNA-miRNA-mRNA (ceRNA) network
The lncRNA-miRNA-mRNA (ceRNA) network, consisting of 53 lncRNA nodes, 16 miRNA nodes, 898 mRNA nodes and 609 edges was reconstructed and visualized by using Cytoscape (Fig. 3A).
3.3 KEGG pathway and GO enrichment analysis of lncRNAs based on ceRNA network
We used DAVID to analysis the biological classification of DEMs according to method 2.5. The results of top 15 significant GO terms and KEGG pathways are shown in Table 2 and Fig. 3B-3E. 60 pathways are significantly enriched through KEGG pathway analysis, including PI3K-Akt signaling pathway, focal adhesion, proteoglycans in cancer, pathway in cancer and most importantly in melanomagenesis. The results of GO-BP analysis reveal that 172 enriched terms, particularly in various regulation of transcription such as positive regulation of transcription from RNA polymerase II promoter, positive regulation of transcription (DNA-templated), transcription from RNA polymerase II promoter and etc.
3.4 Hub gene selection
According to the node degree in ceRNA network, we found that three lncRNAs, MALAT1, LINC00943, LINC00261, had the highest number of lncRNA-miRNA and miRNA-mRNA pairs, suggesting that these three lncRNAs could be chosen as hub nodes, and the results are shown in Table 3. Therefore, these three lncRNAs might play an essential role in melanomagenesis, which might be considered as the key lncRNAs.
3.5 Reconstruction of MALAT1/LINC00943/LINC00261-miRNA-mRNA sub-networks
MALAT, LINC00943,LINC00261 and their paired miRNA and mRNA were used to reconstruct key ceRNA sub-network. The MALAT1 ceRNA network consists of 1 lncRNA nodes, 9 miRNA nodes, 158 mRNA nodes and 209 edges as shown in Fig. 4A. The LINC00943 ceRNA network consists of 1 lncRNA nodes, 7 miRNA nodes, 182 mRNA nodes and 209 edges as shown in Fig. 5A. The LINC00261 ceRNA network consists of 1 lncRNA nodes, 5 miRNA nodes, 123 mRNA nodes and 163 edges as shown in Fig. 6A. The results of functional analysis revealed that 75 GO-BP, 21 GO-CC, 15 GO-MF and 20 pathways were enriched in the MALAT1-miRNA-mRNA sub-network, 67 GO-BP, 14 GO-CC, 17 GO-MF and 13 pathways were enriched in the LINC00943-miRNA-mRNA sub-network, 42 GO-BP, 7 GO-CC, 10 GO-MF and 7 pathways were enriched in the LINC00261-miRNA-mRNA sub-network. The results of top 10 significant GO terms and KEGG pathways of these three lncRNAs were shown in Fig. 4B-5E, Fig. 5B-5E, Fig. 6B-6E, and Table 4-6.
The MALAT1、LINC00943 and LINC00261 sub-network show that GO-BP analysis revealed 75, 69 and 42 enriched terms respectively, which are also significantly enriched in regulation of transcription such as positive regulation of transcription from RNA polymerase II promoter, positive regulation of transcription (DNA-templated), transcription from RNA polymerase II promoter and more. Pathways analysis reveal that 19, 13 and 7 pathways are significantly enriched in tumor-related pathways, including pathway in cancer, focal adhesion, PI3K-Akt signaling pathway and more.
3.6 Expression of MALAT1、LINC00943 and LINC00261 are higher in tumor tissues
To confirm the expression of MALAT1、LINC00943 and LINC00261 in melanoma tissues, we had evaluated MALAT1、LINC00943 and LINC00261 expression level in cancer tissues from 12 melanoma patients and 3 healthy samples via qRT-PCR, as shown in Fig. 7. The results show that the expression of MALAT1、LINC00943 and LINC00261 are significantly higher in tumor tissues compared with healthy samples (p = 0.0243, p = 0.0005, p < 0.0001, respectively). Also, the expression of MALAT1、LINC00943 and LINC00261 are significantly higher in tumor tissues compared with its adjacent normal tissue (p = 0.0002, p < 0.0001, p < 0.0001, respectively). However, no significant difference was observed between healthy samples and adjacent normal skin tissues in expression of MALAT1、LINC00943 and LINC00261 (p = 0.366, p = 0.379, p = 0.262, respectively). The results are consistent with that discussed above. Thus, the expression of MALAT1、LINC00943 and LINC00261 are increased in melanoma and may responsible for the tumorigenesis of melanoma.