DNA Methylation Profiles Distinguish Different Subtypes of Gastroenteropancreatic Neuroendocrine Tumors

Abstract

Aim: Most studies have considered gastroenteropancreatic neuroendocrine tumors (GEP-NETs) as a homogenous group of samples or distinguish only gastrointestinal from pancreatic endocrine tumors. This article investigates if DNA methylation patterns could distinguish subtypes of GEP-NETs. Materials & methods: The DNA methylation level of 807 cancer-related genes was investigated in insulinomas, gastrinomas, non-functioning pancreatic endocrine tumors and small intestine endocrine tumors. Results: DNA methylation patterns were found to be tumor type specific for each of the pancreatic tumor subtypes and identified two distinct methylation-based groups in small intestine endocrine tumors. Differences of DNA methylation levels were validated by pyrosequencing for 20 candidate genes and correlated with differences at the transcriptional level for four candidate genes. Conclusion: The heterogeneity of DNA methylation patterns in the different subtypes of gastroenteropancreatic neuroendocrine tumors suggests different underlying pathways and, therefore, these tumors should be considered as distinct entities in molecular and clinical studies.

Keywords: :

• cancer
• DNA methylation
• endocrine tumors
• growth factors
• insulinomas
• pancreatic endocrine tumors
• tumor classification

Figure 1. Hierarchical clustering of tumoral and peritumoral samples according to DNA methylation patterns.

(A) Dendrogram of all tumoral and all peritumoral samples. (B) Dendrogram of all tumoral samples. Values above each cluster correspond to the approximately unbiased value obtained after 10,000 multiscale bootstrap resampling tests. The numbers under the SIETs in the dendrogram presented in the (B) panel correspond to the two groups of SIETs identified according to the cluster analysis presented in the (A) panel (1 or 2 for the SIETs clustering respectively with the SIET peritumoral samples or the other tumors).

Pan.: Pancreatic; Duo.: Duodenal; NF-PET: Non-functioning pancreatic endocrine tumor; SIET: Small intestine endocrine tumor.

Figure 2. Main diseases, disorders, molecular and cellular functions significantly associated with each tumor type.

(A) Diseases and disorders. (B) Molecular and cellular functions. The threshold represents a p-value of 0.05.

DG: Duodenal gastrinoma; G1S: Group 1 small intestine endocrine tumors; G2S: Group 2 small intestine endocrine tumors; Ins: Insulinomas; NFP: Non-functioning pancreatic endocrine tumor; PG: Pancreatic gastrinoma; SIET: Small intestine endocrine tumor.

Figure 3. Methylation levels of four candidate genes classified according to the WHO criteria for all pancreatic endocrine tumors.

(A) Methylation level of EGF, FGF1, KLK11 and ZIM3 according to GoldenGate. (B) Methylation level of EGF, FGF1, KLK11 and ZIM3 according to pyrosequencing. Statistical significance of DNA methylation differences between the three histological groups using the GoldenGate and pyrosequencing was assessed by the nonparametric Kruskal–Wallis test.

B: Benign; UB: Uncertain behavior; WDEC: Well-differentiated endocrine carcinoma.

Figure 4. Methylation and expression patterns of four selected genes in gastroenteropancreatic neuroendocrine tumors.

(A) The DNA methylation level obtained by the BeadArray analysis in function of the tumor types. (B) Relative expression level according to tumor type.

Gast.: Gastrinomas; Ins.: Insulinomas; NF-PET: Non-functioning pancreatic endocrine tumor, Panc: Pancreatic normal/peritumoral tissue; SIET: Small intestine endocrine tumor.

Gastroenteropancreatic neuroendocrine tumors (GEP-NETs) result from alterations of neuroendocrine cells of the gastrointestinal mucosa or pancreas and secrete specific proteins used for their detection by immunohistochemistry (chromogranin A, synaptophysin, and neuron-specific enolase) [1–4]. GEP-NETs can be divided into gastrointestinal neuroendocrine tumors (GI-NETs) and pancreatic endocrine tumors (PETs), which have a different pathogenesis [5,6], display different chromosomal aberrations [7] and DNA methylation patterns at selected genomic loci [8,9] as well as distinct gene expression profiles [10]. Candidate gene approaches and genome-wide genetic and expression studies have identified genes (RAR-β, hMLH1, RASSF1, HER2, CYCLIN D1, P16^INK4a, P18^INK4c, P27^Kip1 and tyrosine kinase receptors) and pathways (e.g., PI3K/AKT/mTOR, DAXX/ATRX) potentially implicated in the tumorigenesis of GEP-NETs [10–20].

Recent advances in the genetic analysis of GEP-NETs using next-generation sequencing have demonstrated that PETs present recurrent somatic mutations in DAXX, ATRX, PTEN and TSC2 [19]. These alterations are not present in GI-NETs, in which solely recurrent somatic mutations and deletions of CDKN1B have been identified [21], strengthening the idea that these two groups of tumors develop through distinct mechanisms.

PETs can be further divided into functioning and non-functioning pancreatic endocrine tumors (NF-PETs) based on their capacity to produce peptides that give rise to clinical symptoms [6]. Functioning PETs create hormonal syndromes and include gastrinomas, insulinomas, glucagonomas, VIPomas and somatostatinomas named according to the secreted products. NF-PETs are not associated with any hormonal symptoms or syndromes even though they can secrete the same hormones as well as other peptides such as pancreatic polypeptide, neurotensin or ghrelin. PETs arise sporadically or in the context of genetic predisposition syndromes such as multiple endocrine neoplasia 1 or 2 (MEN1, MEN2), Von Hippel–Lindau (VHL), tuberous sclerosis complex (TSC) and neurofibromatosis 1 [22].

For GI-NETs, a very recent integrative study identified three molecular subgroups within small intestine endocrine tumors (SIETs) based on copy number variation analysis presenting with either no copy number variation (CNV), loss-of-heterozygosity (LOH) of only chromosome 18 or multiple CNVs. These CNV-defined groups also presented with different DNA methylation and expression profiles as well as different progression-free survival suggesting that these SIETs could be stratified by molecular features [23].

Alterations of DNA methylation patterns are a well-established and important factor in tumorigenesis [24,25]. The study of aberrant DNA methylation patterns in cancer is of great interest as a better description and comprehension of the DNA methylation patterns, which carry important diagnostic and prognostic information, might prove valuable for tumor classification and potentially individualized treatment decisions in human cancers [26]. Until recently, only few studies investigated the role of DNA methylation in GEP-NET tumorigenesis [8,9]. Gene-specific hypermethylation and global hypomethylation were found to be more pronounced in gastrointestinal endocrine tumors than in PETs [27,28]. Hypermethylation of the promoter region of a few tumor suppressor genes (CDKN2A [p14, p16], COX2, RASSF1A, hMLH1, MGMT, ESR1, RAR-β and PHLDA3) has been detected in GEP-NETs [17,27,29–33] and was sometimes associated with gene repression [17,31]. The discovery of mutations in DAXX and ATRX [19], which are part of a chromatin remodeling complex, associated to the evidence that Menin is part of an HMT complex [34,35] and that its mutation alters significantly the H3K4me3 landscape of NETs [36], show that epigenetic deregulation plays an important role in the tumorigenesis of GEP-NETs. A recently published genome-wide DNA methylation profiling study on PETs with or without loss of DAXX or ATRX expression revealed distinct DNA methylation profiles of DAXX negative and ATRX negative PETs compared with DAXX/ATRX positive PETs [37]. DAXX negative PETs showed more pronounced DNA methylation differences suggesting that the loss of DAXX and not ARTX drives DNA methylation alteration in PETs. The methylome of SIETs has also been published in a recent study identifying thousands of differentially methylated genes between tumors and normal tissue with an enrichment for MAPK, Wnt and PI3K–mTOR signaling pathways [23].

Most of the currently published studies tend, however, to consider GEP-NETs as two groups: GI-NETs and PETs, without considering potential differences within the group of PETs.

We have previously shown that the hypermethylation of the IGF2 DMR2 present in a mice model for insulinomas carrying a conditional knockout of Men1 in pancreatic β cells [38], is also present in human insulinomas [39]. The epigenetic alteration was specific for insulinomas and did not occur in any of the other types of GEP-NETs analyzed indicating the presence of specific epigenetic differences in GEP-NETs.

In the presented study, we performed a large-scale multilocus analysis of DNA methylation patterns in GEP-NETs using the Illumina GoldenGate technology [40]. We report on the DNA methylation patterns of 807 cancer-related genes in 60 GEP-NETs, including GI-NETs, SIETs, three types of PETs (insulinomas, NF-PETs and pancreatic gastrinomas) and duodenal gastrinomas. Comparison of DNA methylation patterns demonstrates that each tumor type exhibits distinct DNA methylation patterns that in some cases subclassify the tumor types further. These specific epigenetic signatures suggest different mechanisms for tumorigenesis in each tumor type, in particular PETs and demonstrate that the different types of PETs should be considered as different molecular entities to advance our understanding of the underlying biology.

Materials & methods

Ethics statement

Tissue specimens, as well as clinical and pathological records, were obtained and handled according to French laws and recommendations of the French National Committee of Ethics (Comité Consultatif National d’Ethique pour les Sciences de la Vie et de la Santé) for the protection of the patient. Written informed consent was obtained from all patients participating in this study.

Study group

A cohort of 62 patients with gastroenteropancreatic neuroendocrine tumors was analyzed to investigate differences in the DNA methylation patterns. All patients of this study were submitted to surgical resections at the Hôpital Cochin (Paris, France). During surgery, fragments of tumors obtained were immediately dissected, snap frozen and stored in liquid nitrogen until DNA and mRNA extraction. 60 tumoral samples derived from primary lesions and 30 peritumoral samples of the 62 patients including 28 paired samples were analyzed, including 18 SIETs and 11 peritumoral ileal mucosa tissues, as well as 42 PETs (Table 1). The PET samples included 11 insulinomas and five peritumoral pancreatic tissues, 12 gastrinomas (five pancreatic and seven duodenal gastrinomas) and seven peritumoral tissues and 19 NF-PETs with seven peritumoral pancreatic tissues. PETs were classified according to the WHO criteria as benign, uncertain behavior and well-differentiated endocrine carcinomas [41,42].

DNA preparation for methylation analysis

Genomic DNA was isolated from frozen tissue samples using standard phenol-chloroform extraction. DNA concentrations were determined using the Quant-iT dsDNA broad range or high sensitivity assay kit (Molecular Probes, Life Technologies) on a SpectraMAX Gemini XPS microplate spectrofluorometer (Molecular Devices, St Grégoire, France) and adjusted to a concentration of 50 ng/µl.

Bisulfite conversion of DNA was performed on 1 µg genomic DNA using the EpiTect 96 Bisulfite Kit (Qiagen) according to the manufacturer’s instruction. Five hundred nanograms of bisulfite-treated DNA was used for GoldenGate methylation assay. The remaining DNA was diluted to a final concentration of 20 ng/µl for subsequent validation studies by pyrosequencing. The absence of unconverted genomic DNA as well as the quality of the bisulfite conversion was assessed by PCR amplification and subsequent pyrosequencing analysis.

Large-scale locus-specific methylation analysis

The methylation status of 1505 CpGs located in the promoter and/or exonic regions of 807 cancer-related genes was analyzed in 91 samples (90 tumoral and peritumoral samples, and one technical replicate). The GoldenGate methylation assay was performed on 500 ng of bisulfite-treated DNA according to the manufacturer’s instructions (Illumina, Inc., CA, USA). Results were obtained as β-values ranging from 0 (unmethylated site) to 1 (methylated site). Quality control of all data was performed using the Beadstudio software (Illumina). Samples and probes with low signal intensity or less than 80% of valid results (p < 0.01 for each CpG) were removed from the analysis and the quality controlled data are provided as Supplementary Data 1.

Data analysis

The GoldenGate methylation data were analyzed using the BeadStudio software (Illumina). For the analysis of DNA methylation patterns between two groups of samples, only differences ≥0.15 (15%) (for each analyzed CpG) were considered and the compiled data are available in Supplementary Data 2. This cut-off value has previously been used in the Illumina GoldenGate BeadArray proof of concept study for the identification of differentially methylated genes as the resolution of this technology was evaluated to range from 2.5–17% depending on the interrogated CpG site [40] and is also close to the reported absolute difference in β-value able to distinguish two samples from each other (0.17, [43]).

Gene networks and gene ontology annotations were identified using Ingenuity Pathway Analysis (IPA, Qiagen). Unsupervised hierarchical clustering analysis and dendrogram representation of the methylation data was done in R 3.2.1 using Euclidean distance and the Ward method as implemented in the R package pvclust, which allows the assessment of the uncertainty in hierarchical clustering using multiscale bootstrap resampling giving approximately unbiased (AU) probabilities [44].

Differences in the distribution of DNA methylation levels between the three histology groups were assessed by a nonparametric Kruskal–Wallis tests with p < 0.05 considered as significant. Possible associations between the two main DNA methylation clusters and the histological type were assessed using a Pearson’s χ² test.

Validation of DNA methylation signatures by pyrosequencing

Validation studies were performed by pyrosequencing of bisulfite-treated DNA. Regions of interest were amplified using 20 ng of bisulfite-treated DNA and analyzed by pyrosequencing as previously described [39,45]. Sequences of amplification and pyrosequencing primers as well as the annealing temperatures used for the PCRs on bisulfite-treated DNA are given in Supplementary Table 1.

A duplicate NF-PET sample was included in the methylation profiling experiments and showed a good reproducibility shown by the R² > 0.93 between the duplicates while an R² < 0.72 was obtained when comparing it to a random NF-PET sample (Supplementary Figure 1). Moreover, four Illumina probes in four different genes were validated by quantitative DNA methylation analysis using pyrosequencing for the very same CpG site. We obtained an R² > 0.71 for all assays indicating reproducible DNA methylation patterns between the two different methods (Supplementary Figure 2). Twenty additional genes were selected from the GoldenGate methylation data for pyrosequencing validation because of either differential DNA methylation patterns between tumoral and peritumoral samples, according to their tissue origin or their specificity for a distinct tumor type (Supplementary Table 2). The methylation patterns were validated by pyrosequencing of a region including or close to the CpG analyzed in the GoldenGate experiment. For 17 of the 20 genes (85%) the initial findings were validated and the results obtained with the two techniques were found to be closely correlated (Supplementary Figure 3).

mRNA expression analysis

RNAs were isolated from frozen tissue samples using Trizol. RNA extraction was successful for four insulinomas, two gastrinomas, five NF-PETs, three SIETs and a normal pancreas. One commercial pancreatic RNA preparation (Ambion, Life Technologies, Courtaboeuf, France) was added as normal/reference sample. RNA concentrations were determined using a NanoDrop 2000c spectrophotometer (Labtech, Palaiseau, France) and RNA integrity was assessed with an Agilent 2100 Bioanalyzer (Agilent, Massy, France). The mean RNA integrity number (RIN) was 7.4 ± 1 with a coefficient of variation (CV) of 14%. One microgram of total RNA was reverse transcribed using Superscript III and oligo(dT) primers (Invitrogen) according to the manufacturer’s instructions.

Real-time quantitative PCR was carried out using the QuantiFast SYBR green PCR kit and QuantiTect Primer Assays (INS: QT01531040, PDX1: QT00201859, APC: QT02407671 and LCN2: QT00028098, Qiagen, Courtaboeuf, France) in a LightCycler 480 II thermocycler (Roche Applied Science, Meylan, France). PCR conditions were 95°C for 10 min, followed by 45 cycles of 95°C for 15 s and 60°C for 1 min. At the end of the amplification reaction, a melting curve analysis was performed to confirm the specificity as well as the integrity of the PCR product by the presence of a single peak. Absence of cross contamination and primer dimers was verified on genomic DNA and a blank water control. The geometric means of two reference genes (β-actin and β2-microglobulin), previously defined as the best combination for this sample collection, was used for normalization [39]. Standard curves were generated from assays made with serial dilutions of reference cDNA to calculate PCR efficiencies (100 ± 15%). C_t values were transformed into quantitative values using the formula (1 + efficiency)^Ct. Only means of duplicates with a CV of <15% were analyzed.

Results

In the current study, we analyzed the methylation status of 1505 CpGs located in the promoter and/or exonic regions of 807 cancer-related genes in 90 tumoral and peritumoral samples from patients with endocrine tumors of the pancreas and the small intestine. Duplicate samples present in the experiment indicated a good technical reproducibility of the method (Supplementary Figure 1) and the validation of a subset of 20 genes by pyrosequencing showed a good concordance between the two technologies and indicated a high level of confidence in the BeadArray results (Supplementary Figures 2 & 3).

DNA methylation analysis reveals hierarchical tissue & tumoral organization of samples

The quality controlled DNA methylation data were first analyzed using unsupervised hierarchical clustering and the robustness of the clusters was assessed using bootstrap resampling analysis which gave approximately unbiased (AU) probability values ranging from 69–99% for the principal subclusters of each dendrogram suggesting that the clusters were very stable and our classification reliable (Figure 1 & Supplementary Figure 4). The peritumoral samples separated into two distinct groups according to their tissue of origin, in other words, the pancreas or the intestinal mucosa as expected (Supplementary Figure 4A).

When all tumoral and peritumoral samples were included into the clustering analysis, samples separated as expected based on their tumoral or peritumoral status demonstrating widespread altered DNA methylation patterns in the tumors (Figure 1A). A subgroup of SIETs and a single gastrinoma sample (Figure 1A) showed a DNA methylation profile similar to the peritumoral samples with which they regrouped. Thus, SIETs were found to form two groups either close to the peritumoral samples (group 1) or creating a distinct subcluster within the tumoral samples (group 2) (Figure 1A). Hierarchical clustering of the tumoral samples revealed a separation of the tumors in three major clusters (Figure 1B). The first two clusters comprised the PET samples with insulinomas and gastrinomas enriched in one of each subcluster respectively, while the third cluster comprised two subclusters regrouping the SIET samples. The PETs clustered mainly together in function of the tumor type, in other words, insulinoma, gastrinoma or NF-PET (Figure 1). The two subclusters of SIETs separate perfectly the SIETs in the two groups of SIET identified in Figure 1A reinforcing the idea that these two groups of SIETs are real. These two groups of SIETs were therefore treated separately for all subsequent analyses.

Gastrinomas of pancreatic and duodenal origin grouped together when analyzing the entire dataset suggesting a similar underlying tumor biology. To reveal potential differences between the two types of gastrinomas, we performed a separate hierarchical clustering on all gastrinomas and corresponding peritumoral samples, in which most of the duodenal tumors were more similar to the peritumoral samples than the pancreatic gastrinomas, which formed a cluster apart (Supplementary Figure 4B).

Variation of DNA methylation patterns in GEP-NETs, gene ontology analysis & identification of differentially methylated gene networks

Globally, promoter regions seemed slightly more prone to methylation variation as a higher percentage of differentially methylated genes was found in promoter regions compared with exonic regions (Supplementary Tables 3 & 5). For each tumor type, the most frequent and statistically significant associated disease was as expected cancer (Figure 2). Endocrine system disorders, gastrointestinal and cardiovascular diseases were other diseases and disorders associated with the differentially methylated genes. The absence of the associated term ‘endocrine system disorders’ for NF-PETs and SIETs (Figure 2) is of note and in agreement with the absence of hormonal syndromes in these two types of tumors. Cellular growth, cell death, cellular movement and cell–cell signaling and interaction were the main molecular and cellular functions of the differentially methylated genes (Figure 2). To confirm our hypothesis of tumor type-specific DNA methylation differences, we compared the DNA methylation profile of each tumor group to their corresponding peritumoral tissues and investigated if the differentially methylated genes regrouped in networks and specific pathways. The detailed analysis of differentially methylated genes in each type of GEP-NETs compared with their respective peritumoral tissue is presented in the Supplementary Material as Supplementary Analyses, Supplementary Figures 4–13 & Supplementary Tables 3–7.

Identification of tumor type specific DNA methylation signatures in PETs

To investigate if the three types of PETs (insulinomas, pancreatic gastrinomas and NF-PETs) had divergent and specific DNA methylation patterns, comparative methylation analysis was performed. The peritumoral tissues of insulinomas, NF-PETs and pancreatic gastrinomas showed only 3–6% of differentially methylated genes (Supplementary Figure 5D & Supplementary Table 5) when compared pairwise indicating homogeneity of the samples, which was expected as they share the same tissue of origin. When compared with their peritumoral counterpart, insulinomas, NF-PETs and pancreatic gastrinomas displayed 15–29% differentially methylated genes with a comparable fraction of hypermethylated genes in the three types of PETs (8–10%) (Supplementary Figure 5C). Differences were mainly due to the number of hypomethylated genes: 7% for insulinomas, 12% for NF-PETs and 19% for pancreatic gastrinomas (Supplementary Figure 5C).

Insulinomas compared with NF-PETs or pancreatic gastrinomas (Supplementary Figure 5D & Supplementary Table 5) yielded 20 and 15% of differentially methylated genes, respectively, with most genes being hypermethylated in the insulinomas (13 and 12%, respectively). Fifteen percent of genes were differentially methylated between pancreatic gastrinomas and NF-PETs (Supplementary Figure 5D & Supplementary Table 5) without any preference for hyper- or hypomethylation.

Among the differentially methylated genes in PETs 63 were shared by all types of PETs, whereas 31, 109 and 47 genes were specific to insulinomas, pancreatic gastrinomas and NF-PETs, respectively (Supplementary Figure 5E & Supplementary Table 8). Genes specific for pancreatic gastrinomas were mostly hypomethylated and included many metalloproteinases and their inhibitors (hypomethylated: MMP1, MMP3, TIMP2, TIMP3, hypermethylated: MMP2, MMP9) as well as members of the serpin family (hypomethylated: SERPINA5, SERPINB5), whereas those specific to insulinomas and NF-PETs were equally hyper- or hypo-methylated. Specific alterations found in NF-PETs included notably three ILs (hypermethylated: IL6, hypomethylated: IL8, IL10), whereas the hypermethylation of Interferon γ (IFNG) and the hypomethylation of HDAC9 and PDX1 were specific for insulinomas.

These results demonstrate therefore that insulinomas, pancreatic gastrinomas and NF-PETs display distinct methylation patterns with specific differentially methylated genes.

DNA methylation profiling of oncogenes & tumor suppressor genes in GEP-NETs

In order to visualize the DNA methylation profiles of oncogenes and tumor suppressor genes in GEP-NETs, each probe present on the GoldenGate BeadArray was annotated using the oncogene and tumor suppressor list provided in [46] or obtained in the Uniprot Knowledge base by a keyword search [47]. Forty-nine tumor suppressor genes and 84 oncogenes were identified corresponding respectively to 81 and 136 GoldenGate probes. Most oncogenes and tumor suppressors were either methylated or unmethylated in all tumoral and peritumoral samples (Supplementary Figure 14). However, some tumor suppressor genes including SMARCB1, CASP8 and NBL1 and oncogenes including IL2, MCF2 and MOS showed hypermethylation and hypomethylation, respectively in PETs and/or SIETs. Interestingly, some oncogenes were hypermethylated in PETs including GAS7, PDGFRA and BCR while some tumor suppressor genes were hypomethylated in PETs and SIETs including APC and DLC1 (Supplementary Figure 14).

Identification of differentially methylated genes in PETs according to histology

PETs are commonly classified based on the WHO criteria into well-differentiated tumors, which can be benign (B) or of uncertain behavior (UB), in well-differentiated endocrine carcinoma (WDEC) and poorly differentiated endocrine carcinoma [41]. No associations were found between the two main clusters of tumors as defined by the DNA methylation profile and the histological classification of the tumors. However, some genes presented variation of their methylation patterns correlating with the three histological groups (Table 2) These genes include notably two growth factors (FGF1, EGF) and two imprinted genes (ZIM3 and SNRPN). Four genes (EGF, FGF1, KLK11 and ZIM3) were selected for further validation by pyrosequencing and confirmed the significant differences for all four genes found in the GoldenGate data (Figure 3). As multiple CpGs were analyzed by pyrosequencing, these analyses further refined the median and spread for DNA methylation values.

Correlation of DNA methylation with expression level of four candidate genes in GEP-NETs

In order to investigate the potential functional impact of the differentially methylated genes on GEP-NET tumorigenesis, we analyzed the gene expression level of four candidate genes (Figure 4).

We selected two genes with a consistent DNA methylation difference in all different types of tumors: APC, which is hypomethylated and LCN2, which is hypermethylated in all GEP-NETs (Figure 4A). For these two genes, an inverse correlation between DNA methylation and RNA expression level was found in most tumor types. Insulinomas exhibited only a slightly decreased expression level of LCN2 despite very high methylation levels in the promoter (Figure 4).

INS and PDX1 showed differential methylation patterns in the tumor types. The INS promoter was found to be hypomethylated in insulinomas and NF-PETs and displayed no variations in the other tumor types, whereas the PDX1 promoter was only strongly hypomethylated in insulinomas (Figure 4A). Only insulinomas overexpressed INS and PDX1 whereas the other tumor types present different level of transcriptional repression (Figure 4) recapitulating the known requirement of PDX1 for insulin expression.

These examples clearly indicate that DNA methylation profiling is a powerful tool to identify novel deregulated genes that have a functional impact on GEP-NET tumorigenesis and might provide novel targets for treatment and diagnosis.

Discussion

In the current study, we performed one of the largest analysis of DNA methylation patterns in endocrine tumors of the pancreas and the small intestine and fill thereby an important gap in the knowledge on epigenetic modifications in endocrine tumors [46]. Technological advances in DNA methylation analysis permit nowadays to obtain much more detailed genome-wide DNA methylation patterns using high-density epigenotyping arrays analyzing 450K CpG positions [48] or whole-genome bisulfite sequencing [49], which could further refine the DNA methylation landscape. However, our study is a first important step in generating knowledge on the DNA methylation differences between different subtypes of GEP-NETs. Hierarchical clustering separated samples into distinct groups mainly in function of the tumoral status (tumor or peritumoral tissue, respectively).

Previously reported studies investigating DNA methylation patterns in GEP-NETs demonstrated that global hypomethylation was more common in GI-NETs than in PETs [39,50] and that PETs display hypermethylation of some tumor suppressor genes (n = 11) at different frequencies [17]. Even if differences between SIETs and PETs demonstrated in this study might have been expected, as they have a different tissue of origin and pathogenesis [51], our analyses demonstrate the presence of different specific methylation signatures within GI-NETs and PETs. A genome-wide DNA methylation study on PETs using the Illumina 450K BeadChip described recently DNA methylation profiles in DAXX negative and ATRX positive compared with DAXX/ATRX negative PETs, but did not demonstrate tumor type-specific DNA methylation as the PETs analyzed were not further characterized into PET subgroups [37]. Therefore, our data add a new level of stratification of DNA methylation of PETs according to their tumor types in addition to the classification of DNA methylation based on the DAXX/ATRX expression status. Menin is another important protein implicated in GEP-NET tumorigenesis, which acts also as a chromatin remodeling factor [34,35] and has shown to be mutated in approximately half of PETs [19,52]. Although MEN1 was not covered by probes of the BeadArray, it is probable that gene mutation rather than DNA methylation is the molecular mechanism for MEN1 inactivation as no changes in DNA methylation patterns of MEN1 have previously been detected in GEP-NETs [27].

In our previous study, which analyzed the DNA methylation patterns of IGF2 in three different types of PETs [39], we suggested that insulinomas tend to be hypermethylated whereas NF-PETs and pancreatic gastrinomas tend to be hypomethylated. This hypothesis was confirmed in the present study at a much larger scale as insulinomas exhibited two- to three-times more hypermethylated genes compared with the other two tumor types.

The main observation of the presented work is the existence of tumor-specific DNA methylation patterns in GEP-NETs and in some cases even within the same tumor type, as exemplified by the two distinct groups of SIETs. After pathological re-examination, the latter showed no histological differences, particularly of the proliferation index or the percentage of tumor cells present in the lesions. Moreover, the clinical data of the patients revealed no association with sex or age at diagnosis. While it is highly probable from our analyses that these two distinct subtypes exist, the number of SIETs in this study is too small to draw any conclusion concerning a potentially different prognostic and requiring potentially different treatment decisions warranting further molecular characterization of these two subtypes. A very recent integrative study identified three molecular subgroups of SIETs based on CNV analyses presenting also with differences in DNA methylation, mRNA expression and progression-free survival suggesting that molecular features are a promising way to stratified further tumor types, notably GI-NETs [23]. However, we lack these molecular and clinical information for the SIETs to allow any comparison of our two subgroups with the three CNV-defined subgroups.

Gastrinomas also were found to display different DNA methylation patterns in function of their origin. Pancreatic gastrinomas were found to display four-times more alterations (mainly hypomethylated genes) of their DNA methylation patterns compared with duodenal gastrinomas. Hypomethylated genes have been associated with bad prognosis, the occurrence of metastases and/or a more aggressive phenotype for several cancer types [53–55]. The preferential hypomethylation observed in the pancreatic gastrinoma would support this idea as indeed pancreatic gastrinomas are more aggressive than duodenal gastrinomas and have a worse survival prognostic [56–58].

Our study discovered a large number of specifically methylated genes for each tumor type as well as genes that seem to be central to gene networks of differentially methylated genes. However, only integrative genome-wide analyses of different gene regulatory and transcriptional levels can yield a comprehensive vision of GEP-NETs tumorigenesis and the results should be taken with caution when only one regulation level is analyzed. Functional follow-up studies are required to determine the importance of each of these genes in the context of a specific tumor type rather than in general in PETs or GI-NETs. We have followed up the functional implications for the four selected candidate genes at the gene expression level and we found in most cases a good inverse correlation of the methylation and the expression status in the different tumor types. This indicates that analysis of variations of DNA methylation is a promising way to identify genes implicated in the GEP-NET tumorigenesis, but of course each gene requires detailed functional analyses.

Many studies have shown the simultaneous occurrence of global hypomethylation associated with tumor suppressor gene hypermethylation in cancer [59], which has also been demonstrated for GEP-NETs [27,28]. In our study, we confirmed the presence of tumor suppressor hypermethylation in PETs and SIETs as well as hypomethylation of some oncogenes although most probes did not present any difference with peritumoral tissue. Interestingly, we could also show the presence of some hypomethylated tumor suppressor genes and hypermethylated oncogene in PETs and SIETs, which suggests that these genes could have a different role in GEP-NET tumorigenesis compared with other cancers. For example, hypomethylation of APC in the different types of GEP-NETS was accompanied by an increased gene expression of this tumor suppressor gene.

The recent discovery of CDKN1B heterozygous frameshift mutations occurring in 8% of SIETs provides the first example for recurrent mutations in this type of tumors [21] as previously no recurrent mutations were found [60]. In PETs, no recurrent mutation nor promoter methylation of CDKN1B was found [19,61–62] and data on the expression status is conflicting. mRNA [62] and p27^Kip1 protein expression was sometimes found increased [13] but also reduced in MEN1 tumors [63] and in 77% of insulinomas that developed in pancreatic β-cell-specific Men1 mutant mice [64]. In our data CNKN1B did not display promoter methylation in SIETs and gastrinomas, while the promoter was strongly methylated in insulinomas and NF-PETs (20 and 35% of methylation difference compared with peritumoral samples respectively). This suggests that p27^Kip1 is predominantly inactivated in SIETs through frameshift mutations and other so far unknown mechanisms whereas in PETs it could, in addition to being regulated by Menin [63], be at least partly mediated, by DNA methylation.

Conclusion

In this study, we have demonstrated the large heterogeneity of molecular alterations in PETs or GI-NETs through the identification of six distinct tumor-specific DNA methylation patterns. Their analyses demonstrate that even pathologically homogeneous tumors like SIETs can be further subdivided into two different subtypes, but these findings require further molecular characterization and validation experiments in independent cohorts. In the light of the different underlying molecular alterations identified in the presented study, it becomes clear that a combined analysis of several tumor types as it has frequently been reported in the literature needs to be revised in future studies.

Future perspective

DNA methylation-based biomarkers could be used in diagnostics to distinguish different types or subtypes of GEP-NETs or to identify the histology of tumors when pathological examination might not be sufficient.

Functional follow-up studies are now required to determine the importance of each of the specifically or commonly altered genes in the context of a specific tumor type rather than in general in PETs or GI-NETs.

Table 1. Clinical characteristics of the patient cohort.

Clinical characteristics	Insulinomas	Gastrinomas		NF-PETs	SIETs	Total
		Pancreatic	Duodenal
Number of patients	12	5	7	20	18	62
Age
– Min	20	25	22	18	44	18
– Max	76	68	63	75	80	80
– Average	48	45	51	50	62	52
Sex
– Male	5	4	6	11	12	38
– Female	7	1	1	9	6	24
Tumors	11	5	7	19	18	60
Peritumoral samples	5	3	4	7	11	30
Number of patients with matched tumor and peritumoral samples available	4	2	4	6	11	28
Metastases
– Positive	0	4	6	14	18	42
– Negative	11	1	1	5	0	18
WHO classification
– B	5	0	0	0	–	5
– UB	6	1	1	5	–	13
– WDEC	0	4	6	14	–	24

B: Benign; NF-PET: Non-functioning pancreatic endocrine tumor; SIET: Small intestine endocrine tumor; UB: Uncertain behavior; WDEC: Well-differentiated endocrine carcinoma.

Table 2. List of genes presenting methylation variation according to the WHO criteria.

Target ID	Description	B (%)	UB (%)	WDEC (%)
CDH17^†	Cadherin 17, LI cadherin	78.39	66.67	46.04
CEACAM^†	Carcinoembryonic antigen-related cell adhesion molecule 1	51.37	27.78	16.93
EGF^†	Epidermal growth factor	69.96	56.00	23.69
FAS^†	Fas (tnf receptor superfamily, member 6)	40.25	31.76	15.98
FGF^†	Fibroblast growth factor 1	82.66	63.07	25.42
KLK11^†	Kallikrein-related peptidase 11	80.48	66.43	48.15
MKRN3^‡	Makorin ring finger protein 3	93.55	87.62	52.06
S100A2^‡	S100 calcium binding protein a2	70.50	61.08	59.95
S100A4^‡	S100 calcium binding protein a4	63.01	41.44	16.69
SNRPN^†	Small nuclear ribonucleoprotein polypeptide n	66.77	47.58	34.88
THBS2^†	Thrombospondin 2	81.75	70.07	27.20
ZIM3^†	Zinc finger, imprinted 3	76.35	73.52	39.61
PDGFRA^‡	Platelet-derived growth factor receptor, alpha polypeptide	24.04	41.04	70.80
PRKCDBP^‡	PKC, delta binding protein	11.68	26.97	62.07
TRIM29^†	Tripartite motif-containing 29	35.62	48.58	80.43

^†Promoter region analyzed.

^‡Exonic region analyzed.

B: Benign; UB: Uncertain behavior; WDEC: Well-differentiated endocrine carcinoma.

Executive summary

• DNA methylation profiling of 807 cancer-related genes in four types of gastroenteropancreatic neuroendocrine tumors (GEP-NETs) including insulinomas, gastrinomas, non-functioning pancreatic endocrine tumors (PETs) and small intestine endocrine tumors revealed DNA methylation patterns specific to each tumor type but no specific pattern according to the histology.

• Insulinomas exhibit preferentially hypermethylated genes compared with the two other tumor types of PETs (gastrinomas and non-functioning PETs) presenting more hypomethylated genes.

• Gastrinomas exhibited two distinct DNA methylation profiles depending on their pancreatic or duodenal origin.

• Two distinct groups of small intestine endocrine tumors were identified based on their DNA methylation profile and did not correlate with any known classification.

• DNA methylation differences correlated inversely with RNA expression levels for selected candidate genes.

Disclaimer

The funders had no role in study design, data collection and analysis, decision to publish or preparation of the manuscript.

Ethical conduct of research

The authors state that they have obtained appropriate institutional review board approval or have followed the principles outlined in the Declaration of Helsinki for all human or animal experimental investigations. In addition, for investigations involving human subjects, informed consent has been obtained from the participants involved.

Acknowledgments

The authors would like to thank Hafida El Abdalaoui for the bisulfite treatment of the samples, Anne Boland-Auge for the measure of DNA concentration and Doris Lechner, Damien Delafoy, Antoine Daunay and Florence Busato for excellent technical assistance.

Supplementary data

To view the supplementary data that accompany this paper please visit the journal website at: www.tandfonline.com/doi/suppl/10.2217/epi.15.85

Financial & competing interests disclosure

This work was supported by the EU framework 6 integrated project Molecular Phenotyping to Accelerate Genomic Epidemiology (MolPAGE, LSHG-CT-2004-512066) and the French Ministry of Research. The authors have no other relevant affiliations or financial involvement with any organization or entity with a financial interest in or financial conflict with the subject matter or materials discussed in the manuscript apart from those disclosed.

No writing assistance was utilized in the production of this manuscript.

Additional information

Funding

This work was supported by the EU framework 6 integrated project Molecular Phenotyping to Accelerate Genomic Epidemiology (MolPAGE, LSHG-CT-2004-512066) and the French Ministry of Research. The authors have no other relevant affiliations or financial involvement with any organization or entity with a financial interest in or financial conflict with the subject matter or materials discussed in the manuscript apart from those disclosed. No writing assistance was utilized in the production of this manuscript.