DNA hypermethylation is associated with invasive phenotype of malignant melanoma

Tumor cell invasion is one of the key processes during cancer progression, leading to life‐threatening metastatic lesions in melanoma. As methylation of cancer‐related genes plays a fundamental role during tumorigenesis and may lead to cellular plasticity which promotes invasion, our aim was to identify novel epigenetic markers on selected invasive melanoma cells. Using Illumina BeadChip assays and Affymetrix Human Gene 1.0 microarrays, we explored the DNA methylation landscape of selected invasive melanoma cells and examined the impact of DNA methylation on gene expression patterns. Our data revealed predominantly hypermethylated genes in the invasive cells affecting the neural crest differentiation pathway and regulation of the actin cytoskeleton. Integrative analysis of the methylation and gene expression profiles resulted in a cohort of hypermethylated genes (IL12RB2, LYPD6B, CHL1, SLC9A3, BAALC, FAM213A, SORCS1, GPR158, FBN1 and ADORA2B) with decreased expression. On the other hand, hypermethylation in the gene body of the EGFR and RBP4 genes was positively correlated with overexpression of the genes. We identified several methylation changes that can have role during melanoma progression, including hypermethylation of the promoter regions of the ARHGAP22 and NAV2 genes that are commonly altered in locally invasive primary melanomas as well as during metastasis. Interestingly, the down‐regulation of the methylcytosine dioxygenase TET2 gene, which regulates DNA methylation, was associated with hypermethylated promoter region of the gene. This can probably lead to the observed global hypermethylation pattern of invasive cells and might be one of the key changes during the development of malignant melanoma cells.


| INTRODUC TI ON
Melanoma is a neural crest-derived tumor that develops from melanocytes originating from a highly migratory embryonic cell population. [1][2][3] The mechanism of migration and invasion, key processes of cancer cell progression which leads to life-threatening metastatic tumors, is poorly understood. [4][5][6] According to recent studies, the cellular plasticity that promotes invasion strategies in malignant melanoma is mainly due to environmental stimuli and is accompanied by transcriptomic reprogramming. [7] Based on a recent study by Verfaillie et al, [7] the proliferative and invasive transcriptomic signatures are highly correlated with permissive and repressive chromatin states underlining the importance of epigenetic regulation in the acquisition of the invasive cellular state.
Indeed, the expression of the MITF and SOX10 transcription factors, which are master regulators of the proliferative gene network, has been confirmed. In the contrary, invasive cells exhibit high expression levels of TEAD and AP1 genes. [8,9] Due to the lack of direct genetic components in transcriptional reprogramming, studying the epigenetic factors that may promote cellular plasticity leading to increased invasion and metastasis is reasonable. [10] Based on the Cancer Genome Atlas Network (TCGA), the well-established mutational classifications of melanomas are not in agreement with gene expression patterns, which could explain not only the low response rate of therapies targeting the aforementioned mutations but also the concerns raised against the durability of such interventions. [11] Nevertheless, the strong association between the mutations of chromatin remodelling genes (ARID2 and IDH1) and the high degree of DNA methylation at several promoter regions described in melanoma (CpG island methylator phenotype; CIMP) suggests that epigenetic factors might play a pivotal role in cellular plasticity leading to increased invasion and metastasis. [11,12] Epigenome-wide (EWAS) DNA methylation studies implemented during the last few years have greatly improved our understanding on the importance of CIMP in the silencing of tumor suppressor and developmental genes. Several genes of the melanocyte lineage differentiation pathway were found to be methylated such as KIT, PAX3, SOX10, different members of the HOX family genes and MITF. [13][14][15][16] Remarkably, comparing matched primary and metastatic melanoma cell lines, Chatterjee et al [17] found EBF3 promoter hypermethylation as a possible epigenetic driver of melanoma metastasis.
Importantly, EWAS on melanomas have more often focused on the metastatic tumors, and therefore, the DNA methylation changes accompanying the early molecular invasion events remain to be elucidated. Only a single study used cell lines derived from primary melanomas, and the authors applied melanocytic markers to distinguish between invasive and less invasive cell lines and found that SOX9 demethylation is associated with melanoma cell invasion and metastasis and decreases patient survival. [18] There is a great need of identifying early metastasis-promoting epigenetic events, and our main goal was to study the DNA methylation landscape of early invasion using a direct, in vitro selection for the invasive melanoma subpopulation derived from primary malignant melanomas.

| In vitro invasion assay and selection
The invasive potential of melanoma cell lines was analysed using BD Biocoat Matrigel invasion chambers (pore size: 8 μm, 24-well; BD Biosciences). The upper chamber of the insert was filled with 500 μL of cell suspension in serum-free media (5 × 10 4 cells/well). Medium containing 10% FBS was added to the lower chamber as a chemoattractant. After the cells were incubated for 24 hours at 37°C, cells in the lower layer were fixed with methanol and stained with haematoxylin-eosin. The invading cells were counted in seven different visual fields under a light microscope at ×200 magnification, and the data are presented as the means ± SD of three independent experiments. To select the invasive subpopulations, the invading cells in the lower layer chamber were treated with 0.5% trypsin/0.2% EDTA (Sigma-Aldrich Inc) for recovery from the membrane and cultured using standard protocols. The selected subpopulations were desig-

| Cell proliferation assay
Cell proliferation rate was determined using the WST-1 assay (Sigma-Aldrich Inc) according to the manufacturer's guidelines. Briefly, cells were seeded in 96-well in triplicate and cultured for 24, 48, 72 and 96 hours. Then, 10 µL of WST-1 was added directly to the culture medium in each well, and the cells were incubated for another 3 hours.
The reference absorbance was set at 700 nm.

| Genome-wide DNA methylation analysis
DNA was extracted using a standard procedure of the G-spin Genomic DNA extraction kit (Intron Biotechnology Inc). DNA quantification was done NanoDrop. For methylation studies, bisulphite modification was performed on 600 ng of DNA using EZ DNA Methylation kit (Zymo Research). The quality of modification was confirmed by PCR (HotStarTaq Master Mix kit; Qiagen GmbH) using modified and unmodified primers for the GAPDH gene. The DNA methylome profiling was performed using an Illumina Infinium Human Methylation 450K (HM450K) BeadChip assay (Illumina), which includes more than 480 000 methylation sites. [19] The array experiments were performed by the Epigenetics Group and the Core Facility of the Genetic Cancer Susceptibility Group, International Agency for Research on Cancer.
In brief, raw methylation data were imported and processed using the "Lumi v2.36.0," "wateRmelon v1.28" and "minfi v1.30" packages. [20][21][22] Probes were filtered for low quality with the "pflter" function, and additionally, known cross-reactive probes were also excluded from further analysis. [21,23] The remaining data set was background-subtracted and normalized using intra-array beta-mixture quantile normalization. [24] Methylation beta values were logarithmically transformed to M values before parametric statistical analyses, as recommended. [25] To define differentially methylated positions (DMPs) and differentially methylated regions (DMRs), first, we modelled the main variables (invasive capacity) as a categorical variable in a linear regression using the "limma v3.40.2" package an empirical Bayesian approach. [26] To infer the detected differentially methylated sites into DMRs, we used the "DMRcate v1.20" package with the recommended proximity-based criteria: if a region harboured at least 3 probes spanning in 1 kb. [27] For the annotations, to obtain information of the nearest gene and transcript of each the detected DMR, we used the FDb.InfiniumMethylation.
hg19 v.2.2.0 package, using hg19 as a reference genome. [28] For the visualizations, we used either the DMRcate or the coMET packages with the functionality of the Gviz package. [30]

| Correlation between gene expression and DNA methylation
RNA was isolated using RNeasy Plus Mini Kit (Qiagen GmbH) and then assessed using NanoDrop and Bioanalyzer (Agilent Technologies). To assess gene expression at genome-wide levels,

| TCGA-SKCM data analysis
We downloaded Illumina Methylation 450K data available for SKCM from the TCGA-GDC data portal (https ://portal.gdc.cancer.gov/) by using the GDCquery and GDCprepare functions of the TCGAbiolinks R package. [31] The latter generated a summarized experiment object that we further analysed by using the TCGAanalyze_DMR function of the TCGAbiolinks package with a mean delta-beta cut-off 10% and a Benjamini-Hochberg adjusted P-value of .05. The rest of the settings were the default options recommended by the developers of the package. We compared the tumors classified as "metastatic" vs "primary" according to the definition column of the clinical data available at the data portal.
Afterwards, we added a variable to the colData data frame of the summarized experiment by using the addAnnotation function of the IntEREst R package, [32] consisting of a merge of any Clark level below stage V into a single category to compare primary tumors by invasiveness (V vs not-V), and finally rerun the TCGAanalyze_DMR function as described above.

| Real-time quantitative PCR analysis
The relative expression level of 20 genes that are related to methyla- LightCycler 480 Real-Time PCR System (Roche Diagnostics GmbH) as previously described. [33] The primer sequences are listed in Table S1.

| Phenotypic characterization of selected invasive cells
To identify the invasion-related gene expression changes and genome-wide DNA methylation patterns in melanoma cells, invasive cell subpopulations were selected from the original cell lines

| Methylation profile of selected invasive cells
To define the methylation patterns of the cell lines, we used robust methylation profiling platform, which allowed to compare epigenome-wide data of selected invasive melanoma cell subpopulations  Figure 2E). Furthermore, we observed significant difference in the GC content between hypermethylated probes and all bead array (HM450) probes ( Figure 2D).
We applied more stringent criteria to determine significant DMRs with increased Δβ mean > 10% between the invasive and the  Table S4).

| Integration of methylation and gene expression profiles
To determine the functional relevance of the DNA methylation changes, we performed integrative analysis of the DNA methylation and gene expression alterations.
We identified a total of 886 significantly correlated CpG sites corresponding to 392 individual genes, of which 220 showed negative, whereas 172 genes exhibited positive correlation between DNA methylation and gene expression (Table S5) Interestingly, CpG island shore hypermethylation was associated with decreased expression in case of four DMRs corresponding to RHOB, ID4, ST8SIA1 and GRIA2 genes ( Figure S1).
Additionally, two hypomethylated genes revealed significant correlation with either up-regulation (NNMT; upper left segment; Figure 3D) or down-regulation (NBPF8; lower left segment; Figure 3D) of gene expression.

| Invasion-related methylation changes in melanoma tumor samples
To validate the importance of methylation changes observed in the selected invasive cells, first, we determined and compared our results to the DNA methylation changes present in the TCGA metastatic melanomas (n = 349) vs tissues of primary sites (n = 88). Altogether, 879 genes (corresponding to 1984 differentially methylated regions) exhibited significant differences between the metastatic and primary melanomas of the TCGA cohort (Table S6) However, this comparison has the limitation, that is later metastatic events are not necessarily characteristic for those that arise at the early stages of invasion. [35,36] For this reason, we aimed to concentrate on the 88 primary melanoma tissues and used Clark staging as the most relevant clinical parameter to differentiate between locally invasive (Clark stage V; n = 20) and early stage (all Clark stages below V; n = 41) referred as less invasive. We identified 448 differentially methylated genes (corresponding to 1269 probes) seem to have a role during early invasion represented by the Clark staging system (Table   S6). Of note, 18 out of the 385 genes in our data set show overlap with the TCGA data, of which several genes (MECOM, CHD5, TRIM55, FZD6, TPBG and TRPC4) were observed in association with invasion. [37][38][39][40][41][42] Comparing our data with the TCGA, the most interesting finding is the hypermethylation of ARHGAP22 and NAV2 genes that were commonly detected in locally invasive primary melanomas as well as during metastasis.
F I G U R E 3 Integration of methylation and gene expression profiles related to melanoma invasiveness. A, Positively and negatively correlated probes relative to the genes (promoter, UTRs or intron/exon) and (B) distance to transcription start sites (TSSs). C, Positively and negatively correlated DMP positions relative to DNAse I hypersensitivity sites. D, A starburst plot corresponding to correlation analysis between DNA methylation and gene expression changes. The filtered mean log expression and methylation data are shown in a correlation plot (1-fold expression differences between the invasive and the original cell lines were correlated to DMPs of Δβmean > 10%). Probes with increased gene expression are shown in red, and probes with decreased expression are highlighted in blue  Figure 4C). Unexpectedly, TET2 down-regulation was associated with hypermethylation at the TET2 gene promoter region in the invasive subpopulation ( Figure 4D).

Invasion is a crucial step in metastasis formation in primary tumors
including malignant melanoma. Recent advances in epigenome-wide DNA methylation methods have allowed for the identification of potential biomarkers that could be exploited in clinical settings. [15,[43][44][45][46][47] However, in the case of early invasion steps in primary melanomas, insufficient data are available regarding the epigenetic mechanisms and especially the functionally relevant DNA methylation changes affecting gene expression patterns.
In the present work, we selected invasive cells in vitro from the original cell lines and analysed their invasion-associated DNA methylation changes, which followed by functional analysis of the observed changes at mRNA expression level. A number of studies have indicated that several tumor suppressor genes are silenced by DNA methylation in malignant melanoma compared with normal melanocytes or nevi, for example MITF gene. [15,[48][49][50] MITF (microphthalmiaassociated transcription factor) has been extensively studied in the context of master regulator of melanin production, suppression of invasion and regulation of the proliferative phenotype in melanoma cells. [7,15,51,52] Its methylation change was also observed in melanoma brain metastases, suggesting its role not only in invasion property, but also in metastasis formation. [14,16] Selected invasive melanoma cells also showed hypermethylation of MITF that may directly affect MITF expression, giving a functional role of the detected epigenetic change.
Several studies have indicated that different biological behaviours of melanoma tumors are associated with distinct methylation subgroups. [14,18,53,54] The methylation changes in the TFI2, HCK, MGMT and TP73 genes have been described in association with advanced clinical stage, shorter overall survival and the presence of metastasis, and it seems that, according to our results, these genes have a potential role in the earlier invasion steps of primary melanoma cells. [54][55][56][57]  In agreement with the widely accepted assumption that increased DNA methylation of certain promoters causes deregulation of the corresponding genes, we observed a negative correlation between the methylation and gene expression for several promoters such as FBN1, ADORA2B and CHL1. Hypermethylation was associated with the deregulation of fibrillin-1 (FBN1) that is a major component of microfibrils, and it can mediate cell adhesion in melanoma cells. [58] ADORA2B has been identified as specific receptor for 5'-methylthioadenosine (MTA) that can affect cell invasiveness in melanoma cells. [59] Neural cell adhesion molecule L1 (CHL1) is frequently downregulated in different types of tumors, and it is verified to inhibit invasive growth and able to suppress further metastatic spread. [60] It seems that down-regulation of CHL1 in association with methylation change was also observed in melanoma cells by Chatterjee et al [61] indicating that differentially methylated CHL1 is a marked alteration in melanoma cells as well.
On the other hand, recent studies have shown that the methylation of the gene body is positively correlated with transcription. [45,62,63] Similar to these observations, hypermethylation in the gene body potentially plays a role in the up-regulation of EGFR and RBP4 genes in the selected invasive cell populations. The role of EGFR in a range of neoplasms including melanoma is well known in association with tumor progression and metastasis. [64][65][66] Based on our results, methylation of EGFR gene body in correlation with up-regulated expression was revealed in the invasive cells. Epigenetic activation of EGFR upon resistant development to BRAF inhibitors has been previously described in melanoma, as well EGFR showed methylation difference in metastatic cell lines compared with its matched primary cell lines. [61,67] Recent studies indicated RBP4 serum levels as biomarker in colorectal cancer, and its overexpression was associated with ovarian cancer cell migration. [68][69][70] However, its function in melanoma has not been observed previously. Hypermethylation of homeobox genes is frequent in several cancers; however, this higher methylation is not consequently repress their downstream genes, as well as differentially methylated homeobox genes are not shown to be down-regulated in our invasive cells. [14,[71][72][73] Differentially methylated HOXA5 and HOXD11 was found as a specific alteration in melanoma brain metastasis, and hypermethylation of HOXD9 was described in lymph node metastasis with poorer overall survival. [14,17,47] It is suggested that methylation pattern of homeobox genes can be specific to melanoma cells, and it is a possible approach to use epigenetic biomarker panels including homeobox genes in diagnosis, prediction and prognosis. [71,74] The most interesting finding, between our results and the TCGA melanomas, is the hypermethylation of ARHGAP22 and NAV2 promoter regions that are commonly presented in locally invasive primary melanomas as well as during metastasis. Both NAV2 (neuron navigator 2) and ARHGAP22 have been identified to be involved in cell migration of different tumor types including melanoma. [75][76][77] NAV2 has several functional domains, which play key roles in the regulation of cytoskeletal remodelling and cell migration facilitating tumor invasion and metastasis. [78,79] Furthermore, a recent study suggested that NAV2 might contribute to melanoma invasion by epithelial-mesenchymal transition through the GSK-3β/β-catenin-SNAI2 pathway. [80] ARGHAP22 is a member of Rho GTPases that regulate the cytoskeleton-dependent processes during migration and invasion. [81] Silencing of ARHGAP22 results in increased number of elongated cells in melanoma cell lines and can regulate the mesenchymal-amoeboid transition. [82] The switch between mesenchymal and amoeboid types of movement allows metastatic tumor cells to adapt their morphology and movement in different microenvironments. [82] Our results indicate the relevance of methylation-mediated gene expression changes in ARHGAP22 and NAV2 during the invasion of primary tumors and also during invasion-related melanoma progression.

Recent studies indicated that increased expression of UHRF1
and/or UHRF2 negatively regulates de novo DNA methylation, and their decreased expression has been observed to correlate with hypermethylation pattern in different tumors. [83,84] Consistent with the recent findings, both UHRF1 and UHRF2 genes showed downregulation in invasive cells, and however, this mechanism need additional investigations. [61] Interestingly, the hypermethylation of

ACK N OWLED G EM ENTS
We thank Florence Le Calvez-Kelm and Geoffroy Durand from IARC's Genetics Platform (GPS) and Fabienne Barbet from the ProfileXpert platform for the processing and scanning of the Illumina 450K arrays.

CO N FLI C T O F I NTE R E S T
None declared.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data sets generated and/or analysed during the current study