Identification METTL18 as a Potential Prognosis Biomarker and Associated With Immune Infiltrates in Hepatocellular Carcinoma

Methyltransferase-like 18 (METTL18), a METTL family member, is abundant in hepatocellular carcinoma (HCC). Studies have indicated the METTL family could regulate the progress of diverse malignancies while the role of METTL18 in HCC remains unclear. Data of HCC patients were acquired from the cancer genome atlas (TCGA) and gene expression omnibus (GEO). The expression level of METTL18 in HCC patients was compared with normal liver tissues by Wilcoxon test. Then, the logistic analysis was used to estimate the correlation between METTL18 and clinicopathological factors. Besides, Gene Ontology (GO), Gene Set Enrichment Analysis (GSEA), and single-sample Gene Set Enrichment Analysis (ssGSEA) were used to explore relevant functions and quantify the degree of immune infiltration for METTL18. Univariate and Multivariate Cox analyses and Kaplan–Meier analysis were used to estimate the association between METTL18 and prognosis. Besides, by cox multivariate analysis, a nomogram was conducted to forecast the influence of METTL18 on survival rates. METTL18-high was associated with Histologic grade, T stage, Pathologic stage, BMI, Adjacent hepatic tissue inflammation, AFP, Vascular invasion, and TP53 status (P < 0.05). HCC patients with METTL18-high had a poor Overall-Survival [OS; hazard ratio (HR): 1.87, P < 0.001), Disease-Specific Survival (DSS, HR: 1.76, P = 0.015), and Progression-Free Interval (PFI, HR: 1.51, P = 0.006). Multivariate analysis demonstrated that METTL18 was an independent factor for OS (HR: 2.093, P < 0.001), DSS (HR: 2.404, P = 0.015), and PFI (HR: 1.133, P = 0.006). Based on multivariate analysis, the calibration plots and C-indexes of nomograms showed an efficacious predictive effect for HCC patients. GSEA demonstrated that METTL18-high could activate G2M checkpoint, E2F targets, KRAS signaling pathway, and Mitotic Spindle. There was a positive association between the METTL18 and abundance of innate immunocytes (T helper 2 cells) and a negative relation to the abundance of adaptive immunocytes (Dendritic cells, Cytotoxic cells etc.). Finally, we uncovered knockdown of METTL18 significantly suppressed the proliferation, invasion, and migration of HCC cells in vitro. This research indicates that METTL18 could be a novel biomarker to evaluate HCC patients’ prognosis and an important regulator of immune responses in HCC.


INTRODUCTION
Hepatocellular Carcinoma (HCC) is the fourth most common cause of death of tumor and the seventh most common malignant tumor in 2018 (1). Due to the hepatitis related virus being endemic, the majority of new patients of HCC occur in China annually (2). Most of HCC patients usually lost the opportunity of surgical treatment for the lack of specific symptoms. The molecular mechanisms of liver cancer initiation and progression are still unclear, which make it hard to attain effective treatment (3). In addition, problems remain in the process of treatment and diagnosis of liver cancer, due to the shortage of efficient biomarkers for cancer type and disease stage. In recent years, several serum biomarkers were used to detect the progression and prognosis of hepatocellular carcinoma, such as Alpha-fetoprotein (AFP) and Carcinoembryonic Antigen (CEA) (4,5). However, both the specificity and sensitivity of these biomarkers in evaluating tumor progression, prognosis, and recurrence are still unsatisfactory. Therefore, identification of the reliable and new biomarkers for the diagnosis of HCC is urgently needed to improve the prognosis.
Methyltransferases contains a S-Adenosyl Methionine (SAM) binding domain, known as a family with a structurally conserved methyltransferase similar domains (6)(7)(8)(9). Research studies have shown that methylation could influence chromatin organization and directly regulate transcription of gene, but modulating the mutations in the gene itself (8,10). In addition, studies have demonstrated that methyltransferases could influence the progression of metabolic diseases, genetic diseases, and cancers (11)(12)(13). Studies have shown that members of the METTL family are involved in a variety of biological functions. For example, METTL2B, METTL3, METTL8, and METTL16 have been shown to be RNA methyltransferases (14,15) and play pivotal roles in tumorigenesis (16). In addition, studies have also shown that the METTL family could play a key role in HCC. Chen et al. have demonstrated that METTL3 could promote the progression of liver cancer, and Ma et al. have indicated the metastatic ability of HCC could be suppressed by METTL14 (17,18). However, the influence of other members of the METTL family in progression of hepatocellular carcinoma has not been thoroughly studied.
Methyltransferase-like 18, also known as METTL18, has been regarded as a protein-coding gene (19). A recent research indicated that METTL18 is an effective candidate protein for patients suffering Fibromyalgia syndrome (5). UniProt, the most authoritative protein database demonstrates that METTL18 performs as a methyltransferase, participating in modification of essential protein and process of binding heat shock protein (20). Nevertheless, little is known about METTL18 in cancer.
To better analyze the role of METTL18 in progression of HCC, we applied RNA-seq data from the TCGA and GEO datasets, with statistical and bioinformatics ways, such as differentially expressed genes (DEG) analysis, Kaplan-Meier (KM) survival analysis, and Cox & Logistic regression analysis, nomogram, Gene Ontology (GO) analysis, Gene Set Enrichment Analysis (GSEA), singlesample Gene Set Enrichment Analysis (ssGSEA). Moreover, we further knock down the expression of METTL18 in vitro to detect the impact on the ability of proliferation, invasion, and migration of hepatocellular carcinoma cells.

Data and Preprocessing
Expression data of RNA with clinical information from HCC patients (included 371 tumor and 50 normal tissues) were acquired from TCGA (https://portal.gdc.cancer.gov/); in addition, 225 tumor and 220 normal samples were acquired from GSE14520. The criteria of exclusion were OS less than 30 days and normal HCC tissues. Then, HTSeq-FPKM information of level 3 has been transformed into Transcripts Per Million (TPM); then TPM information of 371 HCC samples was applied for the next analyses. Unknown and unavailable clinical factors in 371 samples were deemed as missing values, and the information was demonstrated in Supplementary Table 1.

Expression of METTL18 in HCC Samples in the TCGA and GEO Dataset
Applying disease state (normal or tumor) as variable, scatter plots and boxplots were performed to estimate different expression levels of METTL18. Using receiver operating characteristic (ROC) curves, we evaluated the diagnostic performance of METTL18 in HCC patients. By statistical ranking, expression level of METTL18 below or above the median was determined as METTL18-high or METTL18-low. statistically significance. Using GSEA, we have investigated the differences in pathways of biofunction between the METTL18-low and -high patients to analysis METTL18-associated pathways and phenotypes. With 1,000 times, the permutation test was applied to analyze significantly signal pathways. FDR <0.25 and adjusted P <0.01 were recognized as significantly associated genes. Graphical plots and statistically analysis were performed by R package ClusterProfiler (4.0) (21). Using STRING database, the protein-protein interaction (PPI) network was built by the DEGs (22), and the PPI pairs were selected with an interaction score >0.7. To better explore the tumor infiltration levels of immune cells, ssGSEA, interrogating expression data of genes in published gene lists (23), was used to quantify relative infiltration levels of 24 types of immune cell. Comprising 509 genes, the signatures applied containing a multiple set of innate as well as adaptive immune relative cell types (23). To evaluate the association between the infiltration levels of immune cells and METTL18 and the correlation of the different groups of METTL18 expression with infiltration of immune cells, Spearman correlation and Wilcoxon test were used.

Clinical Analysis on Prognostic State, Model Construction and Estimation
Using R package (V3.6.2), we analyzed connection between clinicopathological characteristics and METTL18 by logistic regression and Wilcoxon signed-rank sum test. By Kaplan-Meier method and Cox regression, we analyzed the clinical pathologic factors related to 10-year overall survival (OS), disease-specific survival (DSS), and progression-free interval (PFI) in TCGA samples. Multivariate Cox analysis was applied to analyze the influence of METTL18 expression on prognosis along with other clinicopathologic factors (TNM stage, Histologic grade, Vascular invasion, Residual tumor, Albumin, AFP, TP53 status, Race, Child-Pugh grade, Adjacent hepatic tissue inflammation, Age, Gender, and Prothrombin time). Median of METTL18 expression level was chosen as the cut-off value. In all tests, P < 0.05 was defined statistically significant. Using KM method with a log-rank test, the difference of OS, DSS, as well as PFI between METTL18-low and -high group was analyzed. Acquired from multivariate Cox analysis, we used the independent prognostic factors to construct nomograms, evaluating the prognosis for 1, 3-, and 5-year, respectively. Applying the RMS package (https://cran.r-project.org/web/packages/rms/index. html), we established nomograms that included calibration plots and significant clinical factors. By drafting the nomogram evaluating probability against actually occurrence, the calibration curves were graphical evaluated, and the 45 degrees line indicated the best predictive values. Analyzed by bootstrap way with 1,000 resamples, a concordance index (C-index) was analyzed and applied to estimate the discrimination of the model. The separate prognostic factors and predictive accuracies of the nomogram were evaluated by the Cindex. The statistical test was two tailed with a statistically significant level set at 0.05 in our research.

Cell Culture and Transfection
Hepatocellular cancer cell lines (HepG2, M97H, LM3, Bel7402, SK-HEP1, and Huh7) were acquired from the American Type Culture Collection (Virginia, USA). They were cultured in Dulbecco's Modified Eagle Medium (Hyclone) with 10% fetal bovine serum (Gibco) under incubation at 37°C with 5% carbon dioxide. Transfections were performed applying OPTI-MEM (Invitrogen) and Lipofectamine 3000 by the manufacturer's instructions. The siMETTL18 (5′-GACTTTCCTTAG ACTGTTA-3′) and siNC were bought from RiboBio (Guangzhou, China) and introduced into cells at a concentration of 50 nM. The transfected cells were harvested at 48 h after transfection.

Cell Proliferation, Invasion, and Migration Assays
The Cell Counting Kit-8 (CCK-8) and colony formation assays were applied to explore the ability of proliferation of cancer cells in different groups. In CCK-8 experiment, a total of 2,500 cancer cells were added into each well of 96-well plate. 10 µl of CCK-8 solution (Beijing Solarbio Science & Technology Co., Ltd, Bei, China) was added into 96-well, then the absorbance of each well was analyzed at 450 nm after an incubation at 37°C for 2 h. For colony formation experiment, 1,000 cells of different groups were added into each well of a six-well plate. The culture medium was changed every 72 h. Crystal violet and 4% paraformaldehyde were applied to stain and fix the cells when the appearance of colonies could be recognized. The wound healing and transwell assays were applied to explore the ability of cellular migration and invasion.

Statistical Analysis
Statistical analyses were performed by RStudio software and the R software (version 3.8.0). One-way analysis of variance (ANOVA) and two tailed Student's t-test were applied to analyzed the data. Statistical significance of the difference was set at P-value <0.05.

Identification of DEGs in HCC
We compared 186 HCC METTL18-high samples with 185 METTL18-low samples. Between the two groups, a total of 431 DEGs, including 153 downregulated genes and 278 upregulated genes, were found to be statistically significant (adjusted p-value < 0.05 and absolute Log2-fold change > 1.5) ( Figure 1F and Supplementary Table 2). In addition, applying DESeq2 package, we analyzed DEGs in HTSeq-Counts. Ranked by relative expression, the top 10 differential expressed genes between two groups were demonstrated ( Figure 1G).

Enrichment of Biofunction and Analysis of METTL18 Associated Genes in HCC
To better analyze the enrichment of biological function of METTL18 associated genes, we applied Metascape to explore GO enrichment, which demonstrated that METTL18 associated genes were involved in a number of Biological Processes (BPs), Cellular Compositions (CCs), and Molecular Functions (MFs). For instance, different expression of METTL18 could modulate the cellular transition metal ion homeostasis, stress response to metal ion, detoxification of inorganic compound and cellular zinc ion homeostasis. Moreover, cellular response to cadmium ion, kidney development, response to metal ion, appendage morphogenesis, transmembrane receptor protein threonine/serine kinase pathway, and BMP pathway also have the relationship with METTL18 related genes ( Figure 2A and Supplementary Table 3).

Analysis of Protein-Protein Interaction Network
Using PPI network, we explore the association between the 431 DEGs in HCC group, by the STRING dataset and high confidence (0.70) of interaction score is set.

Higher Expression of METTL18 Was Significantly Associated With Poor Prognosis of Patients With HCC
The OS rates of HCC patients were significantly increased among samples with lower expression of METTL18 than patients with higher expression of METTL18 (P = 0.006; Figure 5A). Furthermore, DSS and PFI in lower METTL18 cohort were significantly longer than higher METTL18 cohort (P = 0.015; P = 0.006; Figures 5B, C). we also have validated the prognostic value of METTL18 in GEO datasets and the results indicated that higher expression of METTL18 was significantly related to poor prognosis (P = 0.021; Figure 5D).  Table  3). In addition, increased expression of METTL18 also connected   Tables 8, 9).

Construction and Evaluation of a Nomogram Related to METTL18
To better predict the survival rates of HCC patients, we performed a nomogram based on METTL18 and other independent clinicopathologic factors ( Figure 6A). The point scale of nomogram was applied to assign points to each variable by the results of multivariate Cox regression. With the adjusted range from 1 to 100, the points of each variable were added up and total scores were calculated. By delineating a direct line down from the total score line to the outcome line, the probable prognosis of each HCC patients at 1-, 3-, and 5-years were defined. For example, an HCC patient with METTL18-high risk (89 points), tumor free (0 points), T3 (100 points), and M0 (0 points) could attain a total point of 189. The survival rates of 1-, 3-, 5-year were about 78.5, 48.5, and 32% ( Figure 6A). Furthermore, the efficiency of the nomogram has been evaluated, and the calibration curve with Hosmer test of the nomogram in the TCGA-LIHC cohort was 0.689, which  indicated that the ability of prediction efficiency of the nomogram is moderately accurate ( Figure 6B).

Knockdown of METTL18 Suppress Malignant Phenotype of Hepatocellular Carcinoma In Vitro
To explore the role of METTL18 in LIHC, the expression of METTL18 in the six hepatocellular carcinoma cell lines (HepG2, M97H, LM3, Bel7402, SK-HEP1, and Huh7) was analyzed respectively. The expression of METTL18 in the LM3 and Huh7 cell lines was higher than other cell lines ( Figure 7A). Then, the LM3 and Huh7 cell lines were selected for functional analysis. Using western blotting, we have detected the efficiency of METTL18-siRNA ( Figure 7B). The results of colony formation and CCK8 experiments demonstrated that lower expression of METTL18 significantly inhibited the ability of proliferation and colony formation of LM3 and Huh7 cells (Figures 7C-E). The results of wound healing experiment indicated that knockdown of METTL18 significantly inhibited the migration of LM3 and Huh7 cells (Figures 7F, G). The results of transwell experiment indicated that the LM3 and Huh7 cells exhibited significantly decreased invasion ability upon METTL18 knockdown ( Figures  7H, I). These results proved that METTL18 could be a promoter of hepatocellular carcinoma. Further research studies are needed to reveal the underlying mechanisms.

DISCUSSION
More than 27 members are involved in the METTL family, but only a few of them have been studied (15,(24)(25)(26). Since the METTL family shares a limited conserved domain, it is found that their functions are diverse (24). Recently research studies have demonstrated that members of METTL play a key role in the progression of a number of tumors via multiple mechanisms. For example, METTL3 could induce m6A modification in coding region of mRNA transcription related to the cell cycle, which is essential for the differentiation of leukemia cells and ultimately promotes the translation of a variety of oncogenic mRNAs (27). In addition, a number of members of METTL family have been demonstrated participating in the process of initiation and development of cancers (11)(12)(13). For instance, METTL16, METTL2B, and METTL8 (14,15) and play significant roles in tumorigenesis (16). However, the impact of METTL18 in tumor development has not been explored. Based on our study, the expression levels and prognostic values of METTL18 were evaluated, we found that expression of METTL18 is abnormal in a number of tumors and significantly high in liver cancer in multiple databases. In addition, we found that METTL18 has a relatively high ROC score with an AUC of 0.948 for HCC in the TCGA database. In general, METTL18 is differentially expressed in tumor and normal samples. However, further prospective research studies are warranted to form the diagnostic accuracy of METTL18 in HCC.
To better explore the METTL18 function, we used GO and GSEA to perform functional analysis. We found that higher METTL18 phenotype was associated with cellular zinc ion homeostasis, detoxification of copper ion, stress response to copper ion, zinc ion homeostasis by GO analysis. Recently, multiple studies have demonstrated that alterations of metal molecular, such as Zinc and Copper, could modulate various molecular targets which playing an important role in progression and development in various cancers, such as cyclic AMP (cAMP)-dependent protein kinase (PKA), protein kinase C (PKC), and transcription factors (including nuclear factor (NF)-kB) (28)(29)(30)(31)(32). In addition, using Cytoscape, the PPI network of METTL18-related genes was performed, and the results indicated that diverse biological processes and signaling pathways have association with these genes. In future research studies, we will explore the relationship between METTL18-related genes and prognosis of HCC. We also indicated that METTL18 was significantly associated with G2M checkpoint, E2F targets, KRAS signaling, and mitotic spindle by GSEA. G2M checkpoint and KRAS signaling pathway have been proved to play an important role in the progression and development of liver cancer (33)(34)(35). A previous study has demonstrated biological processes related to cell proliferation like G2M checkpoint, and cell cycle were enriched in abnormal expression of RNA-methyltransferase NSUN6 (36). Furthermore, research studies have demonstrated that E2Fs are the downstream of cell cycle signaling pathway and have a significant impact in modulation of cell proliferation, differentiation, and apoptosis (37)(38)(39). Also, mitotic spindle has been found enhancing chromosomal instability and liver cancer progression (35). These studies and our results indicated that METTL18 might contribute to HCC initiation and development by modulating cell cycle and KRAS pathway. But the association of METTL18    expression with G2M checkpoint, KRAS signaling, mitotic spindle and E2F targets was the first to be reported, and the regulatory molecular mechanism needs to be further explored.
Recently, research studies have shown that tumor infiltrating immune cells (TIICs) could modulate the process of development as well as progression of tumor (40). In addition, TIICs are potentially clonally expanded and preferentially enriched in HCC (41), and poor prognosis correlates with accumulation of TIICs in HCC (42). Then, our results demonstrated that METTL18 expression in HCC was negatively associated with multiple types of immune cell infiltration, for example, dendritic cells (DCs, iDCs, and pDCs) and cytotoxic cells. DCs, also known as a group of antigen-presenting cells, have a significant impact on the initiation and regulation of cancer immune responses (43). Recently, one study has shown DCs can generate resistance to HCC (44). Immature DCs have the phagocytic ability. However, mature DCs have a significant modulation function as well as produce a large number of cytokines (45). These results demonstrated that the DCs had a potent negatively correlation with METTL18. In addition, cytotoxic cells, also known as CD8 + T lymphocytes with cytotoxic granules, are the important anti-tumor effector cells (46). One research has reported that hepatocellular carcinoma cells could mediate the progression, development and tumor resistance to PD1 by inhibiting the cytotoxic T cell response (47). Therefore, our results reveal the potent modulating role of METTL18 in immune response with liver cancer. Furthermore, significantly connection can be found between expression of METTL18 and the modulation of T helper cells, Th2 cells, and Tfh in HCC. These associations could be indicative of a potent mechanism where METTL18 modulates T cell functions in HCC. In summary, these results demonstrate that the METTL18 plays a significant role in modulation of immune infiltrating cells in HCC.
Although research studies have demonstrated that there was diverse association between the expression of member of METTL and clinical pathological factors in HCC (17,18,48). In our study, results have indicated that elevated expression level of METTL18 was related to poor prognosis and advanced clinical pathologic characteristics in HCC. In a stratified analysis, we found that METTL18 expression remained a powerful factor to forecast the prognosis within these subgroups, such as T2 to 4, stage 2 to 4, N0, M0 and histologic grade G3 and G4 etc, indicating that METTL18 was independent of these significant clinicopathological parameters.
After regulating the conventional clinicopathological factors, our study indicates that METTL18 could perform as an independent prognostic factor of poor OS and DSS for HCC. Then, our METTL18 related nomogram was performed with the expression level of METTL18 and other clinical factors (cancer status, T stage, M stage, etc). The C-index of METTL18-related Cox model for overall survival prediction was 0.689 (95%CI: 0.645−0.732). The calibration plots demonstrated optimal agreement between the prediction by METTL18-related nomogram and actual observation for 1-, 3-, and 5-year OS probability. Our nomogram was performed based on the complementary perspective for individual patients and provided a personalized score for respective tumors. In addition, our model could be a new method to estimate the prognosis of clinicians in the future. Furthermore, the impact of METTL18 on the ability of proliferation, invasion, and migration of hepatocellular carcinoma cells was explored in vitro. We found that malignant phenotype of HCC cells was suppressed when METTL18 is knocked down, indicating the expression of METTL18 is a potent target for hepatocellular cancer therapy. In future research studies, the underlying mechanisms of METTL18 in HCC should be elucidated.
Although this research improved our understanding about the correlation between METTL18 and HCC, there still existed some limitations. Firstly, to better evaluate the significant role of METTL18 in progression of hepatocellular carcinoma, all types of clinical characteristics should be involved, for example, the way of treatment for each patient. However, because the experiments were done in different center, this kind of missing information were inevitable in public databases. Second, the relationship between the expression level of METTL18 and prognosis should be verified using more clinical samples, and using the single cohort from public datasets to predict the prognosis is far from perfect. In addition, because the limitations of this research design, the significant pathways correlated with METTL18 may have been missed, and the related signaling should be explored further. To better examine the mechanism of METTL18 in HCC, we would conduct more experimental studies on METTL18 in the sooner future.

CONCLUSION
In our study, we demonstrated METTL18 as an important molecular biomarker with prognostic value and may have significant impact on the modulation of cell cycle, KRAS signaling and immune infiltration in HCC. This research supported promising visions for subsequent study to clarify the molecular pathogenesis and clinicopathological significance of HCC. Randomized clinical trials and further studies are needed to analyze the underlying mechanism and clinical applications for HCC patients.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding author.

ETHICS STATEMENT
The ethics consent and approval to participate in the current study was approved by the ethics committee of Peking Union Medical College Hospital.

AUTHOR CONTRIBUTIONS
T-HL and CQ contributed equally to this article. T-HL designed the research and participated in manuscript writing. B-BZ, CQ, H-TC, and X-TZ participated in generating and analyzing the information. X-YY, Y-YW, and Z-RL assisted with analysis. All authors contributed to the article and approved the submitted version.