Skip to main content

ORIGINAL RESEARCH article

Front. Genet., 01 April 2021
Sec. Cancer Genetics

Autophagy-Related Long Non-coding RNA Signature as Indicators for the Prognosis of Uveal Melanoma

\r\nYi Cui&#x;Yi Cui1†Mi Zheng&#x;Mi Zheng2†Jing Chen&#x;Jing Chen2†Nuo Xu*Nuo Xu2*
  • 1Department of Ophthalmology, Fujian Medical University Union Hospital, Union Clinical Medical College, Fujian Medical University, Fuzhou, China
  • 2Department of Ophthalmology, Fujian Provincial Hospital, Shengli Clinical Medical College, Fujian Medical University, Fuzhou, China

This study aimed to develop an autophagy-associated long non-coding RNA (lncRNA) signature to predict the prognostic outcomes of uveal melanoma (UM). The data of UM from The Cancer Genome Atlas (TCGA) were enrolled to obtain differentially expressed genes (DEGs) between metastasizing and non-metastasizing UM patients. A total of 13 differentially expressed autophagy genes were identified and validated in Gene Expression Omnibus, and 11 autophagy-related lncRNAs were found to be associated with overall survival. Through performing least absolute shrinkage and selection operator regression analyses, a six-autophagy-related lncRNA signature was built, and its efficacy was confirmed by receiver-operating characteristic, Kaplan–Meier analysis, and univariate and multivariate Cox regression analyses. A comprehensive nomogram was established and its clinical net benefit was validated by decision curve analysis. GSEA revealed that several biological processes and signaling pathways including Toll-like receptor signaling pathway, natural killer cell-mediated cytotoxicity, and B- and T-cell receptor signaling pathway were enriched in the high-risk group. CIBERSORT results showed that the signature was related to the immune response especially HLA expression. This signature could be deployed to assist clinicians to identify high-risk UM patients and help scientists to explore the molecular mechanism of autophagy-related lncRNAs in UM pathogenesis.

Introduction

Uveal melanoma (UM) is the most common and aggressive malignant neoplasm of the eye in adults, and up to 50–62% of UM patients will suffer from metastasis (Kujala et al., 2003; Amaro et al., 2017; Kaliki and Shields, 2017). Local interventions including transpupillary thermotherapy (Mashayekhi et al., 2015), photodynamic thermotherapy (Cerman and Cekic, 2015), brachytherapy (Collaborative Ocular Melanoma Study Group, 2006), and surgical enucleation (Shields and Shields, 2015) are the main strategies to treat the primary tumor, but the overall prognosis for metastatic UM remains poor. Therefore, it is imperative to identify high-risk UM patients with poor prognosis. Clinical features including patient age, tumor apical height, involvement of ciliary body, and largest basal diameter (LBD) of the tumor have been identified as prognostic factors for metastasis (Chew et al., 2015; Correa and Augsburger, 2016). Recently, scientists have proposed certain molecular signatures as predictors of high-risk UM, such as the gene expression profile (GEP) model (Zhang et al., 2014), DNA methylation-based model (Li et al., 2018), microRNA (miRNA) model (Xin et al., 2019), and alternative splicing event model (Wan et al., 2020). However, the prognostic value of long non-coding RNAs (lncRNAs) in UM has not been fully explored.

LncRNAs are a set of RNAs with a transcript length of over 200 nucleotides and have limited or no protein-coding capacity. Previous data has shown that abnormal activation of lncRNAs participates in various biological processes and human tumors including UM. A number of studies have shown that aberrant lncRNA expressions like MALAT1, FTH1P3, ZNNT1, and CANT1 played important roles in the initiation and progression of UM (Xing et al., 2017; Zheng et al., 2017; Li et al., 2020; Wu et al., 2020). Recent evidence highlights the involvement of lncRNA ZNNT1 in the tumorigenesis of UM by regulating key autophagy gene expression (Li et al., 2020). Therefore, there exist some interactions between lncRNAs and autophagy processes that could contribute to the pathogenesis of UM.

Autophagy is a highly conserved intracellular degradation system required for the maintenance of homeostasis and promotion of survival when facing cellular stress. It plays critical roles in age-related diseases including cancers (Mizushima, 2007). According to previous studies, autophagy was proved to be upregulated in UM and to be associated with poor prognosis (Giatromanolaki et al., 2011). Specific autophagy inhibitors might serve as a potential approach to suppress UM (Zhang et al., 2016). On the other hand, increased autophagy could be protective and cooperate with apoptosis to induce cell death in UM (Ambrosini et al., 2013). Hence, exploring autophagy-related gene signatures may deepen the understanding of the underlying molecular mechanism of UM and provide insights into its prognosis.

Based on the above background, we postulated that autophagy-related lncRNAs might help us identify high-risk UM patients and customize their personalized targeted treatment. This study aimed to analyze the correlation between the differentiated expressed lncRNAs and autophagy-related genes in UM and to screen for the most survival-relevant autophagy-related lncRNA signature by least absolute shrinkage and selection operator (LASSO) regression. Furthermore, univariate and multivariate Cox regression analyses were used to confirm its prognostic value in UM. This study may provide new insights into prognostic markers and deepen the understanding of autophagy-related lncRNA in UM.

Materials and Methods

Data Resource and Preprocessing

The RNA expression profiles of 80 UM patients were downloaded in FPKM format from The Cancer Genome Atlas (TCGA) database using the GDC Data Transfer Tool1, and the mRNA and lncRNA expression data were also acquired. The relevant clinicopathological characteristics including sex, age, and cancer stage were downloaded from the UCSC Xena website2. The detailed clinical characteristics are shown in Table 1. The flowchart of the analysis procedure is shown in Figure 1. The gene expression data of GSE22138 that contains 63 tissue samples of UM patients were downloaded from the Gene Expression Omnibus (GEO) as a validation dataset (Laurent et al., 2011).

TABLE 1
www.frontiersin.org

Table 1. Baseline data and correlations between risk score of lncRNA signature and clinicopathological characteristics in TCGA cohort.

FIGURE 1
www.frontiersin.org

Figure 1. The flowchart of the construction of the prognostic model of autophagy long non-coding RNA (lncRNA) in uveal melanoma (UM).

Identification of Differentially Expressed Genes and Autophagy-Related LncRNAs

In order to analyze the metastasis-related gene characteristics of UM patients, the differentially expressed genes (DEGs) between the metastasizing and non-metastasizing groups were calculated using the limma R package (Ritchie et al., 2015). DEGs with absolute log (fold change) > 0.05 and false discovery rate (FDR) < 0.05 were considered to be included for subsequent analysis. The result of differential analysis was presented in a heat map. Then, the autophagy-related gene set was downloaded from the MSigDB database3 (Subramanian et al., 2005), and the metastasis-related autophagy DEGs were extracted by the VennDiagram R package (Chen and Boutros, 2011). The expression correlations between those DEGs were calculated, and the protein–protein interaction (PPI) network was constructed by the Search Tool for the Retrieval of Interacting Genes (STRING) database4. The autophagy-related lncRNAs were identified based on the criteria that the correlation coefficient | R2| > 0.3 and P < 0.05.

Construction of Prognostic Autophagy-Related LncRNA Signature

The LASSO method was used to further analyze the prognostic autophagy-related lncRNAs for dimensionality reduction and feature selection. The coefficients obtained from the LASSO regression were weighted by each normalized gene expression value to construct the prognostic signature model, and the risk score of each UM patient was calculated by the formula as follows:

RiskScore = i Coefficient ( hubgene i ) * lncRNA Expression ( hubgene i )

Based on the median risk score, UM patients in the TCGA cohort were divided into high-risk and low-risk subgroups. Then, univariate and multivariate Cox regression analyses were performed to analyze the predictive ability of the signature model and patients’ clinical–pathological parameters on overall survival (OS) and progression-free interval (PFI). A prognostic nomogram integrating risk scores and clinical–pathological characteristics was constructed to analyze the probable 1-, 3-, and 5-year OS of the UM patients. Harrell’s concordance index (C-index) was applied to evaluate the predictive accuracy of the nomogram. To determine the clinical predictive validity of the nomogram, decision curve analysis (DCA) was used to quantify the net benefits at threshold probabilities in the dataset.

Functional Enrichment Analysis and Gene Set Enrichment Analysis

Gene Ontology (GO) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses were performed by ClusterProfiler R package (Yu et al., 2012). FDR < 0.05 was considered statistically significant. To investigate the differences in biological processes between the metastasizing and non-metastasizing groups, gene set enrichment analysis (GSEA) was performed between the two subgroups. The “c2.cp.kegg.v6.2.symbols.gmt” gene set was chosen as the reference gene list (Subramanian et al., 2005). We performed 1,000 times of permutations each time, and a pathway was considered to be significantly upregulated with FDR < 0.1. Furthermore, we analyzed the correlation between different subgroups and some related biologically processes constructed by Mariathasan et al. (2018), which included (1) immune checkpoints; (2) antigen processing processes; (3) CD8+ T-effector features; (4) epithelial–mesenchymal transition (EMT) markers, including EMT1, EMT2, and EMT3; (5) angiogenesis features; (6) pan-fibroblast TGF-β response signature (Pan-FTBRS); (7) WNT targets; (8) DNA damage repair; (9) mismatch repair; (10) nucleotide excision repair; (11) DNA replication; and (12) antigen processing and presentation.

Comparison of Immune Cell Infiltration Levels and Immune Scores Between the High-Risk and Low-Risk Subgroups of UM Patients

To quantify the proportion of different immune cells in the UM samples, the CIBERSORT algorithm and the LM22 gene set were used to analyze 22 human immune cell phenotypes in the tumor immune microenvironment (TME) (Newman et al., 2015). CIBERSORT is an algorithm for calculating the cell-type proportions in tumor samples, which uses a set of reference gene expression values as a representation (consisting of 547 marker genes). The Mann–Whitney U test was used to compare the infiltration level of immune cells between the two groups.

ESTIMATE (Estimation of STromal and Immune cells in MAlignant Tumor tissues using Expression data) is an algorithm for quantifying the immune infiltration level in a tumor sample based on gene expression data (Yoshihara et al., 2013). We calculated the immune and stromal scores of each UM sample using the ESTIMATE R package and compared the expression differences of human leukocyte antigen (HLA) family genes between the two groups.

Survival Analysis

The effects of the signature model on OS and PFI were compared. The Kaplan–Meier survival curve by survival package was used to perform survival analysis, and the log-rank test was used to assess the significance of the difference in survival time between the two groups.

Statistical Analysis

All data processing and analysis were performed using the R software (version 3.6.2). For the comparison of the two sets of continuous variables, the statistical significance of normally distributed variables was estimated by independent Student’s t-test, and differences between non-normally distributed variables were analyzed by Mann–Whitney U test. Categorical variables in the two groups were compared and analyzed using the chi-squared test or Fisher’s exact test. Pearson correlation analysis was used to calculate correlation coefficients between different genes. Univariate and multivariate Cox regression analyses were performed to identify independent prognostic factors. The pROC R package was used to visualize receiver-operating characteristic (ROC) curves, and the area under the curve (AUC) was calculated to assess the prognostic accuracy of the risk score (Robin et al., 2011). P < 0.05 was considered to be statistically significant.

Results

Identification of Metastasis-Associated Autophagy DEGs

To obtain metastasis-associated signature genes in UM patients, the gene expression profiles of metastasizing and non-metastasizing UM patients were analyzed and 1,107 DEGs were identified [| log (fold change)| > 0.5 and FDR < 0.05 filters, Figure 2A]. Meanwhile, the Venn diagram showed that there were 13 autophagy DEGs (Figure 2B). GO and KEGG analyses of all DEGs are shown in Supplementary Figure 1. GO and KEGG analyses of autophagy DEGs revealed that these DEGs were mainly enriched in cellular responses to starvation and extracellular stimulus and in relevant autophagy processes (Figures 2C,D).

FIGURE 2
www.frontiersin.org

Figure 2. Identification of autophagy differentially expressed genes (DEGs). (A) Heat map of the 1,107 DEGs between metastatic and non-metastatic UM. (B) Venn diagram of 13 autophagy DEGs. (C) GO functional enrichment analysis of autophagy DEGs. (D) KEGG pathway analysis of autophagy DEGs.

We further analyzed the interrelationship between metastasis-associated autophagy DEGs. As shown in Figures 3A,B, all 13 autophagy DEGs were differentially expressed with good concordance both in the TCGA and GEO cohorts. The correlation analysis showed significant correlations between the expression of most of these genes, with the highest correlation coefficient between RAF1 and PIK3R4 (Figure 3C). Meanwhile, the PPI network further revealed the interactions between 13 autophagy DEGs (Figure 3D).

FIGURE 3
www.frontiersin.org

Figure 3. Validation and correlation analysis of autophagy DEGs. (A) Box plot of the expression of the autophagy DEGs between metastatic and non-metastatic UM in TCGA cohort. (B) Box plot of the expression of the autophagy DEGs between metastatic and non-metastatic UM in the GEO cohort. (C) Correlation analysis diagram of autophagy DEGs. (D) Protein–protein interaction network of autophagy DEGs. *P < 0.05, **P < 0.01, ***P < 0.001.

Identification of Prognostic Autophagy-Related LncRNAs in UM

As shown in Figure 4A and Supplementary Figure 2, 27 autophagy-related lncRNAs were identified. The univariate Cox regression analysis showed that 14 lncRNAs had statistical significance, among which five lncRNAs were considered as protective factors and nine lncRNAs were considered as risk factors (Figure 4B). In order to confirm the impacts of these lncRNAs on the UM prognosis, all UM patients were ranked according to each lncRNA expression and subsequently divided into low and high groups (each group n = 40) according to media value. The survival prognostic analysis showed that 11 lncRNAs had statistically significant impacts on the prognosis of UM patient (Figure 4C).

FIGURE 4
www.frontiersin.org

Figure 4. Identification of autophagy-related lncRNAs in UM. (A) Heat map of the autophagy-related lncRNAs and their correlation with clinical characteristics. (B) Univariate Cox regression analysis results show that 14 autophagy-related lncRNAs correlate with overall survival of UM patients. (C) Kaplan–Meier survival curve analysis shows that 11 autophagy-related lncRNAs correlated with overall survival of UM patients.

Construction of Prognostic Autophagy-Related LncRNA Signature

To quantify the impact of autophagy-related lncRNAs on the prognosis of each UM patient, we constructed a prognostic signature by LASSO regression analysis. As shown in Figure 5A, six hub autophagy-related lncRNAs were included in the prognostic signature model. The risk score for each UM patient was calculated using the coefficients of each gene obtained from the LASSO algorithm. The distribution and scatter plot of risk scores, survival curves, and gene expression plots for each patient are shown in Figure 5B.

FIGURE 5
www.frontiersin.org

Figure 5. Establishment and validation of autophagy-related lncRNA signature. (A) LASSO coefficient profiles of lncRNA signature for overall survival (OS) in TCGA cohort. (B) Risk score distribution, survival status scatter plots, and heat map of the six prognostic lncRNAs in the high-risk and low-risk groups. (C) Kaplan–Meier analysis showed that UM patients in the high-risk group had poorer OS and PFI compared with UM patients in the low-risk group. (D) Time-dependent ROC curve of the lncRNA signature at 1, 3, and 5 years. (E) PCA of the lncRNA signature in high-risk and low-risk groups.

Kaplan–Meier analysis showed that the OS and PFI of patients with high-risk scores were significantly shorter than those with low-risk scores (Figure 5C). The ROC results validated the prognostic accuracy of the signature model (Figure 5D). The area under the ROC curves (AUCs) of 1, 3, and 5 years were 0.919, 0.845, and 0.724, which confirmed the accuracy of the signature in predicting survival times of UM patients. The principal component analysis (PCA) showed that the signature model could separate the distribution patterns between high- and low-risk groups in terms of the six autophagy-related lncRNAs and all genes (Figure 5E).

Clinical Correlation Analysis of the Signature and Nomogram Building

The clinical correlation analysis showed that high-risk patients demonstrated a higher risk of metastasis and poorer stage compared with low-risk patients, whereas no significant correlation with patient’s age was found (Figures 6A–C). The univariate and multivariate Cox regression analyses showed that the lncRNA signature was an independent prognostic indicator for UM patients (Table 2 and Figures 6D,E). Then we constructed a nomogram that integrates clinicopathological characteristics and the signature to predict OS of UM patient (Figure 6F). C-index was used to calculate the discriminative power of the nomogram, which showed a high degree of discrimination [0.93 (0.975–0.885)]. The DCA showed a prognostic benefit in 18–95% of patients using the signature model (Figure 6G).

FIGURE 6
www.frontiersin.org

Figure 6. Nomogram for predicting the overall survival probability of patients with UM. (A–C) The relationships between risk score of lncRNA signature and clinical parameters including age, metastasis, and stage. (D,E) Univariate and multivariate Cox regression analysis for risk score of lncRNA signature and other clinicopathological characteristics. (F) The nomogram was established based on lncRNA signature and other clinicopathological characteristics for predicting OS in UM patients at 1, 3, and 5 years. (G) DCA curves of the nomogram in UM patients at 3 years.

TABLE 2
www.frontiersin.org

Table 2. Univariate and multivariate Cox analysis of prognostic lncRNA signature and clinical features for OS in TCGA cohort.

Gene Set Enrichment Analysis and Immune Infiltration Analysis of the Signature Model

GSEA results showed that the proteasome, protein transport, metabolism, and apoptosis pathways were significantly enriched in the high-risk group compared with those in the low-risk group (Table 3 and Figure 7A). Meanwhile, some biologically relevant pathways, such as EMT, DNA replication, and cell cycle regulatory factors (Figure 7B), were significantly different between the high- and low-risk groups (P < 0.05).

TABLE 3
www.frontiersin.org

Table 3. Results of the GSEA base on prognostic lncRNA signature.

FIGURE 7
www.frontiersin.org

Figure 7. Functional enrichment analysis of lncRNA signature. (A) Six representative KEGG pathways in the high-risk group via GSEA. (B) Comparisons of different pathway features between the high-risk and low-risk groups. *P < 0.05, **P < 0.01, ***P < 0.001.

GSEA results also showed that the change in DEGs was significantly enriched in immune-related signaling pathways, including Toll-like signaling pathway, natural killer cell signaling pathway, etc (Figure 8A). Therefore, we further analyzed the immune infiltration between the high- and low-risk groups. As shown in Figures 8B,C, the high-risk patients had significantly higher immune scores and stromal scores compared with the low-risk patients. To further analyze the infiltration of immune cells in UM patients, the gene expression data of UM patients were downloaded from the TCGA cohort and 22 different specific immune cell infiltration in each patient were obtained using the CIBERSORT algorithm (Figure 8D). As shown in Figure 8E, most of the immune cells were significantly different between the two groups (Figure 8E). Furthermore, the gene expression of the HLA family was significantly different between the two groups (Figure 8F).

FIGURE 8
www.frontiersin.org

Figure 8. Correlations between lncRNA signature and immune cell infiltration in UM patients. (A) GSEA showed that various immune-related signaling pathways were enriched in the high-risk group compared with the low-risk group. (B,C) Comparisons of immune score and stromal score between the high-risk and low-risk groups. (D) Histogram of the 22 different specific immune cell infiltration in each UM samples. (E) Comparisons of various cell types in immune infiltration between the high-risk and the low-risk groups. (F) Comparisons of various HLA family genes between the high-risk and low-risk groups. *P < 0.05, **P < 0.01, ***P < 0.001.

Discussion

The mortality rate of UM patients with advanced stage is higher than 95% within 5 years, and the survival of patients with metastasis UM is below 1 year (Kujala et al., 2003), so finding novel and effective prognostic biomarkers is critical for UM. Although some epigenetic signatures like miRNA signature and DNA methylation signature have been established to precisely identify high-risk UM, the lncRNA signature for UM has not been reported. Experimental studies have shown that lncRNAs play a critical role in UM by regulating many cellular processes including autophagy. Therefore, it may be helpful to build an autophagy-related lncRNA signature to predict the prognosis of UM and improve therapeutic approaches.

In the present study, two datasets (TGGA and GSE22138) were enrolled to construct novel prognostic autophagy-related lncRNA signatures for UM. After comparing DEGs between metastasizing and non-metastasizing UM patients, we screened out 13 autophagy genes and further identified 27 autophagy-related lncRNAs through establishing the lncRNA–autophagy gene co-expression network. Moreover, we utilized the LASSO regression mode to build a six-hub autophagy-related lncRNA model that could differentiate high- and low-risk UM patients. The survival analyses, ROC curve, and PCA analyses confirmed the accuracy of this model. Through univariate and multivariate Cox regression analyses, we also revealed that the risk score of this model was an independent factor in predicting prognosis of UM patients. This six-lncRNA model may guide in the selection of more precise and individualized treatment strategies for UM patients.

Through bioinformatics analysis, we identified 13 metastasis-related autophagy DEGs, namely CDKN1A, DLC1, FOXO3, GABARAPL1, GRID1, ITGA6, MTMR14, PIK3R4, PRKCD, RAF1, SPHK1, TUSC1, and WIPI1. Among them, FOXO3, which acts as a transcription factor to correct autophagy perturbations, was reported to attenuate proliferation, invasion, and migration and promote cellular apoptosis of UM (Yan et al., 2016, 2018). PRKCD was shown to suppress autophagy and implicated as a link to activate MAPK in GNAQ mutant UM (Chen et al., 2017; Zhang et al., 2017). Meanwhile, we found that RAF1 and PIK3R4 showed the highest correlation coefficient. RAF1 acts as a regulatory link between the membrane-associated Ras GTPases and the MAPK/ERK pathway, and this critical regulatory link can lead to inhibition of cell proliferation and induction of apoptosis in UM cells (Samadi et al., 2012). PIK3R4 belongs to the regulatory subunit of the PI3K complex that is involved in the maturation of autophagosomes and endocytosis. A previous study revealed that somatic mutation of PIK3R4 contributed to melanoma metastasis (Shull et al., 2012), and its functions in UM progress are worthy of further study. It is worth mentioning that only 13 DEGs out of 1,107 genes in our results were found to be autophagy-related genes, and only 11 autophagy-related lncRNAs were associated with overall survival for UM. However, the autophagy-related lncRNA signature we built was proven to be highly effective for predicting the overall survival for UM. These results may be explained by the fact that autophagy takes part in a certain but crucial phase of the metastatic cascade, but other complicated mechanisms take up more roles in regulating UM metastasis.

Autophagy is a lysosomal catabolic process that is regulated by non-coding RNAs such as miRNAs and lncRNAs. Experimental and clinical studies revealed that autophagy-related lncRNAs had significant diagnostic and therapeutic values in various cancers including UM. For example, a recent research reported that lncRNA ZNNT1 could act as a tumor inhibitor in UM by upregulating ATG12 gene expression, indicating delicate interactions between lncRNAs and autophagy (Li et al., 2020). Our model enrolled six autophagy-related lncRNAs including PVT1 that predicted overall survival in UM. PVT1 has been identified as a prognostic biomarker and played a critical role in the development of pancreatic ductal adenocarcinoma via triggering cytoprotective autophagy (Huang et al., 2018), and it affected EMT and cell proliferation and migration via modulating CDKN1A expression in breast cancer (Wang et al., 2018). Furthermore, PVT1 also took part in drug resistance of pancreatic cancer by interacting with ATG14 and promoting PtdIns3K activity (Zhou et al., 2020). In UM, silencing PVT1 could suppress tumorigenesis by activating the miR-17-3p-dependent p53 signaling pathway (Wu et al., 2019). Further in-depth studies are needed to investigate the molecular mechanism between PVT1 and autophagy in UM.

In order to explore the biological functions of the prognostic six-lncRNA model in UM, GSEA analysis was performed to compare the low- and high-risk groups. The results demonstrated that these lncRNAs might exert influence on UM progression through the immune-related pathway. LncRNAs are involved in immunological processes in cancers by mediating the functional state of immune cells and related pathways. Based on this lncRNA model, we found that the immune scores and stromal scores in high-risk UM were higher than in low-risk UM, and various members of the HLA family showed statistical difference between the two groups of UM. Several studies have shown an abnormal expression of the HLA family in UM and confirmed that a high expression of HLA class I is a risk factor for the development of metastasis (Gezgin et al., 2017; Souri et al., 2019; Wierenga et al., 2019). Recent studies showed that the lncRNA MIR155 host gene (MIR155HG) can be a potential prognostic marker for predicting the curative effect of immune checkpoint inhibitor therapy in multiple cancers including UM (Peng et al., 2019). Our results, combined with these publications, indicated that the six autophagy-related lncRNAs might play crucial roles in the pathogenesis of UM by regulating related genes to affect the immune responses, and there existed close relationships between autophagy and immune regulation in UM.

This study was not without limitations. All results of this study are based on bioinformatics analysis, and the prognostic signature enrolled some autophagy-related lncRNAs that no research confirmed their roles in UM so far. Thus, further experiments combining biochemical experiments and clinical prognostic data are required to confirm how these lncRNAs affect the prognosis of UM patients through autophagy exactly.

Conclusion

In conclusion, this study screened and identified 13 dysregulated autophagy-related genes in UM, and a prognostic risk signature enrolling the six dysregulated autophagy-related lncRNAs was constructed, which was proven to be a useful marker for predicting UM patients’ prognosis with high efficacy. Based on the prognostic lncRNA model, dysregulated immune-related pathway, especially HLA expression, may contribute to the poor prognosis in high-risk UM. This signature model could be deployed to assist clinicians to identify high-risk UM patients and help scientists to explore the molecular mechanism of autophagy-related lncRNA in UM pathogenesis.

Data Availability Statement

All datasets generated for this study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author/s.

Author Contributions

YC and NX: conceptualization. YC: methodology and writing—original draft preparation. YC and MZ: formal analysis and investigation. NX and JC: writing—review and editing. NX, JC, and YC: funding acquisition. NX: supervision. All authors contributed to the article and approved the submitted version.

Funding

This research was funded by the Fujian Provincial Natural Science Foundation (Grant No. 2019J01152), the Joint Funds for the Innovation of Science and Technology of Fujian Province (Grant No. 2018Y9035), the Startup Fund for Scientific Research of Fujian Medical University (Grant Nos. 2019QH1149 and 2019QH1178), and the Fujian Provincial Health Technology Project (Grant Nos. 2019-1-31, 2020GGB003, and 2020QNA002).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2021.625583/full#supplementary-material

Supplementary Figure 1 | GO and KEGG analyses of all DEGs. (A) GO functional enrichment analysis of all DEGs. (B) KEGG pathways analysis of all DEGs.

Supplementary Figure 2 | Correlation analysis diagram between autophagy-related genes and lncRNAs.

Footnotes

  1. ^ https://portal.gdc.cancer.gov/
  2. ^ http://xena.ucsc.edu/
  3. ^ https://www.gsea-msigdb.org/gsea/msigdb
  4. ^ https://www.string-db.org/

References

Amaro, A., Gangemi, R., Piaggio, F., Angelini, G., Barisione, G., Ferrini, S., et al. (2017). The biology of uveal melanoma. Cancer Metastasis Rev. 36, 109–140.

Google Scholar

Ambrosini, G., Musi, E., Ho, A. L., De Stanchina, E., and Schwartz, G. K. (2013). Inhibition of mutant GNAQ signaling in uveal melanoma induces AMPK-dependent autophagic cell death. Mol. Cancer Ther. 12, 768–776. doi: 10.1158/1535-7163.mct-12-1020

PubMed Abstract | CrossRef Full Text | Google Scholar

Cerman, E., and Cekic, O. (2015). Clinical use of photodynamic therapy in ocular tumors. Surv. Ophthalmol. 60, 557–574. doi: 10.1016/j.survophthal.2015.05.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, H., and Boutros, P. C. (2011). VennDiagram: a package for the generation of highly-customizable Venn and Euler diagrams in R. BMC Bioinformatics 12:35. doi: 10.1186/1471-2105-12-35

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, X., Wu, Q., Depeille, P., Chen, P., Thornton, S., Kalirai, H., et al. (2017). RasGRP3 mediates MAPK pathway activation in GNAQ mutant uveal melanoma. Cancer Cell 31, 685–696.e6.

Google Scholar

Chew, A. L., Spilsbury, K., and Isaacs, T. W. (2015). Survival from uveal melanoma in Western Australia 1981-2005. Clin. Exp. Ophthalmol. 43, 422–428. doi: 10.1111/ceo.12490

PubMed Abstract | CrossRef Full Text | Google Scholar

Collaborative Ocular Melanoma Study Group (2006). The COMS randomized trial of iodine 125 brachytherapy for choroidal melanoma: V. Twelve-year mortality rates and prognostic factors: COMS report No. 28. Arch. Ophthalmol. 124, 1684–1693. doi: 10.1001/archopht.124.12.1684

PubMed Abstract | CrossRef Full Text | Google Scholar

Correa, Z. M., and Augsburger, J. J. (2016). Independent prognostic significance of gene expression profile class and largest basal diameter of posterior uveal melanomas. Am. J. Ophthalmol. 162, 20–27.e1.

Google Scholar

Gezgin, G., Luk, S. J., Cao, J., Dogrusoz, M., Van Der Steen, D. M., Hagedoorn, R. S., et al. (2017). PRAME as a potential target for immunotherapy in metastatic uveal melanoma. JAMA Ophthalmol. 135, 541–549. doi: 10.1001/jamaophthalmol.2017.0729

PubMed Abstract | CrossRef Full Text | Google Scholar

Giatromanolaki, A. N., Charitoudis, G. S., Bechrakis, N. E., Kozobolis, V. P., Koukourakis, M. I., Foerster, M. H., et al. (2011). Autophagy patterns and prognosis in uveal melanomas. Mod. Pathol. 24, 1036–1045. doi: 10.1038/modpathol.2011.63

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, F., Chen, W., Peng, J., Li, Y., Zhuang, Y., Zhu, Z., et al. (2018). LncRNA PVT1 triggers Cyto-protective autophagy and promotes pancreatic ductal adenocarcinoma development via the miR-20a-5p/ULK1 Axis. Mol. Cancer 17:98.

Google Scholar

Kaliki, S., and Shields, C. L. (2017). Uveal melanoma: relatively rare but deadly cancer. Eye (Lond.) 31, 241–257. doi: 10.1038/eye.2016.275

PubMed Abstract | CrossRef Full Text | Google Scholar

Kujala, E., Makitie, T., and Kivela, T. (2003). Very long-term prognosis of patients with malignant uveal melanoma. Invest. Ophthalmol. Vis. Sci. 44, 4651–4659. doi: 10.1167/iovs.03-0538

PubMed Abstract | CrossRef Full Text | Google Scholar

Laurent, C., Valet, F., Planque, N., Silveri, L., Maacha, S., Anezo, O., et al. (2011). High PTP4A3 phosphatase expression correlates with metastatic risk in uveal melanoma patients. Cancer Res. 71, 666–674. doi: 10.1158/0008-5472.can-10-0605

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, P., He, J., Yang, Z., Ge, S., Zhang, H., Zhong, Q., et al. (2020). ZNNT1 long noncoding RNA induces autophagy to inhibit tumorigenesis of uveal melanoma by regulating key autophagy gene expression. Autophagy 16, 1186–1199. doi: 10.1080/15548627.2019.1659614

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, Y., Yang, X., Yang, J., Wang, H., and Wei, W. (2018). An 11-gene-based prognostic signature for uveal melanoma metastasis based on gene expression and DNA methylation profile. J. Cell. Biochem. doi: 10.1002/jcb.28151 [Epub ahead of print]

PubMed Abstract | CrossRef Full Text | Google Scholar

Mariathasan, S., Turley, S. J., Nickles, D., Castiglioni, A., Yuen, K., Wang, Y., et al. (2018). TGFbeta attenuates tumour response to PD-L1 blockade by contributing to exclusion of T cells. Nature 554, 544–548.

Google Scholar

Mashayekhi, A., Shields, C. L., Rishi, P., Atalay, H. T., Pellegrini, M., Mclaughlin, J. P., et al. (2015). Primary transpupillary thermotherapy for choroidal melanoma in 391 cases: importance of risk factors in tumor control. Ophthalmology 122, 600–609.

Google Scholar

Mizushima, N. (2007). Autophagy: process and function. Genes Dev. 21, 2861–2873. doi: 10.1101/gad.1599207

PubMed Abstract | CrossRef Full Text | Google Scholar

Newman, A. M., Liu, C. L., Green, M. R., Gentles, A. J., Feng, W., Xu, Y., et al. (2015). Robust enumeration of cell subsets from tissue expression profiles. Nat. Methods 12, 453–457. doi: 10.1038/nmeth.3337

PubMed Abstract | CrossRef Full Text | Google Scholar

Peng, L., Chen, Z., Chen, Y., Wang, X., and Tang, N. (2019). MIR155HG is a prognostic biomarker and associated with immune infiltration and immune checkpoint molecules expression in multiple cancers. Cancer Med. 8, 7161–7173. doi: 10.1002/cam4.2583

PubMed Abstract | CrossRef Full Text | Google Scholar

Ritchie, M. E., Phipson, B., Wu, D., Hu, Y., Law, C. W., Shi, W., et al. (2015). limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 43:e47. doi: 10.1093/nar/gkv007

PubMed Abstract | CrossRef Full Text | Google Scholar

Robin, X., Turck, N., Hainard, A., Tiberti, N., Lisacek, F., Sanchez, J. C., et al. (2011). pROC: an open-source package for R and S+ to analyze and compare ROC curves. BMC Bioinformatics 12:77. doi: 10.1186/1471-2105-12-77

PubMed Abstract | CrossRef Full Text | Google Scholar

Samadi, A. K., Cohen, S. M., Mukerji, R., Chaguturu, V., Zhang, X., Timmermann, B. N., et al. (2012). Natural withanolide withaferin a induces apoptosis in uveal melanoma cells by suppression of Akt and c-MET activation. Tumour Biol. 33, 1179–1189. doi: 10.1007/s13277-012-0363-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Shields, J. A., and Shields, C. L. (2015). Management of posterior uveal melanoma: past, present, and future: the 2014 Charles L. Schepens lecture. Ophthalmology 122, 414–428. doi: 10.1016/j.ophtha.2014.08.046

PubMed Abstract | CrossRef Full Text | Google Scholar

Shull, A. Y., Latham-Schwark, A., Ramasamy, P., Leskoske, K., Oroian, D., Birtwistle, M. R., et al. (2012). Novel somatic mutations to PI3K pathway genes in metastatic melanoma. PLoS One 7:e43369. doi: 10.1371/journal.pone.0043369

PubMed Abstract | CrossRef Full Text | Google Scholar

Souri, Z., Wierenga, A. P. A., Mulder, A., Jochemsen, A. G., and Jager, M. J. (2019). HLA expression in uveal melanoma: an indicator of malignancy and a modifiable immunological target. Cancers (Basel) 11:1132. doi: 10.3390/cancers11081132

PubMed Abstract | CrossRef Full Text | Google Scholar

Subramanian, A., Tamayo, P., Mootha, V. K., Mukherjee, S., Ebert, B. L., Gillette, M. A., et al. (2005). Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc. Natl. Acad. Sci. U.S.A. 102, 15545–15550. doi: 10.1073/pnas.0506580102

PubMed Abstract | CrossRef Full Text | Google Scholar

Wan, Q., Sang, X., Jin, L., and Wang, Z. (2020). Alternative splicing events as indicators for the prognosis of uveal melanoma. Genes (Basel) 11:227. doi: 10.3390/genes11020227

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, L., Wang, R., Ye, Z., Wang, Y., Li, X., Chen, W., et al. (2018). PVT1 affects EMT and cell proliferation and migration via regulating p21 in triple-negative breast cancer cells cultured with mature adipogenic medium. Acta Biochim. Biophys. Sin. (Shanghai) 50, 1211–1218. doi: 10.1093/abbs/gmy129

PubMed Abstract | CrossRef Full Text | Google Scholar

Wierenga, A. P. A., Gezgin, G., Van Beelen, E., Eikmans, M., Spruyt-Gerritse, M., Brouwer, N. J., et al. (2019). Soluble HLA in the aqueous humour of uveal melanoma is associated with unfavourable tumour characteristics. Cancers (Basel) 11:1202. doi: 10.3390/cancers11081202

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, S., Chen, H., Han, N., Zhang, C., and Yan, H. (2019). Long noncoding RNA PVT1 silencing prevents the development of uveal melanoma by impairing MicroRNA-17-3p-dependent MDM2 upregulation. Invest. Ophthalmol. Vis. Sci. 60, 4904–4914. doi: 10.1167/iovs.19-27704

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, S., Chen, H., Zuo, L., Jiang, H., and Yan, H. (2020). Suppression of long noncoding RNA MALAT1 inhibits the development of uveal melanoma via microRNA-608-mediated inhibition of HOXC4. Am. J. Physiol. Cell Physiol. 318, C903–C912.

Google Scholar

Xin, X., Zhang, Y., Ling, F., Wang, L., Sheng, X., Qin, L., et al. (2019). Identification of a nine-miRNA signature for the prognosis of Uveal Melanoma. Exp. Eye Res. 180, 242–249. doi: 10.1016/j.exer.2019.01.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Xing, Y., Wen, X., Ding, X., Fan, J., Chai, P., Jia, R., et al. (2017). CANT1 lncRNA triggers efficient therapeutic efficacy by correcting aberrant lncing cascade in malignant uveal melanoma. Mol. Ther. 25, 1209–1221. doi: 10.1016/j.ymthe.2017.02.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Yan, F., Liao, R., Farhan, M., Wang, T., Chen, J., Wang, Z., et al. (2016). Elucidating the role of the FoxO3a transcription factor in the IGF-1-induced migration and invasion of uveal melanoma cancer cells. Biomed. Pharmacother. 84, 1538–1550. doi: 10.1016/j.biopha.2016.11.027

PubMed Abstract | CrossRef Full Text | Google Scholar

Yan, F., Liao, R., Lin, S., Deng, X., Little, P. J., and Zheng, W. (2018). Forkhead box protein O3 suppresses uveal melanoma development by increasing the expression of Bcl2like protein 11 and cyclindependent kinase inhibitor 1B. Mol. Med. Rep. 17, 3109–3114.

Google Scholar

Yoshihara, K., Shahmoradgoli, M., Martinez, E., Vegesna, R., Kim, H., Torres-Garcia, W., et al. (2013). Inferring tumour purity and stromal and immune cell admixture from expression data. Nat. Commun. 4: 2612.

Google Scholar

Yu, G., Wang, L. G., Han, Y., and He, Q. Y. (2012). clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS 16, 284–287. doi: 10.1089/omi.2011.0118

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, D., Xu, X., and Dong, Z. (2017). PRKCD/PKCdelta contributes to nephrotoxicity during cisplatin chemotherapy by suppressing autophagy. Autophagy 13, 631–632. doi: 10.1080/15548627.2016.1269990

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, Y., Song, H., Guo, T., Zhu, Y., Tang, H., Qi, Z., et al. (2016). Overexpression of Annexin II receptor-induced autophagy protects against apoptosis in uveal melanoma cells. Cancer Biother. Radiopharm. 31, 145–151. doi: 10.1089/cbr.2016.1991

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, Y., Yang, Y., Chen, L., and Zhang, J. (2014). Expression analysis of genes and pathways associated with liver metastases of the uveal melanoma. BMC Med. Genet. 15:29. doi: 10.1186/1471-2350-15-29

PubMed Abstract | CrossRef Full Text | Google Scholar

Zheng, X., Tang, H., Zhao, X., Sun, Y., Jiang, Y., and Liu, Y. (2017). Long non-coding RNA FTH1P3 facilitates uveal melanoma cell growth and invasion through miR-224-5p. PLoS One 12:e0184746. doi: 10.1371/journal.pone.0184746

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, C., Yi, C., Yi, Y., Qin, W., Yan, Y., Dong, X., et al. (2020). LncRNA PVT1 promotes gemcitabine resistance of pancreatic cancer via activating Wnt/beta-catenin and autophagy pathway through modulating the miR-619-5p/Pygo2 and miR-619-5p/ATG14 axes. Mol. Cancer 19:118.

Google Scholar

Keywords: long non-coding RNAs, autophagy, uveal melanoma, prognostic signature, The Cancer Genome Atlas, Gene Expression Omnibus

Citation: Cui Y, Zheng M, Chen J and Xu N (2021) Autophagy-Related Long Non-coding RNA Signature as Indicators for the Prognosis of Uveal Melanoma. Front. Genet. 12:625583. doi: 10.3389/fgene.2021.625583

Received: 03 November 2020; Accepted: 25 February 2021;
Published: 01 April 2021.

Edited by:

Shunbin Xiong, University of Texas MD Anderson Cancer Center, United States

Reviewed by:

Xiaojie Yu, University of Texas MD Anderson Cancer Center, United States
Xin Hu, Fudan University, China

Copyright © 2021 Cui, Zheng, Chen and Xu. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Nuo Xu, ophnuox@tmu.edu.cn

These authors have contributed equally to this work

Disclaimer: All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.