Development and validation of an m6A RNA methylation regulator-based signature for the prediction of prognosis and immunotherapy in cutaneous melanoma.

To the Editor: The incidence of cutaneous melanoma, which has an extremely poor prognosis and is responsible for the majority of skin cancer-related deaths, has been increasing rapidly.[1] m6A (N6-methyladenosine), which is methylated at the N6 position of adenosine, is the most prominent modification in messenger ribonucleic acid (mRNAs) and long non-coding ribonucleic acid (RNAs) in eukaryotic cells and modifies the RNA's stability, translocation, RNA splicing and translation into protein. Multiple studies have revealed that tumorigenesis is cogently associated with m6A RNA methylation.[2] However, the relationship between m6A methylation regulators and the prognosis, as well as the risk of immune evasion of cutaneous melanoma (CM), remains unknown. Our study aims to identify the prognostic value of m6A regulators in CM. We used data from the The Cancer Genome Atlas (TCGA) and Genotype-Tissue Expression databases and included samples of 471 CM patients and 813 normal control samples for further analysis. Differential expression of 21 m6A regulators between melanoma tissues and normal skin tissues and the interactions of these 21 m6A regulators were analyzed [Supplementary Figure1, https://links.lww.com/CM9/B430]. In the tumor samples, the correlation between the expression levels of m6A regulators and overall survival (OS) was revealed by implementing univariate Cox regression analysis. The results demonstrated that RNA binding motif protein 15B (RBM15B), methyltransferase-like protein 16 (METTL16), and Wilms tumor 1 associated protein (WTAP) were significantly related to OS. Except for WTAP, both RBM15B and METTL16 were risk regulators with hazard rates (HRs) > 1[Figure 1A]. Then, these three m6A regulators were chosen to build the prognostic signature, and the least absolute shrinkage and selection operator (LASSO) algorithm was used to determine the coefficients of these three regulators [Figure 1B, C]. The risk scores for each CM patient were generated by the following formula: risk score = 0.4295 × RBM15B + 0.3191 × METTL16 – 0.3317 × WTAP. Based on the median cutoff value of the risk scores, a total of 471 melanoma patients were grouped into the low-risk (green) and high-risk groups (red) [Figure 1D]. Finally, the survival analysis showed that the melanoma patients in the high-risk group had a significantly shorter OS than those in the low-risk group [Figure 1E]. This signature for predicting OS in CM patients was further successfully validated in the GEO database, which contained 55 CM cases [Figure 1F].Figure 1: Construction and validation of an m6A methylation regulator-based signature for prediction of prognosis and immunotherapy in cutaneous melanoma.Then, the univariate [Figure 1G] and multivariate [Figure 1H] Cox regression analyses indicated that clinicopathological features, such as T stage, N stage, MKI67 and HMB45 were identified as independent risk indicators for CM in the TCGA database. Surprisingly, the predictive power of our three-gene risk score appears to be more sensitive. By using the TIMER algorithm and ggpubr package, we also found a significant negative correlation between the risk score and infiltration levels of dendritic cells, neutrophils, macrophages, CD8+ T cells, CD4+ T cells and B cells [Supplementary Figure 2, https://links.lww.com/CM9/B430]. The high-risk group was associated with higher expression of several immune checkpoints, including programmed death 1 (PD-1), programmed death-ligand 1 (PD-L1), T cell immunoglobulin and ITIM domain (TIGIT), cytotoxic T-lymphocyte–associated antigen 4 (CTLA-4), V-set immunoregulatory receptor (VSIR) and lymphocyte activation gene 3 (LAG3) [Supplementary Figure 3, https://links.lww.com/CM9/B430]. The results above indicated that a high m6A risk score correlated with the immune microenvironment, which urged us to consider whether the m6A risk score could act as a biomarker to predict the response rate to immunotherapy. To test our hypothesis, tumor patients from the GSE78220 and GSE100797 cohorts who received anti-PD-1 monotherapy were divided into the high and low m6A risk score groups using the median value as the cutoff value. In the GSE78220 cohort, the boxplots [Figure 1I] and bars [Figure 1J] suggested that patients with a lower m6A risk score tended to have a higher response to anti-PD-1 monotherapy. The receiver operating characteristic (ROC) curve also indicated that the m6A risk score had high accuracy in predicting the response to immunotherapy [Figure 1K]. Similar results were observed in the GSE100797 cohort [Figure 1L–N]. In addition, the m6A risk score is more predictive of the response to immunotherapy than PD-L1[Supplementary Figure 4, https://links.lww.com/CM9/B430]. Collectively, the results demonstrated that the m6A risk score could be a biomarker for predicting the response to immunotherapy. In this study, according to the TCGA database, we found that all 21 m6A regulators were abnormally expressed in CM patients. A three-gene risk signature including RBM15B, METTL16 and WTAP was generated and served as a credible factor for predicting the survival time of CM patients, which was successfully validated in the GEO database. WTAP was involved in promoting m6A installation by recruiting METTL3 and METTL14 into nuclear speckles. The level of m6A can be decreased significantly by silencing WTAP. It has been reported that WTAP is upregulated in tumors such as acute myelocytic leukemia. However, there is little information available on the role of WTAP in CM. According to our prognostic model, WTAP had a positive correlation with OS in CM, serving as a protective prognostic factor, which is consistent with the results of a previous study.[3] RBM15B, a paralog of RBM15, was reported to combine METTL3 and WTAP and guide them into specific RNA sites for m6A modification. Wang T's2 study showed that RBM15B could act as an independent prognostic biomarker and that it had clinical significance in uveal melanoma. METTL16, a homolog of METTL3, can N6-methylate a conserved hairpin of MAT2A mRNA and U6 spliceosomal small nuclear RNA with stimulation of S-adenosylmethionine.[4] However, the relationship between these two regulators and CM is still unclear. In our research, the expression of these two regulators had a negative association with the prognosis of CM, indicating that RBM15B and METTL16 might be promoters of CM tumorigenesis. Multiple studies have revealed the relationship between m6A and tumor immunity, especially regarding immunotherapy. Our findings revealed that the CM patients in the high-risk group had lower infiltration levels of dendritic cells, neutrophils, macrophages, CD8+ T cells, CD4+ T cells and B cells and higher expression levels of PD-1, PD-L1, TIGIT, CTLA-4, VSIR, and LAG3. Notably, patients with lower m6A risk scores were found to respond better to anti-PD-1 monotherapy than those with higher m6A risk scores. PD-L1 was approved by the United States Food and Drug Administration as a complementary diagnostic test for patients who are candidates for anti-PD-1 therapy,[5] but our data showed that the m6A risk score signature that we first established was more predictive of the response to anti-PD-1 therapy than PD-L1. All these findings revealed that our prognostic model can better predict the clinical outcome and the response to immunotherapy of CM patients and may provide significant guidance for immunotherapy strategies in CM. Conflicts of interest None.