Abstract

   DNA methyltransferase inhibitors (DNMTi) decitabine and azacytidine are
   approved therapies for myelodysplastic syndrome and acute myeloid
   leukemia, and their combinations with other anticancer agents are being
   tested as therapeutic options for multiple solid cancers such as colon,
   ovarian, and lung cancer. However, the current therapeutic challenges
   of DNMTis include development of resistance, severe side effects and no
   or partial treatment responses, as observed in more than half of the
   patients. Therefore, there is a critical need to better understand the
   mechanisms of action of these drugs. In order to discover molecular
   targets of DNMTi therapy, we identified 638 novel CpGs with an
   increased methylation in response to decitabine treatment in HCT116
   cell lines and validated the findings in multiple cancer types (e.g.,
   bladder, ovarian, breast, and lymphoma) cell lines, bone marrow
   mononuclear cells from primary leukemia patients, as well as peripheral
   blood mononuclear cells and ascites from platinum resistance epithelial
   ovarian cancer patients. Azacytidine treatment also increased
   methylation of these CpGs in colon, ovarian, breast, and lymphoma
   cancer cell lines. Methylation at 166 identified CpGs strongly
   correlated (|r|≥ 0.80) with corresponding gene expression in HCT116
   cell line. Differences in methylation at some of the identified CpGs
   and expression changes of the corresponding genes was observed in TCGA
   colon cancer tissue as compared to adjacent healthy tissue. Our
   analysis revealed that hypermethylated CpGs are involved in cancer cell
   proliferation and apoptosis by P53 and olfactory receptor pathways,
   hence influencing DNMTi responses. In conclusion, we showed
   hypermethylation of CpGs as a novel mechanism of action for DNMTi
   agents and identified 638 hypermethylated molecular targets (CpGs)
   common to decitabine and azacytidine therapy. These novel results
   suggest that hypermethylation of CpGs should be considered when
   predicting the DNMTi responses and side effects in cancer patients.

   Keywords: DNA methyltransferase inhibitors, decitabine, azacytidine,
   anticancer treatment, olfactory receptor pathway, alternative splicing

Introduction

   DNA methyltransferase inhibitors (DNMTi) are widely used as chemical
   tools for hypomethylating the genome, with an aim to understand the
   role of DNA methylation in multiple processes (e.g., X-chromosome
   inactivation and DNA imprinting) and as an anti-cancer therapy
   ([25]Minkovsky et al., 2015; [26]Ramos M.P. et al., 2015; [27]Bohl et
   al., 2018). At present, two DNMTi drugs, decitabine and azacytidine,
   have been approved for treating patients with myelodysplastic syndrome
   (MDS) and acute myeloid leukemia (AML) ([28]Döhner et al., 2017;
   [29]Bohl et al., 2018), and are also being tested as therapeutic
   options in multiple solid cancers ([30]Fu et al., 2011; [31]Singal et
   al., 2015; [32]Lee et al., 2018). However, treatment response rates
   have remained very low in both solid and hematological cancers patients
   ([33]Derissen et al., 2013; [34]Nervi et al., 2015; [35]Koch et al.,
   2018). At the same time, DNMTi treatment incurs considerable side
   effects (e.g., bleeding, anemia, and joint pain). Attempts are being
   made to improve treatment efficacy (e.g., by testing combinations of
   DNMTis with other anticancer agents) and safety (e.g., by stratifying
   cancer patients into DNMTi responders and non-responders) ([36]Nervi et
   al., 2015) of the drugs. However, the precise prediction of drugs that
   show synergy with DNMTi in combination as well as distinguishing the
   drug responders from the non-responders are challenging tasks and are
   currently compromised by our incomplete knowledge of the drug
   mechanisms of action in various cancer-types ([37]Ramos F. et al.,
   2015; [38]Wang et al., 2018).

   Hypomethylation of the genome is one of the principal mechanisms behind
   the therapeutic benefit of DNMTi treatment and hypomethylated targets
   of the drugs (e.g., CpGs, genes, genomic region, pathways) have been
   functionally well-characterized in multiple cancer cell lines and other
   model systems ([39]Sarkar et al., 2013; [40]Tobiasson et al., 2017;
   [41]Koch et al., 2018). However, DNMTi induced hypomethylation of the
   genome is not sufficient to clearly understand DNMTis’ response and
   side effects ([42]Ramos F. et al., 2015; [43]Wang et al., 2018),
   suggesting presence of additional mechanisms by which DNMTis’ may
   modulate their responses. Therefore, there is a critical need to better
   understand the mechanism of action of DNMTis’ in more details,
   including their molecular targets in various cancer cells. Only a few
   studies indicated that DNMTi treatment can also lead to DNA
   hypermethylations in specific cell types, probably due to an increased
   expression of DNA methylating enzymes ([44]Broday et al., 1999;
   [45]Kastl et al., 2010; [46]Chowdhury et al., 2015). For example,
   [47]Kastl et al. (2010) reported an increase in the mRNA level of
   DNMT1, DNMT3a, and DNMT3b genes in docetaxel-resistant MCF7 cells, as
   compared to drug-sensitive cells when treated with decitabine. However,
   we currently lack the information of the genomic location and function
   of molecular targets in the cancer genome that can resist the DNA
   demethylation and show hypermethylation in response to DNMTi treatment.

   In the present work, we hypothesized that DNMTi treatment causes
   hypermethylation in the genome, and therefore we systematically
   investigated the extent of hypermethylated CpGs, their location, and
   biological role by analyzing cell line and primary tumor data. Since
   colon cancer tissue has revealed silencing of various tumor suppressor
   genes due to hypermethylation ([48]Huidobro et al., 2012), and
   especially the HCT116 cell line has widely been utilized to study DNA
   methylation and its role in regulating gene expression in colon cancer
   ([49]De Carvalho et al., 2012; [50]Yang et al., 2014), we selected
   HCT116 cell line as our primary disease model to discover novel
   hypermethylated CpGs. We validated the findings in a panel of 53 other
   cancer cell lines (e.g., lymphoma, colon, ovarian, and breast cancer),
   as well as in bone marrow cells of 2 primary leukemia patients, 2 PBMCs
   and 2 ascites samples of ovarian cancer patients. There was a
   significant correlation in the methylation of a fraction of identified
   CpGs and corresponding genes’ expression in HCT116 cell line indicating
   that the identified hypermethylated CpGs are clinically functional. Our
   work aims at initiating the search of hypermethylation in the response
   to DNMTi treatment in cancer cells and tissues as a novel mechanism by
   which one could predict in future the therapeutic and adverse effects
   of DNMTis in multiple cancer sub-types (i.e., stratified medicine), or
   even in individual cancer patients (i.e., personalized medicine).

Materials and Methods

Processing of Published Methylation Data

   To identify CpGs with increased methylation after decitabine treatment,
   we re-analyzed the DNA methylation (Illumina 450K platform,
   [51]GSE51810) data from the study of [52]Yang et al. (2014), where
   colon cancer HCT116 cells were treated with decitabine (0.3 mM) for 72
   h. Cells were maintained in McCoy’s 5A medium, supplemented with 10%
   fetal bovine serum along with 1% penicillin/streptomycin after drug
   treatment, and followed through 5, 14, 24, 42, and 68 days. The
   increase in DNA methylation in HCT116 cells were validated using
   methylation data from the study of [53]Han et al. (2013) (Illumina
   450K, [54]GSE41525), where HCT116 and T24 bladder cancer cell lines
   were treated with 0.3 and 1 μM of decitabine, respectively, for 24 h,
   and Illumina 450K assay was performed for both untreated and
   decitabine-treated cells. We also tested the increase in DNA
   methylation of identified CpGs using DMSO (as mock treatment) and
   decitabine-treated (0.06 μM for 72 h) MCF7 breast cancer cells using
   data generated by [55]Leadem et al. (2018) (Illumina 450K platform,
   [56]GSE97483). These cells were cultured in Minimum Essential Medium
   (MEM) with 10% fetal bovine serum).

   Further, we validated our findings identified in cell lines using DNA
   methylation data from primary bone marrow mononuclear cells (Illumina
   450K platform, [57]GSE20945-[58]GPL13534) from two freshly diagnosed
   acute myelogenous leukemia patients (IDs 1107 and 1307) and one healthy
   control (BM) ([59]Tsai et al., 2012). The frozen primary cells were
   thawed and cultured in Poietics HPGM supplemented with 50 ng/ML human
   thrombopoietin, 25 ng/mL stem cell factor, 50 ng/mL Flt-3 ligand, 10
   ng/mL recombinant human IL-3, 10 ng/mL recombinant human IL-6, 10 ng/mL
   recombinant human GM-CSF and 1 ng/mL recombinant human G-CSF, as
   specified in [60]Tsai et al. (2012). These cells have been treated with
   10 nM of decitabine or PBS (control) for 72 h with daily medium change.
   DNA methylation for patient 1307 was assessed at day 3 as well as day
   14 and for patient 1107 at day 7. We also analyzed genome wide DNA
   methylation (Illumina Infinium 27K Beadchip, [61]GSE31826) profile of
   14 peripheral blood mononuclear cells (PBMC), 8 tumor biopsies and 6
   ascites samples from recurrent, platinum resistant epithelial ovarian
   cancer patients, before and after treatment with decitabine ([62]Matei
   et al., 2012). Patients have been treated with decitabine (Eisai) at 10
   mg/m^2 given intravenously daily for 5 days. The methylation for each
   sample was measured before treatment at day 1 (control) and after
   decitabine treatment at day 8 (treated). The comparison of treated and
   control samples was done using the paired Wilcoxon non-parametric test
   on the basis of 16 common CpGs (out of 638) with hypermethylation.

   We extended our discovery analyses also to another DNMTi inhibitor,
   azacytidine, by re-analyzing DNA methylation data (Illumina 450K,
   [63]GSE45707) for untreated and azacytidine-treated (5 mM for 72 h)
   lymphoma U937 cell line. We further re-analyzed additional methylation
   (Illumina 450K platform, [64]GSE57342) and gene expression data
   ([65]GSE57343) from 26 breast cancer cell lines (MDA231, SKBR3, HCC38,
   ZR7530, HCC1937, CAMA1, MDA415, HCC1500, BT474, EFM192A, MDA175,
   MDA468, MDA361, HCC1954, BT20, ZR751, HCC1569, EFM19, T47D, MDA453,
   MCF7, HCC1187, HCC1419, MDA436, SUM149, and SUM159), 12 colorectal
   cancer cell lines (SW48, HCT116, HT29, RKO, SW480, Colo320, Colo205,
   SW620, SNUC-1, CACO-2, SK-CO1, and Colo201), and 13 ovarian cancer cell
   lines (TykNu, CAOV3, OAW28, OV2008, ES2, EF27, Kuramochi, OVKATE, Hey,
   A2780, OVCAR3, OVCAR5, and SKOV3) measured after mock treatment and 0.5
   μM azacytidine treatment for 72 h ([66]Li et al., 2014). The cells were
   cultured and maintained under recommended conditions for each cell line
   ([67]Li et al., 2014).

   To investigate the alteration in methylation status of identified
   probes in cancerous tissue and their role in gene expression
   regulation, TCGA level 3 HumanMethylation 450K data and normalized
   RNA-seq gene expression profiles for colon adenocarcinoma (COAD)
   patients (both tumor and adjacent healthy tissue samples) were
   downloaded using the FireBrowse tool^[68]1.

   In all the analysis, methylation status at a CpG site was measured as
   beta value (β), which is the ratio of the methylated probe intensity
   and the overall intensity (sum of methylated and unmethylated probe
   intensities designed for a particular CpG in 450K beadchip). β ranges
   from 0 to 1, indicating no methylation (β = 0) to complete methylation
   of the CpGs (β = 1). We performed appropriate quality control of the
   published data before their downstream analysis as described previously
   ([69]Giri et al., 2017). As a part of quality control, we removed all
   the CpGs with missing values and CpGs assessed by probes that have a
   tendency of cross-hybridization, as specified in the supplementary file
   of [70]Chen et al. (2013). To remove any possible bias due to design
   differences in the type of probes (the type I and type II probes)
   present in the Illumina 450K platform, we performed BMIQ normalization
   ([71]Teschendorff et al., 2013) to the DNA methylation data from the
   TCGA samples before correlation and differential methylation analysis.
   All other data processing was done using local inbuilt commands in R.

Identification of Probes With Increased Methylation in HCT116 Cell Line

   We calculated the difference in methylation level of CpGs before and
   after treatment with decitabine in HCT116 cell line at day 5 in data
   from [72]Yang et al. (2014) ([73]GSE51810). CpG that showed an increase
   in β-value of greater than or equal to 0.10 (Δβ ≥ 0.10) between
   untreated control and decitabine treated HCT116 cells were identified
   as hypermethylated CpGs. The population doubling time reported by
   [74]Yang et al. (2014) was used for correlation analysis with
   methylation level of identified CpGs. We estimated population doubling
   time as 24 h for untreated cells and 34, 38, 31, 28, 26, and 25.5 h for
   treated cells at time points of 5, 14, 24, 42, and 68 days.

Gene Expression Data Analysis

   We downloaded gene expression profiles of untreated control and
   decitabine treated HCT116 cell lines at 5, 14, 24, and 42 days
   ([75]GSE51811, Illumina HumanHT-12 V4.0 expression beadchip) from the
   same study of [76]Yang et al. (2014). The downloaded data were
   log[2]-transformed and normalized using robust spline normalization
   (RSN) method as implemented in the lumi package in R ([77]Du et al.,
   2008). Pearson correlation between methylation and gene expression
   across different time points in HCT116 cell line was assessed using the
   cor.test function of R.

   Normalized RNA-seq data of TCGA colon adenocarcinoma patients
   downloaded from FireBrowse tool^1 were further normalized using voom
   function in the limma package ([78]Ritchie et al., 2015), and these
   data were Z-transformed before the differential and correlation
   analyses. Wilcoxon non-parametric test was used to identify
   differentially methylated and differentially expressed genes between
   TCGA adenocarcinoma tumors and adjacent healthy tissue samples (FDR <
   0.05). Only those adenocarcinoma samples that had both DNA methylation
   and gene expression information were used for the correlation analyses.

   In order to find the possible reasons for absence of hypermethylation
   by low dose treatment with azacytidine in breast cancer cell lines, we
   extracted and compared the gene expression DNMT3A, DNMT3B, DNMT1, and
   DNMTL across 26 breast cancer, 12 colorectal cancer, and 13 ovarian
   cancer cell lines at day 3 (after azacytidine treatment) using gene
   expression data from [79]Li et al. (2014) ([80]GSE57343). We downloaded
   the processed gene expression data of azacytidine treated cells
   normalized to mock treated respective controls and compared across cell
   lines of different tissue types.

Gene Annotation and Pathway Enrichment Analysis

   Identified CpGs were annotated for their location in the genome based
   on annotation file provided by Illumina^[81]2. Gene ontology and
   pathway enrichment analysis of the genes corresponding to CpGs with
   hypermethylation were done using GeneCodis ([82]Tabas-Madrid et al.,
   2012). Statistical enrichment was assessed by FDR corrected p-values
   from hypergeometric test for separate ontology terms and pathways.
   GENEMANIA ([83]Zuberi et al., 2013) was used to construct and visualize
   the interaction network between genes and the transcription factor
   regulating them. Key term enrichment analysis for the genes
   corresponding to CpGs was done using DAVID ([84]Huang Da et al., 2009).

Results

DNMTi Treatment Increases Methylation of a Small Portion of the CpGs

   After quality control, we re-analyzed a total of 369886 CpGs across the
   genome of HCT116 colon cancer cells from [85]Yang et al. (2014), and
   identified hypermethylation (Δβ ≥ 0.10) of 638 CpGs in 393 unique genes
   after 5 days of decitabine treatment (0.3 mM for 72 h), as compared to
   untreated cells ([86]Figure 1A). Most of the loci were hypomethylated
   in the genome (median β = 0.18), and after decitabine treatment, a
   median increase of 0.12 (Δβ = 0.12, p < 0.0005) in methylation level
   was observed for these sites. The detailed list of the identified CpGs
   is provided in [87]Supplementary Table 1. Change in methylation at 34
   of the identified CpGs were strongly correlated (Pearson correlation
   coefficient > 0.80, FDR <0.05) with the population doubling time of
   HCT116 cell lines after decitabine treatment ([88]Supplementary Table
   2), indicating that methylation at a fraction of identified CpGs
   affects proliferation and growth of cancer cells. However, most of the
   identified sites loss their hypermethylation by day 10 ([89]Figure 1)
   suggesting that the observed hypermethylation is transient. Re-analysis
   of another methylation data for HCT116 cell line from the [90]Han et
   al. (2013) study validated our finding, as we found a corresponding
   increase in methylation level (median Δβ = 0.09) at 583 common CpGs
   after decitabine treatment (0.3 μM for 24 h) ([91]Figure 1B). These
   results indicate that the increase in DNA methylation at most of the
   identified sites starts as early as 24 h after the DNMTi treatment and
   lasts up to at least day 5. The result suggests that there are CpGs
   that not only resist the demethylation in response to DNMTi but also
   show transient hypermethylation.

FIGURE 1.

   [92]FIGURE 1
   [93]Open in a new tab

   Decitabine treatment increases DNA methylation levels of a subset of
   CpGs. (A) Scatter plots showing DNA methylation patterns of 638
   differentially methylated CpGs between untreated control cells (x-axis)
   and decitabine-treated cells (y-axis) at various time points in the
   study of [94]Yang et al. (2014). (B) Violin plot showing the median
   methylation level (horizontal line) and distribution patterns (density
   and IQR) of the identified 583 CpGs in untreated and 0.3 μM
   decitabine-treated HCT116 cells after 24 h in the study of [95]Han et
   al. (2013). The statistical significance was assessed using the
   non-parametric Wilcoxon test. ^∗∗∗p < 0.0005.

   Further, we also tested the effect of decitabine on identified
   differentially methylated CpGs in a bladder cancer cell line (T24). An
   increase in median DNA methylation levels (median Δβ = 0.14, p <
   0.0005) at 616 common CpGs was observed after the drug treatment (1 μM
   for 24 h) in T24 cells ([96]Figure 2A) in contrast to a significant
   decrease in the methylation level of other CpGs present in the 450K
   beadchip (median Δβ = −0.14) as shown in [97]Supplementary Figure 1.
   However, we did not observe any increase in methylation level of 590
   common identified CpGs (median Δβ = −0.01, p < 0.0005) in breast cancer
   MCF7 cell line treated with 0.06 μM of decitabine for 72 h ([98]Figure
   2B). Replication of the results in multiple cancer cell lines indicates
   that hypermethylation in the cancer genome is a common effect of
   decitabine treatment that can contribute to DNMTis’ response.

FIGURE 2.

   [99]FIGURE 2
   [100]Open in a new tab

   Increase in methylation of identified CpGs is cell line-specific. (A)
   Methylation level of 616 identified probes common in untreated and
   decitabine-treated (1 μM for 24 h) bladder cancer T24 cell line. An
   increase in median DNA methylation levels (median Δβ = 0.14) at 616
   common CpGs was observed after the drug treatment. (B) Methylation
   level of 590 identified probes common in mock (DMSO) and
   decitabine-treated (0.06 μM for 72 h) breast cancer MCF7 cell line. We
   observed a smaller decrease in methylation level (median Δβ = –0.01, p
   < 0.005) when the cells were treated with lower dose (0.06 μM) of
   decitabine in MCF7 cell lines. The statistical significance was
   assessed using the non-parametric Wilcoxon test. ^∗∗∗p < 0.0005.

   To study the effect of another DNMT inhibitor, azacytidine, on
   methylation of identified CpGs, we next re-analyzed 450K data from 52
   pan-cancer cell lines. We found a significant increase in methylation
   of identified CpGs in 9 out of 13 ovarian cancer cell lines (median Δβ
   > 0.03), 3 out of 26 breast cancer cell lines (median Δβ > 0.01), and
   in 9 out of 12 colon cancer cell lines (median Δβ > 0.02), and in 1
   lymphoma cell line (U937, median Δβ > 0.06) after azacytidine treatment
   ([101]Figure 3A–D). Validation of the observed hypermethylation at
   identified target CpGs in multiple cancer types (bladder, ovarian,
   breast, and lymphoma cancer cell lines) suggests that hypermethylation
   in the genome is a shared phenomenon in response to multiple DNMTi
   agents and that the degree of the hypermethylation is cancer- and cell
   line-specific.

FIGURE 3.

   [102]FIGURE 3
   [103]Open in a new tab

   Azacytidine treatment increases methylation of identified sites in a
   subset of cell lines. Change in median methylation level of identified
   CpGs in (A) 13 ovarian cancer cell lines (B) 26 breast cancer cell
   lines, and (C) 12 colon cancer cell lines. The bar plot represents the
   difference in the median methylation level of identified CpGs between
   cells treated with 0. 5 μM azacytidine (test group) or carboplatin
   (mock group) after 72 h. (D) Scatter plot showing the distribution of
   methylation level of identified probes in untreated (control) and
   treated cells (5 mM of azacytidine for 72 h) in U937 lymphoma cell
   lines. The statistical significance was assessed using the
   non-parametric Wilcoxon test. ^∗p < 0.05, ^∗∗p < 0.005.

Identified CpGs Show Trend of Hypermethylation in AML and Ovarian Cancer
Primary Cells

   To investigate the CpGs identified in the HCT116 colon cancer cells in
   primary patient samples, we analyzed hypermethylation of identified
   CpGs in mononuclear cells from bone marrow sample of two primary AML
   patients ([104]Figure 4). Increase in DNA methylation (Δβ = 0.06, p =
   4.68 × 10^−7 at day 7) was observed at the identified CpGs after
   decitabine treatment in the bone marrow mononuclear cells of one
   patient (ID 1107) but not of the healthy control (median Δβ = −0.02, p
   = 0.26 at day 3), suggesting that hypermethylation in response to DNMTi
   might be specific to cancerous cells but not to healthy cells
   ([105]Figure 4). We also observed trend of hypermethylation at
   identified CpGs (median Δβ = 0.06, p = 0.76) at day 3 in another
   primary leukemia patient (ID 1307) that was lost by day 14 (median Δβ =
   −0.01, p = 2.13 × 10^−35) ([106]Figure 4). This result again suggests
   that hypermethylation at identified sites after DNMTi treatment is
   transient also in primary leukemia cancer cells, similarly as was
   observed in the HCT116 cell line model. Further, significant
   hypermethylation was observed in 2 out of 14 PBMCs (P2 and P3) and 2
   out of 6 ascites (A1 and A2) samples of platinum resistant epithelial
   ovarian cancer patients treated with decitabine as compared to
   untreated paired samples ([107]Supplementary Figure 2). Validation of
   the hypermethylation at identified CpGs in the PBMCs and ascites
   samples treated with clinically used dose of decitabine (10 nM) further
   confirm our findings. We did not either observe hypermethylation in any
   of the eight tumor biopsy samples ([108]Supplementary Figure 2) from
   the ovarian cancer patients which might be because of the requirement
   of higher dose of drug to result in the same effect in solid cancer as
   compared to PBMC and ascites samples.

FIGURE 4.

   [109]FIGURE 4
   [110]Open in a new tab

   Pirate plot show the median methylation level (horizontal colored line)
   and distribution pattern (mean ± standard error as vertical black line)
   of the identified 638 CpGs in untreated and treated (10 nM of
   decitabine) mononuclear cells from a healthy bone marrow (BM) sample
   and two primary leukemia patient samples (IDs 1307 and 1107), based on
   the data from the study of [111]Tsai et al. (2012). DNA methylation
   assay for the patient 1307 was performed at days 3 and 14. For the
   patient 1107, the assay was performed only at day 7. The median
   methylation level is indicated in each plot. The statistical
   significance was assessed using the paired non-parametric Wilcoxon
   test. NS, not significant. ^∗∗∗p < 0.0005.

The Methylation of Identified CpGs Correlates With Expression of
Corresponding Genes

   To discover the functional role of methylation in differentially
   methylated CpGs, we correlated DNA methylation level of identified CpGs
   and expression levels of the corresponding genes at multiple time
   points (0, 5, 14, 24, and 42 days) measured in HCT116 cell line after
   decitabine treatment. The distribution of correlation coefficients for
   130 CpGs in promoter and 200 CpGs in gene body area has been shown in
   [112]Figure 5A. A strong correlation (|r|≥ 0.80) between the gene
   expression and DNA methylation profiles was observed at 27% of loci
   (166 CpGs in 153 genes) in HCT116 cells, out of them, 47 CpGs
   corresponding to 44 genes were nominally significant (p < 0.05)
   ([113]Supplementary Table 3) indicating that observed hypermethylations
   affect cellular function by altering gene expression in cancer cell
   lines. The correlation between gene expression and DNA methylation
   suggests that the increase in DNA methylation is functional and causes
   transcriptional changes that affects cellular processes (e.g., cell
   proliferation).

FIGURE 5.

   FIGURE 5
   [114]Open in a new tab

   DNA methylation at hypermethylated sites affects gene expression in
   HCT116 cell lines and TCGA colon adenocarcinoma samples. (A) Histogram
   depicts density of Pearson correlation coefficients (PCCs) between DNA
   methylation and gene expression for 130 CpGs in promoter region (left
   panel) and 200 CpGs in gene body area (right panel). The PCCs were
   calculated using DNA methylation of a CpGs and the expression level of
   corresponding genes at day 0 (untreated cells), 5, 14, 24, and 42 days
   after decitabine treatment. Further, 53 promoter CpGs and 78 gene body
   CpGs showing a strong methylation-expression correlation (|r| ≥ 0.8) in
   HCT116 cell lines were selected and their methylation levels were
   compared between cancerous (n = 273) and adjacent healthy tissues (n =
   19) samples in TCGA colon cancer data. A significant difference in
   methylation level (FDR < 0.05) was observed at 38 promoter and 73 gene
   body CpGs. (B) Heatmap showing the methylation level of these 38
   promoters (left panel) and 73 gene body (right panel) CpGs in healthy
   tissue and colon cancer tissue. Subsequently, the expression level of
   genes corresponding to CpGs showing significant methylation differences
   in TCGA colon cancer data were compared and 19 genes corresponding to
   significant CpGs in promoter region, and 38 genes corresponding to the
   significant CpGs in gene body region were found differentially
   expressed (FDR < 0.05). (C) Heatmap showing the expression profile of
   genes with significant differences in adjacent healthy tissue (n = 19)
   and colon cancer samples (n = 273).

Genes Corresponding to the Identified CpGs Show Differential Expression in
Cancerous Tissue

   We further investigated the pathological relevance of the increased
   methylation at identified CpGs using DNA methylation and gene
   expression profile of normal and cancerous tissue from TCGA colon
   adenocarcinoma patients. Out of 638 CpGs discovered in HCT116 cell
   line, we extensively investigated the methylation of CpGs in promoter
   and gene body area in TCGA data as methylations at these locations have
   a major effect on gene expression ([115]Yang et al., 2014; [116]Koch et
   al., 2018). We selected 53 promoter CpGs and 78 gene body CpGs showing
   a strong methylation-expression correlation (| r|≥ 0.8) in HCT116 cell
   lines and their methylation levels were compared between cancerous (n =
   273) and adjacent healthy tissues (n = 19) samples in TCGA colon cancer
   data. This analysis revealed differential methylation of 38 promoter
   CpGs and 73 gene body CpGs in cancer tissue ([117]Figure 5B and
   [118]Supplementary Table 4). Subsequently, the expression level of
   genes corresponding to CpGs showing significant methylation differences
   in TCGA colon cancer data were compared and the analysis revealed
   significant difference (FDR < 0.05) in expression of 19 and 38 genes
   corresponding to differentially methylated CpGs at promoter and gene
   body area in colon cancer tissue respectively ([119]Figure 5C and
   [120]Supplementary Table 5). The analysis suggests that DNMTi treatment
   increases methylation at CpGs that involved in the pathology of cancer
   and hence such alteration can contribute to the therapeutic response of
   the drugs.

Genes With Increased Methylation Are Enriched in Cancer-Related Pathways and
Are NFAT, LEF1, and MAZ-Regulated

   We next investigated the functions of the genes corresponding to the
   identified CpGs using the GeneCodis (v2) gene set enrichment analysis
   tool. Gene ontology enrichment analysis revealed that five out of the
   10 (50%) most significant GO processes (FDR = 0.05) were related to
   transcription regulation, which is one of the key functional role of
   DNA methylation in order to control gene expression ([121]Figure 6A).
   Notably, the list of enriched genes included well-known oncogenes, such
   as AFF3, CTNND2, ELK4, ESR1, PAX3, TRRAP, and WHSC1L1. The pathway
   enrichment analysis revealed that olfactory transduction and p53
   signaling pathway were overrepresented (FDR = 0.05) in the gene set
   with increased methylation ([122]Figure 6B), and alternating splicing
   was enriched as the major keyword (fold enrichment = 1.27, FDR =
   0.000132, [123]Figure 6C). Identified CpGs were also enriched in the
   enhancer region (p = 3.36 × 10^−4) of the genomes ([124]Supplementary
   Figure 3).

FIGURE 6.

   [125]FIGURE 6
   [126]Open in a new tab

   Gene ontology (GO) and pathway enrichment analyses of genes
   corresponding to differentially methylated CpGs. (A) Top-10 most
   significantly enriched cellular processes are shown as a bar plot. The
   lengths of the bars denote the number of genes present in each of the
   top GO categories. (B) Pie-chart showing the significantly enriched
   pathways for the genes. The number of genes present in each pathway
   group is shown along with the hypergeometric test p-value corrected for
   multiple testing as implemented in the DAVID tool. (C) The keywords
   enrichment analysis for the genes is shown as a bar chart. The length
   of the bar represents the number of genes enriched for each keyword,
   the FDR-corrected p-value is shown at the top. ^∗FDR < 0.05; ^∗∗FDR <
   0.005; ^∗∗∗FDR < 0.0005.

   Furthermore, enrichment analysis among the transcription factors
   regulating these identified genes revealed that 59 genes were regulated
   by lymphoid enhancer-binding factor 1 (LEF1) (p = 2.1 × 10^−8,
   [127]Table 1 and [128]Supplementary Figure 4A), 51 genes by MYC
   Associated Zinc Finger Protein (MAZ) (p = 4.07 × 10^−8, [129]Table 1
   and [130]Supplementary Figure 4B) and 47 of genes were regulated by the
   nuclear factor of activated T-cells (NFATC1, p = 1.48 × 10^−8,
   [131]Supplementary Figure 4C). Enrichment of genes in cancer-related
   pathways, especially in the p53 tumor suppressor pathway and the
   olfactory pathway, suggests that hypermethylation in the identified
   sites affects cancer cell proliferation and apoptosis.

Table 1.

   Top five transcription factor enriched in gene set corresponding to
   identified 638 CpGs and its regulated genes.
   Transcription factor Genes (N) FDR Regulated genes
   NFATC1 47 1.48 × 10^−8 SLIT3, RORA, SYT10, CNTNAP2, SCN3A, MRPL28,
   ADAMTSL1, CTNND2, ESR1, PPM1B, PDE4D, RARB, OPCML, SNX15, PRKD2, ID3,
   PTPRO, ADCY2, CUL3, DMD, ITM2C, KLF12, BCOR, CTNND1, SGCD, ACACA,
   HDAC6, POGK, AUTS2, PAX3, DLG2, SLC6A5, SOX5, DLC1, ANTXR1, NGFRAP1,
   LSAMP, GRM8, CACNA2D3, ETS1, S100A10, ADAMTS17, KCNH5, ARHGAP6, KCNMA1,
   MAP7, KCNN2
   LEF1 59 2.10 × 10^−8 PDCD10, NRXN1, SLIT3, RORA, YWHAZ, COX7B, SYNPR,
   SCN3A, SORCS1, TMSB4X, ADAMTSL1, NXN, CLSTN2, ZNF8, CNTN6, MAGED2,
   WHSC1L1, GMPR2, PDE4D, ABCF2, RARB, CNKSR2, TIA1, SMARCA1, SFRP2,
   OPCML, WDFY3, MBTD1, CACNA1E, PTPRO, MCTS1, DMD, KLF12, CTNND1, SGCD,
   ACACA, POGK, OXCT1, TAF1, PAX3, SLC6A5, SOX5, 1, SIX4, GPC6, DLC1,
   GTF2A2, TLE3, GAB2, TMSL3, CACNA2D3, ETS1, BZW1, ADAMTS12, CD160,
   TCERG1L, KCNH5, ARHGAP6, NR2F1
   MAZ 51 4.07 × 10^−8 PDCD10, SLIT3, SV2B, RORA, YWHAZ, CNTNAP2, SORCS1,
   PRKCI, CTNND2, ACCN2, ESR1, MAGED2, HRK, FKBP2, RARB, CNKSR2, SSR1,
   SMARCA1, SFRP2, PRKD2, UBE2L3, PTPRF, ID3, DMD, ITM2C, KLF12, P4HA1,
   BCOR, RGS7, CTNND1, DUSP6, POGK, TAF1, LYPLA2, AUTS2, DLG2, SOX5, SIX4,
   DLC1, TLE3, PRDM16, NGFRAP1, POLR1D, ARVCF, BZW1, THRAP3, TRRAP, ZNRF1,
   KCNH5, KCNMA1, NR2F1
   OCT1 16 1.49 × 10^−7 NRXN1, NR2C2, SCN3A, PDE4D, RARB, DMD, KLF12,
   DUSP6, SOX5, DLC1, TLE3, GAB2, TMSL3, TCERG1L, IRX2, NR2F1
   MEF2 27 2.32 × 10^−7 HIPK1, CTNNA1, NRXN3, ESR1, PPM1B, SLC12A1,
   SMARCA1, OPCML, WBP5, ADCY2, CUL3, DMD, KLF12, P4HA1, CTNND1, SGCD,
   DLG2, SOX5, GLG1, NELL2, GRM8, CACNA2D3, ETS1, CIAO1, ADAMTS12,
   AGTPBP1, NR2F1
   [132]Open in a new tab

   FDR-corrected p-values from the hypergeometric test are shown in the
   table.

Discussion

   Our study indicates that a clinically feasible dose of decitabine (0.3
   to 300 μM) treatment causes a transient increase in DNA methylation
   level of a small fraction of CpGs in the colon cancer genome. The use
   of 3-day exposure with such doses in vitro produces a quick increase in
   DNA methylation that may reflect the immediate response of cells to
   external stimuli. The transient increase in methylation (that vanishes
   10 days after treatment) at identified targets also suggests that the
   response is mainly a temporary alteration of transcript level of
   certain genes, rather than permanent shut down or enhancement of their
   expression.

   The further study of the hypermethylation in primary bone marrow
   mononuclear cells from AML patients, peripheral blood cells and ascites
   from epithelial ovarian cancer patients indicated that hypermethylation
   of subsets of the identified CpGs is an important mechanism behind
   DNMTi effects against primary cancer cells, which, as expected, varies
   across patients and cancer types. Additionally, the absence of
   hypermethylation in normal bone marrow mononuclear cells indicates that
   the hypermethylation is limited to cancerous tissues and cell lines.
   However, a future study exploring hypermethylation at identified CpGs
   in pre- and post-treated cancerous tissue types (e.g., kidney, skin,
   brain), in comparison to corresponding healthy cells, is needed to
   detail the tissue-specificity of the hypermethylation at the identified
   CpGs.

   We did not observe hypermethylation in case of breast cancer MCF7 cell
   lines treated with low dose (0.06 μM) of decitabine, nor in 20 out of
   26 breast cancer cell lines treated with 5 mM of azacytidine (low
   dose), even though the same dose of azacytidine was sufficient to
   hypermethylate the CpGs in colon and ovarian cancer cell lines. It
   suggests that lower doses of DNMTi agents are not sufficient to
   hypermethylate the identified CpGs in all the cancer cell lines, which
   could be because of their molecular peculiarities. For example, lower
   change in expression level of DNA methylating enzymes DNMT3A was
   observed in breast cancer cell lines 3 days after azacytidine
   treatment, as compared to colon and ovarian cancer cell lines
   ([133]Supplementary Figure 5). DNMT3A is a de novo methylating enzyme
   and have been associated with CpGs hypermethylation phenotype in AML
   ([134]Spencer et al., 2017). Further, the expression level of DNMT3A
   and DNMT3B have been associated with sensitivity to decitabine
   cytotoxic response in embryonic cells ([135]Oka et al., 2005),
   testicular germ cell tumor ([136]Beyrouthy et al., 2009) and have also
   shown to be predictive of decitabine treatment response in
   patient-derived xenograft model for triple negative breast cancer
   ([137]Yu et al., 2018). Differences in expression of DNMT enzymes,
   together other molecular differences (e.g., genetic mutations or
   epigenetic differences) may resists breast cancer cells against the
   hypermethylation in response to DNMTi at a lower dose. However, a
   detailed mechanistic study will be needed to understand the absence of
   hypermethylation in breast cancer cells in relation with DNMT3A or
   DNMT3B levels and other molecular heterogeneity.

   Additionally, genes associated with identified CpGs were enriched among
   olfactory receptor and p53 pathways. P53 pathway regulates DNA repair
   and apoptosis in HCT116 cells ([138]Li et al., 2015), and there is
   growing evidence for the involvement of olfactory receptor pathway in
   cancers of non-olfactory tissues such as prostate ([139]Neuhaus et al.,
   2009), lung ([140]Giandomenico et al., 2013; [141]Kalbe et al., 2017),
   and colon cancers ([142]Morita et al., 2016). For instance, olfactory
   receptor OR51B4 is highly expressed in colon cancer HCT116 cells and
   the activation of the receptor by its ligand (Troenan) inhibited cell
   proliferation and apoptosis by phosphorylation of p38, mTOR, and Akt
   kinases ([143]Weber et al., 2017). The available literature therefore
   suggests that identified hypermethylation can affect colon cancer cells
   proliferation and apoptosis through olfactory receptor pathway.
   However, a detailed mechanistic study is warranted to decipher the
   exact role of hypermethylation at identified CpGs in cellular growth
   through p53 and olfactory receptor pathways.

   Further, enriched transcription factor LEF1 is implicated in
   tumorigenesis and cancer cell proliferation, migration, invasion, and
   stemness in multiple cancer (e.g., colorectal cancer, AML, oral
   squamous cell carcinoma) ([144]Santiago et al., 2017). Moreover, the
   methylation at a fraction of CpGs was positively correlated with the
   population doubling time of HCT116 cells indicating the role of
   identified CpGs in cell division, and proliferation. Enrichment of
   genes among the cancer-related pathways involved in cell proliferation
   and apoptosis suggests a non-random selection of CpGs for
   hypermethylation in order to alter the cellular dynamics (e.g.,
   population doubling time). We consider that increased in methylation at
   identified CpGs can contribute to the drug response in multiple cancer
   types partly by altering the cell division rate. The methylation levels
   at identified CpGs could be useful to design better treatment
   strategies (e.g., predicting the DNMTi responders and non-responders),
   and to understand side effects due to higher drug dose (e.g., the
   increase in DNA methylation at identified CpGs could be an indicator
   for degree of possible toxic effect). These novel findings therefore
   warrant further analyses in a number of cancer types using either
   pre-clinical animal studies or clinical treatment data from cancer
   patients with baseline molecular profiles available.

   There are multiple reasons to suggest that at least one key mechanism
   underlying the anti-tumor responses to DNMTi treatment may involve
   hypermethylation of specific genes in cancer cells. First, we showed an
   increase in DNA methylation at identified CpGs in more than one type of
   cancer cells. Second, as defined for an epigenetic change, these
   sustained changes persist for significant periods of time (at least
   more than 5 days) after a transient, subsequently withdrawn, drug
   exposure (in this case 72 h). Third, the expression patterns for a
   subset of the genes are different between cancer and normal tissue
   types. Importantly, these changes are induced by drug doses that do not
   acutely kill cells and, thus, allow the transient alterations in gene
   methylation patterns to act on emerging molecular changes to cells
   after DNMTi therapy.

   In summary, our findings suggest that DNMTi mediates its therapeutic
   effect through a complex mechanism of action, and therefore a
   generalized pattern for the activity is challenging to find as DNMTi
   treatment causes hypermethylation at certain loci and hypomethylation
   at others. Hence, the effects of DNMTi on cancer tissues should be
   analyzed at the individual gene level, rather than at the entire
   genomic level, and separately for each tissue type and even for each
   cancer patient. Further, any attempt to predict DNMTi response should
   also incorporate the hypermethylated targets of the drugs along with
   the hypomethylated targets of the drugs. However, a functional
   mechanistic study in more advanced model systems and several human
   tissue types will be required for further revealing how the increase in
   methylation at identified loci alters treatment response and the
   pathological burden of disease. We hope that these identified targets
   associated with hypermethylation have broad implications for further
   research on DNMTi response mechanisms in multiple cancer types.

Author Contributions

   AG conceptualized, analyzed the data, and wrote the manuscript. TA
   critically revised and edited the manuscript.

Conflict of Interest Statement

   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.

Acknowledgments