Systematic analysis of the expression profile and prognostic significance of m6A regulators and PD-L1 in hepatocellular carcinoma

DOI: https://doi.org/10.21203/rs.3.rs-1992254/v1

Abstract

Background

Hepatocellular carcinoma (HCC) is one of the major contributors to cancer burden worldwide, and its incidence has increased in recent decades. N6-methyladenosine (m6A) modification plays an important role in the tumorigenesis and metastasis of HCC by changing gene expression and function at multiple levels, such as ribonucleic acid (RNA) splicing, stability, translocation and translation.

Methods

The ONCOMINE, UALCAN, GEPIA, Kaplan–Meier plotter, cBioPortal for Cancer Genomics, STRING and TIMER2 databases were used for bioinformatic analyses. Quantitative polymerase chain reaction and western blotting were used to detect the expression of m6A RNA methylation regulators in HCC tissues.

Results

The transcription of m6A RNA methylation regulators was upregulated in patients with HCC, and overexpression of YTHDF1/2, YTHDC1, RBM15 and METTL3 was significantly correlated with clinical stages of HCC. In addition, downregulation of ZC3H13 and METTL14 and upregulation of other m6A RNA methylation regulators were associated with a poor prognosis. A high mutation rate (89%) of m6A RNA methylation regulators was also observed in patients with HCC, and mutations in methylation regulators were associated with poor overall survival and disease-free survival. Finally, the expression of the YTHDF family was significantly associated with immune infiltration in the HCC microenvironment.

Conclusion

m6A RNA methylation regulators and programmed death-ligand 1 may play an important role in the tumorigenesis and immune invasion and escape of HCC and may be risk factors affecting the survival of patients with HCC.

Introduction

Primary liver cancer is the seventh most common cancer and the second most common cause of cancer mortality worldwide[1, 2], with Asia and Africa having the highest incidence rates[3, 4]. Hepatocellular carcinoma (HCC) is the major type of liver cancer worldwide, accounting for approximately 75% of all liver cancer cases[5, 6], with a 5-year survival rate of only 18%, making it the second most life-threatening cancer after pancreatic cancer[2, 7]. However, biomarkers for early detection and prognostication of HCC and for the prediction and monitoring of treatment responses should be identified[8, 9]. In recent years, better immunotherapy outcomes have been observed in patients with advanced HCC owing to its safety and fewer adverse events. Therefore, immunotherapy can be considered a new strategy for clinical management. However, the overall disease control rate of HCC and its treatment strategies should be improved because only a few patients benefit from immunotherapy. Immune system imbalance plays a critical role in HCC development and is characterised by immunosuppression[10]. Therefore, regulatory mechanisms of the tumour immune microenvironment (TIME) should be further examined to identify effective biomarkers that can accurately predict the prognosis and improve and optimise individualised immunotherapy management.

As the main form of ribonucleic acid (RNA) modification, N6-methyladenosine (m6A) plays a critical role in various biological processes and has been increasingly recognised in recent years[11, 12]. m6A modification is dynamic and reversible and is mainly mediated by methyltransferases and demethylases[12, 13]. It can be identified according to different reader proteins and performed based on its unique biological functions[14]. The biological effects of m6A are mainly mediated by writer, eraser and reader proteins. Because m6A modification is dynamic and reversible, RNA can be methylated by writers and demethylated by erasers[15]. In addition to writers and erasers, another group known as readers can determine m6A modifications and perform different biological functions[16]. ‘Writers’ mainly include KIAA1429, METTL3, RBM15, ZC3H13, METTL16 and METTL14 and their cofactor WTAP[17]. ‘Readers’ mainly include HNRNPC, YTHDC1, YTHDC2, YTHDF1, YTHDF2 and YTHDF3[18]. ‘Erasers’ primarily include ALKBH5 and FTO[12]. Abnormal expression of m6A regulatory factors plays an important role in tumour progression, prognosis and radiation resistance[19].

To date, most studies on the mechanisms of immune evasion in HCC have focused on the programmed death receptor 1 (PD-1/PDCD1)/programmed death ligand 1 (PD-L1) pathway. PDCD1 is an immunoinhibitory receptor expressed on the surface of activated T and B cells and natural killer (NK) T cells. PD-L1, the major ligand of PDCD1, is expressed on the surface of various immune cells such as antigen-presenting and endothelial cells. PDCD1/PD-L1 ligation suppresses the activation of immune cells and the production of certain cytokines such as IFN-γ, thereby inducing immune suppression and peripheral tolerance[20].

Although the study of m6A methylation in HCC is continuously progressing, its potential role in TIME remains unclear. Han et al. demonstrated that m6A methylation prolonged novel antigen-specific immunity through YTHDF1, with significantly higher CD8 + T and NK cell levels in YTHDF1-deficient mice than that in wild-type mice, to enhance antitumour responses. This finding suggests that YTHDF1 is an important mediator of tumour immune escape[21]. Therefore, m6A regulatory factors involved in the tumour immune response pathway may be a promising target for enhancing the clinical immunotherapeutic response.

In this study, bioinformatic methods were used to systematically evaluate the relationship between the expression of m6A RNA methylation regulators and programmed death-ligand 1 (PD-L1) with the prognosis, mutations and tumour stages of HCC. In addition, related characteristics of m6A regulatory factors in immune cell infiltration were analysed. The function and pathways of m6A regulators and their frequently mutated genes were also investigated. The results showed that abnormal expression of m6A modification regulators was related to the tumour stage and prognosis of patients with HCC and was involved in the regulation of immune invasion in the HCC microenvironment.

Materials And Methods

ONCOMINE (http://www.oncomine.org):

The ONCOMINE database was used to analyse the transcription levels of m6A RNA methylation and PD-L1 in HCC tumour tissues and the corresponding adjacent healthy control tissues. The cell colour was determined by the best gene rank percentile for analysis within the cell, and a Q–Q graph and histogram were used to detect whether the sample data followed a normal distribution. Subsequently, the Student’s t-test was used to generate p-values. A p-value of 0.05, a fold change of 2 and a gene rank in the top 10% were set as the thresholds of significance.

UALCAN (http://ualcan.path.uab.edu/):

UALCAN is a clinical database that is used to analyse gene expression among tumour subgroups, survival and cancer genome-based mapping (TCGA) of level 3 RNA-seq data across 31 cancer types [22]. In this study, UALCAN was used to examine the distinct expression levels of tumour and healthy tissues. TCGA database was used for analysing the expression of m6A regulatory factors and PD-L1 in HCC, and the data were downloaded. The Student’s t-test was used to generate p-values, with a p-value cutoff of 0.05.

GEPIA (http://gepia.cancer-pku.cn/):

GEPIA is a web server for cancer and normal gene expression profiling and interaction analyses [23]. In this study, the ‘Multiple Gene Comparison’ module was used to evaluate the expression of multiple m6A RNA methylation regulators. The ‘Expression DIY’ module was used to examine the relationship between 15 m6A RNA methylation regulators and clinicopathological parameters. The relationship between the expression of m6A regulators and tumour stages was analysed using the ‘Single Gene Analysis’ module. The Student’s t-test was used to generate p-values, with a p-value cutoff of 0.05.

Kaplan–Meier Plotter (http://kmplot.com/analysis/):

The Kaplan–Meier (K-M) plotter can be used to assess the effects of approximately 54,000 genes on survival across 21 cancer types [24]. In this study, K-M curves were used to analyse the relationship between gene mutations in m6A RNA methylation regulators and PD-L1 and overall survival (OS), progression-free survival (PFS), recurrence-free survival (RFS) and disease-specific survival (DSS) in patients with HCC. The log-rank test was used to evaluate significant differences in the survival curves. A p-value of <0.05 was considered statistically significant. In addition, the K-M plotter was used to evaluate the prognostic value of m6A RNA methylation regulators and PD-L1 mRNA expression; patients with HCC were divided into the high- and low-expression groups based on the median values of mRNA expression, and the results were validated using K-M survival curves based on the hazard ratio (HR) with 95% confidence intervals (CIs) and log-rank p-values. Patients were grouped to select the best cutoff, and survival analysis was performed, including OS, FP and PPS analyses. A p-value of <0.05 was considered statistically significant.

TCGA Data and cBioPortal for Cancer Genomics (http://www.cbioportal.org):

TCGA database was used to analyse the genomic profiles of m6A RNA methylation regulators and PD-L1, including mutations, putative copy number alterations from GISTIC and mRNA expression z-scores (RNA-Seq V2 RSEM) with a z-score threshold of ±1.8. The co-expression of m6A regulators and PD-L1 was analysed using the ‘co-expression’ module of cBioPortal. Pearson correlation coefficient was used to examine the correlation among co-expressed genes of m6A regulators and PD-L1. The top 10 co-expressed genes with the highest Pearson correlation coefficients of each m6A RNA methylation regulator and PD-L1 were screened. In the cBioPortal database, ‘Liver Hepatocellular Carcinoma (TCGA, Firehose Legacy)’ was selected in the Liver section, and the ‘Query By Gene’ function was used for analysis. The ‘mRNA expression z-scores relative to diploid samples (RNA Seq V2 RSEM)’ and ‘Samples with mRNA data (RNA Seq V2)’ modules were selected for genomic profiling. The correlation between the expression of m6A regulators and PD-L1 was analysed using the ‘Mutual Exclusivity’ module.

STRING (http://www.string-db.org):

STRING is a database of known and predicted protein–protein interactions (PPIs)[25]. In this study, STRING was used to analyse PPI networks of m6A RNA methylation regulators and PD-L1. Homo sapiens was selected as the species, and a combined score of >0.7 was considered statistically significant. The nodes indicated proteins, the edges indicated the interaction of proteins and disconnected nodes in the network were hidden.

DAVID (https://david.ncifcrf.gov/summary.jsp):

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses were used to examine the functions of mutated m6A RNA methylation regulators and 170 genes significantly related to the mutations of m6A regulators using the DAVID database. GO analysis was used to predict the functions of mutated m6A RNA methylation regulators and 170 genes significantly related to the mutations of m6A regulators, whereas KEGG analysis was used to predict the related pathways of mutated m6A RNA methylation regulators and 170 co-expressed genes significantly related to the mutations of m6A RNA methylation regulators. A p-value of <0.05 was considered statistically significant.

TIMER2 (http://timer.cistrome.org/):

The relationship between the expression of the YTHDF gene family and immune infiltration across all TCGA tumour types was explored using the ‘Immune Genes’ module of the TIMER2 Web server. Immune cells such as cancer-associated fibroblasts and regulatory T cells (Tregs) were selected. The TIMER, CIBERSORT, CIBERSORT-ABS, QUANTISEQ, XCELL, MCPCOUNTER, TIDE and EPIC algorithms were used for assessing immune infiltration. p-values and partial correlation (cor) values were evaluated using the purity-adjusted Spearman’s rank correlation test. The data were visualised on a heatmap and scatter plot.

RNA extraction and qRT-PCR:

Tumour and para-cancerous tissues were treated with the TRIzol reagent (Invitrogen, USA) for 10 mins and centrifuged at 12,000 g at 4℃ for 15 min. Thereafter, RNA to be suppressed was collected and mixed with isopropanol for its isolation. After RNA was obtained, its purity and concentration were analysed using a NanoDrop 1000 spectrophotometer (Thermo Fisher, USA). cDNA was synthesised using a high-capacity cDNA reverse transcription kit (Life Tec, America). The primers used for genes are listed in Table S1. qRT-PCR was performed using 2X Universal SYBR Green Fast qPCR mix (Abclonal, China) on a LightCycler 96 system (Roche, America). All experiments were performed in triplicate.

Western blotting (WB):

HCC specimens were collected during surgery. The proteins separated on gels after electrophoresis are transferred to PVDF membranes and subsequently incubated with primary and secondary antibodies. The following primary antibodies were incubated for 12 h at 4°C: GAPDH (1:10000, ABclonal), YTHDF1 (1:1000, ABclonal), YTHDF2 (1:1000, ABclonal) and YTHDF3 (1:1000, ABclonal). The secondary antibody (1:5000, ABclonal) was incubated for 2 h at room temperature. Immunoreactive bands were developed using an enhanced chemiluminescence detection kit (Genview Scientific Inc., USA). All experiments were performed in triplicate.

RNA m6A quantification

A total of 10 pairs of surgical specimens (tumour and para-cancerous tissues) were collected from patients with HCC. Total RNA was extracted using the TRIzol reagent (Invitrogen, CA, USA) as described above, and its purity and concentration were analysed using a NanoDrop 1000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA). The m6A modification level of total RNA was examined using the EpiQuik m6A RNA Methylation Quantification Kit (p-9005; Epigentek Group Inc., Farmingdale, NY, USA) according to the manufacturer’s instructions. Briefly, 200 ng of RNA along with the m6A standard were coated on assay wells, followed by the addition of capture antibody and detection antibody solutions. m6A levels were quantified colourimetrically by reading the absorbance of each well at a wavelength of 450 nm (OD450), and subsequent evaluation was performed based on the standard curve.

Ethics Statement:

This study was approved by the Medical Ethics Committee of Taizhou Central Hospital. All experiments were conducted following the study protocol and the research was carried out following the guidelines of the ethics committee listed in the ethics statement.

Statistical analysis:

All descriptive data were expressed as mean ± standard deviation. Statistical analyses were performed using the SPSS Statistics version 21.0 (IBM, USA) software. Paired t-test was used to analyse differences between two groups. A p-value of <0.05 indicated significant differences.

Results

1. Abnormal expression of m6A RNA methylation regulators is associated with tumour stages in patients with HCC

The expression of m6A RNA methylation regulators and PD-L1 in HCC was analysed using ONCOMINE and UALCAN. A total of 50 healthy tissues and 371 HCC tissues were examined using the UALCAN database. As shown in Figure 1, the expression of all genes, except PD-L1, was significantly upregulated in HCC tissues. In addition, in the ONCOMINE database, analysis of multiple HCC datasets showed that the expression of m6A RNA methylation regulators was significantly higher in HCC tissues than in healthy tissues (Table S2), suggesting that m6A RNA methylation regulators play an important role in the occurrence and development of HCC and may serve as molecular markers for early diagnosis and prognosis of HCC. Furthermore, the relative expression levels of each m6A RNA methylation regulator and PD-L1 were compared in HCC tissues. As shown in Figure 2, m6A RNA methylation regulators were generally significantly upregulated in HCC tissues. 

Furthermore, the relationship between the expression of m6A RNA methylation regulators and clinical stages of HCC was evaluated using the GEPIA database. As shown in Figure 3, the expression of YTHDF1, YTHDF2, YTHDC1, RBM15 and METTL3 was found to be significantly associated with HCC staging, whereas other genes did not have a significant relationship with the clinical stages of HCC. Moreover, the expression of m6A methylation regulators was downregulated in the late stage of HCC (stage 4), suggesting that m6A modification mainly occurs in the early stages of HCC. This finding may also be related to the small sample size of patients with stage 4 disease; therefore, the sample size should be expanded for further verification. In conclusion, the abovementioned results suggest that the expression level of m6A methylation regulators is significantly correlated with individual tumour stage and is higher in patients with early-stage disease. 

2. Prognostic characteristics of m6A RNA methylation regulators in patients with HCC

To investigate the effects of m6A RNA methylation regulators on the prognosis of HCC, the relationship between the expression of m6A regulators and PD-L1 and OS, relapse-free survival (RFS), PFS and DSS was examined using the K-M plotter. As shown in Table S3, patients in each cohort were divided into the low- and high-risk groups based on the cut-off values. The cut-off values for OS associated with m6A RNA methylation regulators are shown in Figure S2. As shown in Figure 4, upregulated mRNA expression of KIAA1429 was associated with poor OS (p = 0.045) but was not associated with RFS, PFS or DSS. Upregulated mRNA expression of METTL3 was associated with poor OS (p = 0.003), RFS (p = 0.0014), PFS (p = 0.0072) and DSS (p = 0.0068), suggesting that patients with HCC with upregulated METTL3 expression had a poor prognosis. In addition, upregulated mRNA expression of RBM15 was associated with PFS (p = 0.027), and upregulated mRNA expression of WTAP was correlated with RFS (p = 0.0092) and PFS (p = 0.0047). Furthermore, downregulated mRNA expression of ZC3H13 was associated with poor OS (p = 0.00033), RFS (p = 0.025), PFS (p = 0.019) and DSS (p = 0.00082), whereas downregulated mRNA expression of METTL14 was associated with poor OS (p = 0.00039), RFS (p = 0.021) and DSS (p = 0.0039) (Figure 5). Low expression of ZC3H13 and METTL14 as m6A ‘writer’ genes has been associated with a poor prognosis in breast and colon cancers[26-28] and may lead to downregulation of m6A RNA modification in tumours, thereby reducing the level of immune cell infiltration and resulting in a poor prognosis, which is consistent with the results of this study. Upregulated mRNA expression of HNRNPC was associated with poor OS (p = 0.028) and DSS (p = 0.012) (Figure 5). Upregulated mRNA expression of YTHDC1 was associated with poor PFS (p = 0.04) and good DSS (p = 0.047). In addition, upregulated mRNA expression of YTHDF1 was associated with poor OS (p = 0.0042), RFS (p = 0.00012), PFS (p = 0.0017) and DSS (p = 0.0072), whereas mRNA expression of YTHDF2 was associated with poor OS (p = 0.017) and RFS (p = 0.015). Other m6A RNA methylation regulators associated with the prognosis of HCC are shown in Figure S1. In conclusion, abnormal expression of m6A RNA methylation regulators plays an important regulatory role in the prognosis of HCC, with clinical significance. For example, upregulated expression of METTL3, WTAP, YTHDF1 and YTHDF2 indicates a poor prognosis of HCC. However, low expression of ZC3H13 and METTL14 as tumour suppressor genes can also lead to a poor prognosis of HCC, indicating that m6A RNA methylation regulators play an extremely important role in the prognosis of HCC.

3. Mutations in m6A RNA methylation regulators and their relationship with survival in patients with HCC

Epigenetic changes play a crucial role in the occurrence and development of tumours. The effects of mutated m6A RNA methylation regulators on the progression and prognosis of HCC remain unknown. In this study, mutations in m6A RNA methylation regulators in patients with HCC were analysed using the online tool cBioPortal for Cancer Genomics (TCGA, Firehose Legacy). Mutated m6A RNA methylation regulators were found in 331 (89%) of the 371 patients with HCC (Figure 6A). Among these regulators, the mutation rates of KIAA1429, YTHDF3, ALKBH5, WTAP, YTHDF1, YTHDF2 and FTO were 43%, 29%, 27%, 23%, 21%, 18% and 17%, respectively. In addition, the mRNA expression of these regulators and PD-L1 (RNA Seq V2 RSEM) in LIHC (TCGA, Firehose Legacy) was analysed using cBioPortal, followed by Pearson correlation analysis. The correlation between RNA methylation regulators and PD-L1 was observed (Figure 6B), and the results showed that KIAA1429 was positively correlated with PD-L1 and YTHDF3; METTL3 was positively correlated with PD-L1, YTHDF1, HNRNPC, METTL14, MTEEL16 and RBM15; RBM15 was positively correlated with YTHDC2, YTHDC1, HNRNPC and METTL14; ZC3H13 was positively correlated with ALKBH5; METTL16 was positively correlated with YTHDC2 and METTL14; METTL14 was positively correlated with FTO and YTHDC2 and HNRNPC was positively correlated with PDCD1 and YTHDF1. In addition, YTHDC1 was positively correlated with PD-L1; YTHDF1 was positively correlated with PDCD1; YTHDF2 was positively correlated with YTHDF3 and PD-L1 was positively correlated with PDCD1. In addition, the relationship between mutated m6A RNA methylation regulators and OS and disease-free survival (DFS) in patients with HCC was investigated. The K-M plot and log-rank test results showed that mutated m6A RNA methylation regulators were associated with shorter OS (Figure 6C, p = 5.339E-3) and DFS (Figure 6C, p = 4.411E-3) in patients with HCC. These results suggest that mutations in m6A RNA methylation regulators and the PD-L1 gene play an important role in the prognosis of HCC.

4. GO and KEGG enrichment analyses of m6A RNA methylation regulators, PD-L1 and 170 co-expressed genes in patients with HCC

After analysing the gene mutations of m6A RNA methylation regulators and their prognostic value in patients with HCC, the ‘co-expression’ module of cBioPortal was used to analyse the top 10 co-expressed genes (among 170 genes) significantly associated with the mutations of each m6A RNA methylation regulator (see Table S4 for specific co-expressed genes). Subsequently, a PPI network was established using the STRING database. As shown in Figure 7A, the HNRNP family (HNRNPA1/3, HNRNPL and HNRNPU), SNRPD1/3, SRSF3/7, PHF5A and PRPF38A were found to be closely associated with mutated m6A RNA methylation regulators. The HNRNP family is a class of RNA-binding proteins with various key cellular functions. Cell dysfunction, including selective splicing, translation and RNA processing, plays an important role in tumorigenesis. Therefore, the HNRNP family has attracted increasing attention owing to its association with cancer progression. In addition, SNRPD1/3 plays a regulatory role in breast cancer through cell cycle regulation, and SRSF3/7 is an RNA-binding protein associated with metastasis and recurrence. PHF5A and PRPF38A play an important role in regulating tumour proliferation and migration. Therefore, in this study, GO and KEGG analyses were performed using DAVID to analyse the potential role of m6A RNA methylation regulators and their 170 co-expressed genes in HCC (Table S5). As shown in Figure 7B–E, mutated m6A RNA methylation regulators were found to significantly regulate biological processes (BPs), such as regulation of signal transduction by a p53 class mediator, cell proliferation, T-cell co-stimulation, mRNA processing and positive regulation of T-cell proliferation (Figure 7B). Furthermore, molecular functions included protein, poly(A) RNA and protein kinase binding (Figure 7C), whereas cellular components included the cytoplasm, nucleus and nucleoplasm (Figure 7D). KEGG analysis indicated that mutated m6A RNA methylation regulators were associated with the activation of the T-cell receptor signalling pathway, RNA transport and spliceosomes (Figure 7E). Moreover, the relationship between mutated m6A RNA methylation regulators and the clinical characteristics of patients with HCC was also analysed (Figure S3, Table S6). The results suggested that mutated m6A RNA methylation regulators were associated with vascular invasion. In addition, the cancer type, neohistological grade and primary tumour site were significantly correlated with the mutated regulators, suggesting that the regulators were associated with a poor prognosis.

5. Immune infiltration analysis of m6A RNA methylation regulators in the tumour microenvironment of patients with HCC 

   Tumour-infiltrating immune cells, as an important part of the tumour microenvironment, are closely related to the occurrence, progression or metastasis of tumours. Previous studies have shown that differences in tumour tissue expression, survival and gene mutations were highly significant among genes in the YTHDF family. In addition, studies have shown that m6A RNA methylation regulators play an important role in the regulation of tumour immunity. Therefore, the underlying mechanisms and role of the YTHDF family in immune infiltration within the tumour microenvironment should be further investigated. In this study, the CIBERSORT, CIBERSORT-ABS, QUANTISEQ, XCELL, MCPCOUNTER, TIDE and EPIC algorithms were used to investigate the potential relationship between the infiltration levels of various immune cells and the expression of the YTHDF family in TCGA-HCC cohort (371 patients). The role of the YTHDF family in the infiltration of tumour-associated fibroblasts and Tregs was determined using the MCPCOUNTER, TIDE and EPIC algorithms, which showed that the expression of the YTHDF family was significantly positively correlated with the infiltration of tumour-associated fibroblasts in HCC (Figure 8A). In addition, analyses performed using the CIBERSORT, CIBERSORT-ABS, QUANTISEQ and XCELL algorithms (Figure 8B) showed that the expression of the YTHDF family was significantly correlated with Treg infiltration in HCC, suggesting its important role in the regulation of tumour immunity and immune escape.

6. qRT-PCR and WB validated the upregulation of m6A RNA methylation regulators in HCC

Analyses performed using UALCAN, TCGA and other databases revealed that the expression of m6A RNA methylation regulators was upregulated in HCC. To further verify these results, surgical specimens were collected from 6 patients clinically and pathologically diagnosed with HCC. qRT-PCR and WB were used to verify changes in the expression of m6A RNA methylation regulators. As shown in Figure 9A, the results of qRT-PCR showed that the expression of m6A RNA methylation regulators was upregulated in HCC, and the results of WB (Figure 9B) showed that the expression of the YTHDF family was upregulated in HCC. These results are consistent with those mentioned above. In addition, to examine the potential role of m6A modification in HCC, the EpiQuik m6A RNA Methylation Quantification Kit was used to detect m6A modification levels in the total RNA of tumour and para-cancerous tissues. The results revealed that m6A modification levels were lower in HCC tissues than in the corresponding adjacent tissues (Figure 9C). These results suggest that the overall m6A modification levels in HCC are mainly mediated by eraser and reader proteins, which warrants further study.

Discussion

m6A methylation is the most common form of mRNA modification and plays a key role in post-transcriptional regulation[29]. Dysregulation of m6A methylation regulatory proteins triggers dysregulation of downstream RNA metabolism and is involved in the progression of various tumours[30]. However, the modification level of specific m6A regulators acting as tumour inhibitors or promoters remains unknown[15]. Therefore, their role varies in different tumour types, suggesting an extremely complex regulation of the modification level of m6A methylation[17, 31]. To date, most studies have focused on intrinsic carcinogenic pathways of tumours; therefore, regulatory factors of m6A should be investigated further to elucidate potential regulatory mechanisms of m6A methylation. At present, the effects of m6A RNA methylation on the occurrence of HCC remain unclear.

In this study, the expression profile and prognostic value of m6A regulators, the role of mutated m6A regulators and immune infiltration in HCC were analysed. The expression of m6A methylation regulators was found to be upregulated in HCC tissues, and upregulated expression of YTHDF1/2, YTHDC1, RBM15 and METTL3 was significantly correlated with tumour stages of patients with HCC. In addition, abnormal expression of m6A methylation regulators reduced OS, RFS, PFS and DSS in patients with HCC, which is consistent with the findings of previous studies. For example, Chen et al. found that upregulation of METTL3 was associated with a poor prognosis in patients with HCC[32]. In another study, YTHDF2 expression was upregulated in PCA tissues and cell lines (DU-145 and PC3), and YTHDF2 knockout significantly increased m6A levels in DU-145 and PC3 cell lines and inhibited cell proliferation and migration[33]. Furthermore, FTO and YTHDF1 were associated with a poor prognosis in patients with HCC and enhanced the translation of PKM2 mRNA by reducing m6A modification and accelerating liver cancer invasion. YTHDF1 is strongly expressed in HCC and is closely related to poor prognosis[34]. Therefore, m6A regulators are more likely to regulate mRNA stability, translation and splicing and can be considered potential biomarkers and prospective therapeutic targets for HCC. Furthermore, PDCD1 expression was upregulated in HCC tissues, whereas PD-L1 expression was not significantly different between HCC and healthy tissues. Survival analysis suggested that low expression of PDCD1 was associated with poorer prognosis; however, PD-L1 did not show a significant association. Although several studies have confirmed that abnormal expression of PD-L1 in tumours promotes tumour immune escape and blocking of anti-tumour responses is enhanced by PD-L1[35], the prognostic effects of PD-L1 in HCC remain inconsistent, with contradictory results reported in some studies. For instance, Sideras et al. [36]reported that high PD-L1 expression indicated significantly better survival in HCC. Adaptive immune resistance was used to explain the inconsistency. PD-L1 upregulation represents the presence of immunosurveillance and may be associated with a better prognosis. These results indicate that the two mechanisms may coexist, and the main mechanism may change at different time points according to the varying immunogenicity of tumours [37].

In recent years, oncogenic mechanisms underlying tumour prognosis have emerged as a prime focus area for research. In this study, upregulated expression of m6A regulatory factors, such as the YTHDF family, METTL3, RBM15, WTAP and KIAA1429, was found to be closely associated with a poor prognosis of HCC. In addition, low expression of ZC3H13 and METTL14 was associated with a poor prognosis. Gong et al. found that low expression of ZC3H13 and METTL14 was associated with a poor prognosis in breast and colon cancers[2628], a finding consistent with that of the present study. In addition, Ma et al. [38] found that downregulation of METTL14 is associated with HCC metastasis and can be used as a prognostic factor for HCC, and METTL14 deletion can enhance the ability of HCC metastasis. Chen et al. [39] found that ALKBH5 was downregulated in HCC, and reduced ALKBH5 expression was an independent prognostic factor for identifying worsening survival of patients with HCC. This finding is different from that obtained using the ONCOMINE and UALCAN databases in this study. This inconsistency may be because HCC is mostly complicated by hepatitis B virus infection in Eastern countries. Moreover, the results of the two studies on FTO remain unclear; FTO may serve as an oncogene or a tumour suppressor gene in HCC. Previous studies have also indicated that the effects of FTO on the proliferation ability of different HCC cells are controversial [40, 41]. In addition, Wang et al. [42] demonstrated that the absence of METTL3 and METTL14 reduced the viability of HeLa cells. In this study, upregulation of METTL3 and downregulation of METTL14 were associated with the poor prognosis of HCC. Although METTL3, METTL14 and WTAP can form RNA methyltransferase complexes, which regulate the fate of tumour development, studies have shown that METTL3 and METTL14 also have independent regulatory functions in regulating transcription. For example, Liu et al.[43] established the function of the METTL3–MeTTL14 complex independent of m6A. The reassignment of METTL3 to the promoter of SASP gene and the reassignment of METTL14 to the enhancer of SASP gene suggest that the METTL3–METTL14 complex plays an important role in regulating transcription independent of its m6A function, thereby playing an independent regulatory role. Liu et al. demonstrated that METTL3 and METTL14 play a contradicting regulatory role in HCC, and their expression and prognostic value are also opposite, which indicates the synergistic role of METTL3 and METTL14 in the catalytic modification of m6A[44]. These studies indicate that there is a complex regulatory network of m6A modification in different tissues, cell lines and spatiotemporal models and also verify the dual and complex regulatory mechanisms of m6A modification in tumours. Therefore, although METTL14 and ALKBH5 play an important role in inhibiting HCC metastasis, they may regulate the progression of HCC through different functional mechanisms dependent on or independent of m6A modification, which requires further investigation and validation.

Epigenetic changes play an important role in early malignant tumours[45]. In this study, m6A regulators were found to have an extremely high mutation rate (89%) in HCC, and the high mutation rates were significantly associated with shorter OS and DFS. Studies have shown that the mutation of m6A regulatory factors is significantly correlated with clinicopathological characteristics, and that of TP53 is significantly associated with poor OS and DFS[46]. Zhu et al. found that the mutation and differential expression of m6A resulted in a significant relationship between increased m6A levels and poor survival, especially in patients with HCC with TP53 mutations[47]. Genetic alterations in the m6A gene may cooperate with TP53 and its regulatory targets in the pathogenesis of HCC[47]. This finding is consistent with that of GO and KEGG enrichment analyses in this study (activation of the P53 signal). Therefore, mutated m6A regulatory factors may serve as new clinical targets for HCC treatment.

From an epigenetic point of view, tissue specificity and uneven distribution of m6A modification provide a new direction for understanding the pathogenesis of different diseases, especially tumours. However, whether m6A modification is related to the tumour microenvironment has been rarely reported. The tumour microenvironment plays an important regulatory role in tumorigenesis, and its heterogeneity can reflect patient prognosis and treatment response[48]. Further investigation of the mechanisms of m6A methylation in tumour-infiltrating lymphocytes and evaluation of immune scores can help to promote clinical diagnosis and targeted tumour therapy. Li et al.[29] reported that METTL3 or METT14 deletion induces T-cell deproliferation and differentiation, thereby reducing IL-7 sensitivity in vivo. Han et al.[21] found that the infiltration levels of CD8 + T and NK cells were increased in mouse models of YTHDF1-deficient tumours, which enhanced the cross-expression of tumour antigens and cross-primers of CD8 + T cells in vivo. These studies suggest that m6A RNA methylation regulators are, to some extent, involved in time regulation and tumour immune cell infiltration. In this study, the most significant differences in the expression pattern, prognostic value and gene mutations were observed among the members of the YTHDF family in HCC. Therefore, the immune infiltration status of the YTHDF family in HCC was further analysed, which showed that the YTHDF family was significantly associated with tumour-associated fibroblasts and Tregs in the HCC tumour microenvironment, which may promote immune escape. Moreover, Tregs play an important role in regulating tumour immunotherapy and immune escape. Therefore, mechanisms underlying m6A regulator-mediated Treg immune escape in HCC should be investigated further.

In conclusion, the expression, prognostic value and immune invasion levels of m6A RNA methylation regulators in HCC were elucidated, and potential mechanisms of m6A RNA methylation regulators were further analysed. m6A RNA methylation regulators were found to be upregulated in HCC and were significantly associated with the clinical stage of patients with HCC. In addition, higher m6A mutations were associated with shorter OS and DFS. Immune infiltration analysis showed that the YTHDF family was significantly associated with tumour immune infiltration and escape. Therefore, m6A RNA methylation regulators play an important role in the tumorigenesis of HCC and may be a risk factor for the survival of patients with HCC, especially the YTHDF family. Targeting the YTHDF family may be a promising uncharted area for future research into HCC treatment.

Declarations

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.

Acknowledgements:

None.

Funding:

None.

Author Contributions

All authors contributed to the study conception and design. Material preparation, data collection and analysis were performed by Kunpeng Wang. The first draft of the manuscript was written by Fanhua Kong and all authors commented on previous versions of the manuscript. All authors read and approved the final manuscript.

Data Availability Statements

  The datasets generated during and/or analysed during the current study are available in the data repository:1. ONCOMINE (http://www.oncomine.org). 2. UALCAN (http://ualcan.path.uab.edu/). 3. GEPIA (http://gepia.cancer-pku.cn/). 4. Kaplan–Meier Plotter (http://kmplot.com/analysis/). 5. TCGA Data and cBioPortal for Cancer Genomics (http://www.cbioportal.org). 6. STRING (http://www.string-db.org). 7. DAVID (https://david.ncifcrf.gov/summary.jsp). 8. TIMER2 (http://timer.cistrome.org/).

References

  1. Bray F, Ferlay J, Soerjomataram I, Siegel RL, Torre LA, Jemal A. Global cancer statistics 2018: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin. 2018;68(6):394-424.
  2. Chen W, Zhang Z, Fang X, et al. Prognostic value of the ALBI grade among patients with single hepatocellular carcinoma without macrovascular invasion. Medicine (Baltimore). 2021;100(24):e26265.
  3. Petrick JL, Florio AA, Znaor A, et al. International trends in hepatocellular carcinoma incidence, 1978-2012. Int J Cancer. 2020;147(2):317-30.
  4. Kong F, Zou H, Liu X, et al. miR-7112-3p targets PERK to regulate the endoplasmic reticulum stress pathway and apoptosis induced by photodynamic therapy in colorectal cancer CX-1 cells. Photodiagnosis Photodyn Ther. 2020;29:101663.
  5. McGlynn KA, Petrick JL, El-Serag HB. Epidemiology of Hepatocellular Carcinoma. Hepatology. 2021;73 Suppl 1(Suppl 1):4-13.
  6. Kong FH, Ye QF, Miao XY, et al. Current status of sorafenib nanoparticle delivery systems in the treatment of hepatocellular carcinoma. Theranostics. 2021;11(11):5464-90.
  7. Jemal A, Ward EM, Johnson CJ, et al. Annual Report to the Nation on the Status of Cancer, 1975-2014, Featuring Survival. J Natl Cancer Inst. 2017;109(9).
  8. Ahn JC, Teng PC, Chen PJ, et al. Detection of Circulating Tumor Cells and Their Implications as a Biomarker for Diagnosis, Prognostication, and Therapeutic Monitoring in Hepatocellular Carcinoma. Hepatology. 2021;73(1):422-36.
  9. Kong FH, Miao XY, Zou H, et al. End-stage liver disease score and future liver remnant volume predict post-hepatectomy liver failure in hepatocellular carcinoma. World J Clin Cases. 2019;7(22):3734-41.
  10. Cao J, Kong FH, Liu X, Wang XB. Immunotherapy with dendritic cells and cytokine-induced killer cells for hepatocellular carcinoma: A meta-analysis. World J Gastroenterol. 2019;25(27):3649-63.
  11. Fang Z, Hu Y, Hu J, Huang Y, Zheng S, Guo C. The crucial roles of N(6)-methyladenosine (m(6)A) modification in the carcinogenesis and progression of colorectal cancer. Cell Biosci. 2021;11(1):72.
  12. Zheng F, Du F, Zhao J, et al. The emerging role of RNA N6-methyladenosine methylation in breast cancer. Biomark Res. 2021;9(1):39.
  13. Kong F, Liu X, Zhou Y, et al. Downregulation of METTL14 increases apoptosis and autophagy induced by cisplatin in pancreatic cancer cells. Int J Biochem Cell Biol. 2020;122:105731.
  14. Zhong X, Yu J, Frazier K, et al. Circadian Clock Regulation of Hepatic Lipid Metabolism by Modulation of m(6)A mRNA Methylation. Cell Rep. 2018;25(7):1816-28.e4.
  15. Huo FC, Zhu ZM, Pei DS. N(6) -methyladenosine (m(6) A) RNA modification in human cancer. Cell Prolif. 2020;53(11):e12921.
  16. Li D, Zhu X, Li Y, Zeng X. Novel insights into the roles of RNA N(6)-methyladenosine modification in regulating gene expression during environmental exposures. Chemosphere. 2020;261:127757.
  17. Sun T, Wu R, Ming L. The role of m6A RNA methylation in cancer. Biomed Pharmacother. 2019;112:108613.
  18. Huang J, Chen Z, Chen X, Chen J, Cheng Z, Wang Z. The role of RNA N (6)-methyladenosine methyltransferase in cancers. Mol Ther Nucleic Acids. 2021;23:887-96.
  19. Yi L, Wu G, Guo L, Zou X, Huang P. Comprehensive Analysis of the PD-L1 and Immune Infiltrates of m(6)A RNA Methylation Regulators in Head and Neck Squamous Cell Carcinoma. Mol Ther Nucleic Acids. 2020;21:299-314.
  20. Keir ME, Butte MJ, Freeman GJ, Sharpe AH. PD-1 and its ligands in tolerance and immunity. Annu Rev Immunol. 2008;26:677-704.
  21. Han D, Liu J, Chen C, et al. Anti-tumour immunity controlled through mRNA m(6)A methylation and YTHDF1 in dendritic cells. Nature. 2019;566(7743):270-4.
  22. Chandrashekar DS, Bashel B, Balasubramanya SAH, et al. UALCAN: A Portal for Facilitating Tumor Subgroup Gene Expression and Survival Analyses. Neoplasia. 2017;19(8):649-58.
  23. Tang Z, Li C, Kang B, Gao G, Li C, Zhang Z. GEPIA: a web server for cancer and normal gene expression profiling and interactive analyses. Nucleic Acids Res. 2017;45(W1):W98-w102.
  24. Nagy Á, Lánczky A, Menyhárt O, Győrffy B. Validation of miRNA prognostic power in hepatocellular carcinoma using expression data of independent datasets. Sci Rep. 2018;8(1):9227.
  25. Szklarczyk D, Gable AL, Lyon D, et al. STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets. Nucleic Acids Res. 2019;47(D1):D607-d13.
  26. Gong PJ, Shao YC, Yang Y, et al. Analysis of N6-Methyladenosine Methyltransferase Reveals METTL14 and ZC3H13 as Tumor Suppressor Genes in Breast Cancer. Front Oncol. 2020;10:578963.
  27. Zhu D, Zhou J, Zhao J, et al. ZC3H13 suppresses colorectal cancer proliferation and invasion via inactivating Ras-ERK signaling. J Cell Physiol. 2019;234(6):8899-907.
  28. Cai C, Long J, Huang Q, et al. M6A "Writer" Gene METTL14: A Favorable Prognostic Biomarker and Correlated With Immune Infiltrates in Rectal Cancer. Front Oncol. 2021;11:615296.
  29. Li HB, Tong J, Zhu S, et al. m(6)A mRNA methylation controls T cell homeostasis by targeting the IL-7/STAT5/SOCS pathways. Nature. 2017;548(7667):338-42.
  30. Roignant JY, Soller M. m(6)A in mRNA: An Ancient Mechanism for Fine-Tuning Gene Expression. Trends Genet. 2017;33(6):380-90.
  31. Lan Q, Liu PY, Haase J, Bell JL, Hüttelmaier S, Liu T. The Critical Role of RNA m(6)A Methylation in Cancer. Cancer Res. 2019;79(7):1285-92.
  32. Chen M, Wei L, Law CT, et al. RNA N6-methyladenosine methyltransferase-like 3 promotes liver cancer progression through YTHDF2-dependent posttranscriptional silencing of SOCS2. Hepatology. 2018;67(6):2254-70.
  33. Li J, Meng S, Xu M, et al. Downregulation of N(6)-methyladenosine binding YTHDF2 protein mediated by miR-493-3p suppresses prostate cancer by elevating N(6)-methyladenosine levels. Oncotarget. 2018;9(3):3752-64.
  34. Zhao X, Chen Y, Mao Q, et al. Overexpression of YTHDF1 is associated with poor prognosis in patients with hepatocellular carcinoma. Cancer Biomark. 2018;21(4):859-68.
  35. Kudo M. Immune Checkpoint Inhibition in Hepatocellular Carcinoma: Basics and Ongoing Clinical Trials. Oncology. 2017;92 Suppl 1:50-62.
  36. Sideras K, Biermann K, Verheij J, et al. PD-L1, Galectin-9 and CD8(+) tumor-infiltrating lymphocytes are associated with survival in hepatocellular carcinoma. Oncoimmunology. 2017;6(2):e1273309.
  37. Schumacher TN, Schreiber RD. Neoantigens in cancer immunotherapy. Science. 2015;348(6230):69-74.
  38. Ma JZ, Yang F, Zhou CC, et al. METTL14 suppresses the metastatic potential of hepatocellular carcinoma by modulating N(6) -methyladenosine-dependent primary MicroRNA processing. Hepatology. 2017;65(2):529-43.
  39. Chen Y, Zhao Y, Chen J, et al. ALKBH5 suppresses malignancy of hepatocellular carcinoma via m(6)A-guided epigenetic inhibition of LYPD1. Mol Cancer. 2020;19(1):123.
  40. Li J, Zhu L, Shi Y, Liu J, Lin L, Chen X. m6A demethylase FTO promotes hepatocellular carcinoma tumorigenesis via mediating PKM2 demethylation. Am J Transl Res. 2019;11(9):6084-92.
  41. Liu X, Liu J, Xiao W, et al. SIRT1 Regulates N(6) -Methyladenosine RNA Modification in Hepatocarcinogenesis by Inducing RANBP2-Dependent FTO SUMOylation. Hepatology. 2020;72(6):2029-50.
  42. Wang Y, Li Y, Toth JI, Petroski MD, Zhang Z, Zhao JC. N6-methyladenosine modification destabilizes developmental regulators in embryonic stem cells. Nat Cell Biol. 2014;16(2):191-8.
  43. Liu P, Li F, Lin J, et al. m(6)A-independent genome-wide METTL3 and METTL14 redistribution drives the senescence-associated secretory phenotype. Nat Cell Biol. 2021;23(4):355-65.
  44. Liu X, Qin J, Gao T, et al. Analysis of METTL3 and METTL14 in hepatocellular carcinoma. Aging (Albany NY). 2020;12(21):21638-59.
  45. Network CGaA. Comprehensive molecular portraits of human breast tumours. Nature. 2012;490(7418):61-70.
  46. Wang P, Wang X, Zheng L, Zhuang C. Gene Signatures and Prognostic Values of m6A Regulators in Hepatocellular Carcinoma. Front Genet. 2020;11:540186.
  47. Zhu GQ, Yu L, Zhou YJ, et al. Genetic Alterations and Transcriptional Expression of m(6)A RNA Methylation Regulators Drive a Malignant Phenotype and Have Clinical Prognostic Impact in Hepatocellular Carcinoma. Front Oncol. 2020;10:900.
  48. Hui L, Chen Y. Tumor microenvironment: Sanctuary of the devil. Cancer Lett. 2015;368(1):7-13.