Aberrant DNA methylation of WNT pathway genes in the development and progression of CIMP-negative colorectal cancer

ABSTRACT

The WNT signaling pathway has an essential role in colorectal carcinogenesis and progression, which involves a cascade of genetic and epigenetic changes. We aimed to analyze DNA methylation affecting the WNT pathway genes in colorectal carcinogenesis in promoter and gene body regions using whole methylome analysis in 9 colorectal cancer, 15 adenoma, and 6 normal tumor adjacent tissue (NAT) samples by methyl capture sequencing. Functional methylation was confirmed on 5-aza-2′-deoxycytidine-treated colorectal cancer cell line datasets. In parallel with the DNA methylation analysis, mutations of WNT pathway genes (APC, β-catenin/CTNNB1) were analyzed by 454 sequencing on GS Junior platform. Most differentially methylated CpG sites were localized in gene body regions (95% of WNT pathway genes). In the promoter regions, 33 of the 160 analyzed WNT pathway genes were differentially methylated in colorectal cancer vs. normal, including hypermethylated AXIN2, CHP1, PRICKLE1, SFRP1, SFRP2, SOX17, and hypomethylated CACYBP, CTNNB1, MYC; 44 genes in adenoma vs. NAT; and 41 genes in colorectal cancer vs. adenoma comparisons. Hypermethylation of AXIN2, DKK1, VANGL1, and WNT5A gene promoters was higher, while those of SOX17, PRICKLE1, DAAM2, and MYC was lower in colon carcinoma compared to adenoma. Inverse correlation between expression and methylation was confirmed in 23 genes, including APC, CHP1, PRICKLE1, PSEN1, and SFRP1. Differential methylation affected both canonical and noncanonical WNT pathway genes in colorectal normal-adenoma-carcinoma sequence. Aberrant DNA methylation appears already in adenomas as an early event of colorectal carcinogenesis.

KEYWORDS:

• Adenoma
• β-catenin mutation
• APC
• colorectal cancer
• promoter and gene body methylation
• gene expression
• whole methylome analysis
• WNT signaling pathway

Introduction

Colorectal cancer (CRC) is one of the most frequent cancers world-wide with an annual incidence of approximately 1 300 000 newly diagnosed cases, and a very high global mortality of more than 690 000/year.¹ Although more and more CRC-associated molecular alterations are determined, epidemiological data emphasize the necessity of further research of this disease. CRC can develop through several distinct molecular pathways (adenoma-carcinoma sequence, serrated neoplastic pathway or on the basis of long-standing inflammatory bowel disease) involving a cascade of genetic, epigenetic events affecting gene expression changes.^2-4 The analysis of these pathways can improve our understanding of colorectal carcinogenesis and can reveal potential diagnostic, prognostic, and therapeutic predictive markers.

Abnormalities of epigenetic regulation such as DNA methylation alterations can contribute to malignant transition.^5-6 It is well known that global hypomethylation and promoter hypermethylation of several tumor suppressor genes are characteristic in cancers.⁷ Furthermore, with the advent of whole methylome analysis, several studies allude to the important role of methylation alterations in the transcribed regions of genes (gene body region)⁷ with a remarkable proportion of differentially methylated CpG sites located there^5-6 and influencing the regulation of gene expression.^8-9

The WNT signaling pathway (Table 1, Supplementary Table 1) has an essential role in colorectal carcinogenesis and progression, which is characterized by accumulation of genetic and epigenetic changes. The most known WNT pathway alterations, such as APC and AXIN2 inactivation by mutations or loss of heterozygosity, affect the APC-GSK3-Axin complex leading to β-catenin stabilization, thereby resulting in the constitutive activation of canonical WNT signaling in CRC.¹⁰ Promoter hypermethylation of negative regulators of canonical WNT pathway like secreted frizzled-related proteins (SFRPs),^11-18 dickkopf family proteins (DKKs),^11,17,19-21 and WIF1 WNT inhibitory factor^11,17,21 has also been described in malignancies, including CRC. Methylation patterns of CRC and normal adjacent tissue samples were compared in most of the previous global methylation studies focusing on WNT signaling pathway.^20,22 Colorectal adenoma samples were involved in methylation analyses for selected WNT pathway genes^{11,15,17,23,24} focusing mainly on promoter methylation level alterations. Serrated adenoma samples—representing another neoplastic pathway of CRC development—were also profiled in several previous methylation projects.^23,25-27

Table 1. Analyzed WNT pathway genes.

WNT pathway	Gene symbol
Canonical pathway	AES, AKT1^*, AKT2^*, AKT3^*, APC, APC2, AXIN1, AXIN2, BAMBI, BTRC, CACYBP, CCND1, CCND2, CCND3, CER1, CHD8, CREBBP, CTBP1, CTBP2, CTNNB1, CTNNBIP1, CUL1, CXXC4, CSNK1A1, CSNK1A1L, CSNK1E, CSNK2A1, CSNK2A2, CSNK2B, DKK1, DKK2, DKK4, DVL1, DVL2, DVL3, EP300, FBXW11, FOSL1, FRAT1, FRAT2, FZD1, FZD10, FZD2, FZD3, FZD4, FZD5, FZD6, FZD7, FZD8, FZD9, GSK3B, JUN, LEF1, LRP5, LRP6, MAP3K7^*, MMP7, MYC, NKD1, NKD2, NLK, PORCN, PPARD, PPP2CA, PPP2CB, PPP2R1A, PPP2R1B, PPP2R5A, PPP2R5B, PPP2R5C, PPP2R5D, PPP2R5E, PRKACA, PRKACB, PRKACG, PRKX, PSEN1, RBX1, RUVBL1, SENP2, SFRP1, SFRP2, SFRP4, SFRP5, SIAH1, SKP1, SMAD2^*, SMAD3^*, SMAD4^*, SOX17, TBL1X, TBL1XR1, TBL1Y, TCF7, TCF7L1, TCF7L2, TLE1, TLE2, TLE3, TLE4, TLE6, TP53^*, WIF1, WNT1, WNT10A, WNT10B, WNT11, WNT16, WNT2, WNT2B, WNT3, WNT4, WNT5A, WNT5B, WNT6, WNT7A, WNT7B, WNT8A, WNT8B, WNT9A, WNT9B
Planar cell polarity pathway	DAAM1, DAAM2, DVL1, DVL2, DVL3, FZD1, FZD10, FZD2, FZD3, FZD4, FZD5, FZD6, FZD7, FZD8, FZD9, GPC4, MAPK10^*, MAPK8^*, MAPK9^*, NKD1, NKD2, PRICKLE1, PRICKLE2, RAC1, RAC2, RAC3, RHOA, ROCK1, ROCK2, VANGL1, VANGL2, WNT11, WNT5A
WNT/Ca^2+ pathway	CAMK2A, CAMK2B, CAMK2D, CAMK2D, CHP1, CHP2, FZD1, FZD10, FZD2, FZD3, FZD4, FZD5, FZD6, FZD7, FZD8, FZD9, NFAT5, NFATC1, NFATC2, NFATC3, NFATC4, PLCB1, PLCB2, PLCB3, PLCB4, PPP3CA, PPP3CB, PPP3CC, PPP3R1, PPP3R2, PRKCA, PRKCB, PRKCG, WNT5A, WNT5B

Official gene symbols of WNT pathway genes according to KEGG WNT signaling map. Detailed data for these genes can be seen in Supplementary Table 1.

* WNT signaling pathway genes in other intracellular pathways, such as PI3K-AKT, TGFB, P53, and MAPK signaling. Genes with differentially methylated promoters in CRC vs. NAT comparison are indicated in bold; genes with differentially methylated promoters in adenoma vs. NAT comparison are underlined.

Illumina BeadChip technology is a frequently used method in global methylation studies^22,28,29 and allows the determination of the methylation status of >485,000 CpG sites at single-nucleotide resolution. Next generation sequencing of the methylated DNA regions, such as methyl capture sequencing (MethylCap-seq), is another option for genome-wide methylation analysis for revealing novel differentially methylated regions (DMRs).^30,31 MethylCap-seq provides a more extensive overview of whole-genome methylation than BeadChip methylation arrays; however, it has lower resolution and sensitivity.³¹

In this study, we aimed to analyze the most frequent genetic and global DNA methylation alterations of WNT pathway genes in parallel during colorectal normal-adenoma-carcinoma sequence progression. We tested both promoter and gene body DNA methylation alterations using MethylCap-seq and studied the potential regulatory role of the identified DMRs on mRNA expression of WNT signaling pathway genes. Approximately 15-20% of CRCs belong to the distinct molecular subtype of cancers called CpG island methylator phenotype (CIMP) with high degree of DNA methylation of certain genes (CIMP-high).⁷ However, the DNA methylation changes also affect CIMP-negative (CIMP-zero or CIMP-low) CRCs, which represent the majority (approximately 75%) of sporadic CRCs.^7,17 In this study, we focused on CIMP-negative CRCs and their precancerous lesions in order to find common DNA methylation alterations of WNT pathway genes typical in development and progression of sporadic colorectal cancer.

Results

Determination of microsatellite instability (MSI) and CpG island methylation phenotype (CIMP) status

According to the results of immunohistochemical staining for MMR genes (MLH1, MSH2, MSH6, PMS2), all samples were found to be microsatellite stable (MSS). Using dichotomized CIMP status determination [CIMP-negative (CIMP-zero or CIMP-low) or CIMP-positive (CIMP-high)], all tumor samples were proven to be CIMP-negative, according to both classic^7,32 and Weisenberger CIMP status panel.³³

DNA methylator analysis of WNT signaling pathway genes

Methylation status of 160 WNT pathway-related genes (Gene symbols shown in Table 1 and detailed gene information shown in Supplementary Table 1) was analyzed at 100 bp resolution in 30 colorectal tissue samples [9 colorectal cancer (CRC), 15 adenoma (AD), 6 normal adjacent tissue (NAT)] by MethylCap-seq (Table 2 and Table 3).

Table 2. Clinical data of colorectal tissue samples used in MethylCap-seq study.

Sample ID	Patient ID	Sex	Age	Localization	Histology	TNM	Grade	Dukes' stage	Adenoma size
Normal adjacent tissue
MCS_20090930_0001_01_100	1	M	74	rectum	normal adjacent tissue (from adenoma patient 1)
MCS_20110000_0019_01_100	19	F	76	rectum	normal adjacent tissue (from CRC patient 19)	T3N0M0	G2	B2
MCS_20061114 _0020_01_100	20	M	72	coecum	normal adjacent tissue (from CRC patient 20)	T3N0M0	unknown	B2
MCS_20061113 _0021_01_100	21	M	68	sigma	normal adjacent tissue (from CRC patient 21)	T3N0M0	unknown	B2
MCS_20110000_0022_01_100	22	M	58	sigma	normal adjacent tissue (from CRC patient 22)	T3N0M0	G2	B2
MCS_20110000_0023_01_100	23	F	68	sigma and rectum	normal adjacent tissue (from CRC patient 23)	T4N0M0	G1	B3
Adenoma
MCS_20100222_0004_01_210	4	M	78	rectum	adenoma with low-grade dysplasia				25 mm
MCS_20100427_0005_01_210	5	M	65	rectum	adenoma with low-grade dysplasia				25 mm
MCS_20100609_0007_01_210	7	M	66	rectum	adenoma with low-grade dysplasia				20 mm
MCS_20110310_0012_01_210	12	M	83	sigma	adenoma with low-grade dysplasia				35mm
MCS_20110914_0013_01_210	13	M	55	colon desc, sigma	adenoma with low-grade dysplasia				10 mm
MCS_20120109_0014_01_211	14	M	81	sigma	adenoma with low-grade dysplasia				5 mm
MCS_20120308_0015_01_211	15	M	48	colon transv	adenoma with low-grade dysplasia				5 mm
MCS_20120320_0016_01_211	16	M	58	sigma	adenoma with low-grade dysplasia				15 mm
MCS_20120801_0017_01_211	17	F	56	colon asc, transv	adenoma with low-grade dysplasia				8 mm
MCS_20091104_0002_01_220	2	M	60	rectum-sigma	adenoma with high-grade dysplasia				15 mm
MCS_20100610_0006_01_220	6	M	68	colon desc	adenoma with high-grade dysplasia				5 mm
MCS_20100810_0008_01_220	8	M	73	sigma	adenoma with high-grade dysplasia				8 mm
MCS_20100923_0009_01_220	9	M	74	rectum	adenoma with high-grade dysplasia				8 mm
MCS_20110111_0010_01_220	10	M	68	colon desc	adenoma with high-grade dysplasia				15 mm
MCS_20110203_0011_01_220	11	F	76	colon desc	adenoma with high-grade dysplasia				20 mm
Colorectal cancer
MCS_20110203_0011_02_312	11	F	76	rectum	adenocarcinoma	T2N0M0	G2	B1
MCS_20110310_0012_02_312	12	M	83	coecum-colon asc	adenocarcinoma	T2N0M0	G2	B1
MCS_20061114 _0020_02_310	20	M	72	coecum	adenocarcinoma	T3N0M0	unknown	B2
MCS_20061113 _0021_02_310	21	M	68	sigma	adenocarcinoma	T3N0M0	unknown	B2
MCS_20110000_0023_02_311	23	F	68	sigma and rectum	adenocarcinoma	T4N0M0	G1	B3
MCS_20100201_0003_01_320	3	M	68	sigma-rectum	inoperablof CpG sitee CRC with multiplex metastases	unknown	unknown	D
MCS_20100610_0006_02_322	6	M	68	colon desc	adenocarcinoma	T3N2Mx	G2	D
MCS_20110111_0010_02_322	10	M	68	rectum	adenocarcinoma	T2N2M0	G2	C2
MCS_20110914_0013_02_322	13	M	55	rectum	adenocarcinoma	T3N1M0	G2	C2

Table 3. Samples used in the study.

		DNA METHYLATION ANALYSIS			GENE EXPRESSION ANALYSIS
		In silico analysis			In silico analysis
					Gene expression analysis					Demethylation study
		MethylCap-seq	Bisulfite sequencing validation^18	Methylation array analysis^28	GSE37364	GSE8671	GSE18105/GSE32323	GSE22242	GSE9348	GSE29060	GSE32323	GSE14526	GSE41588
TISSUE	NAT/N	6	15	24	38	32	17	4	12
	AD	15	15	32	29	32		3
	CRC	9	15	64	27		17	7	70
CELL LINE	HT29									+	+		+
	RKO										+
	SW480										+	+
	HCT116										+	+
	COLO320										+

Methylated CpG containing DNA fragments were effectively isolated using MethylCap protein (Supplementary Figure 1A). Recovery was 12.51 ± 10.14% in methyl-captured DNA samples, while the unmethylated fraction was slightly detected (0.41 ± 0.16%). Sequencing results of all samples fulfilled the quality criteria, according to the MEDIPS quality control protocol.³⁴ Pearson correlation coefficient of saturation analysis exceeded the 0.5 threshold for every sample (0.71 ± 0.07). Coverage of at least 5% of the methylation sites (11.43 ± 2.56%) was more than 5 times for each sample during the CpG coverage analysis. According to the CpG enrichment test, the efficacy of the methyl capture was proven to be appropriate as all samples showed value well above the 1.4 threshold (enrichment score relH: 2.87 ± 0.42) (Supplementary Figure 1B, C, D).

During the determination of DMRs both in promoter and in gene body regions, methylation changes were considered significant if P-values were under 0.05 and at least a 10% methylation difference could be measured (absolute value of Δβ was > 0.1).

The promoter regions of 33 out of the 160 analyzed WNT pathway genes were differentially methylated in CRC compared to NAT, including hypermethylated AXIN2, CHP1, PRICKLE1, SFRP1, SFRP2, SOX17, and hypomethylated CACYBP, CTNNB1, MYC (Fig. 1A, Table 1). In AD vs. NAT comparison, altered promoter methylation of 44 WNT signaling genes was detected, including hypermethylated APC, AXIN2, DAAM2, DKK4, PRICKLE1, SOX17, SFRP1, SFRP2 and SFRP4, and hypomethylated CACYBP, FZD3 (Fig. 1B, Table 1). Forty-one genes were identified showing different promoter methylation levels between CRC and AD samples. Hypermethylation of AXIN2, DKK1, VANGL1, and WNT5A gene promoters was increased in CRC compared to AD, while promoter methylation of SOX17, PRICKLE1, DAAM2, and MYC genes was higher in AD than in CRC (Fig. 1C). The whole gene list, including genomic positions and P- and β-values for DMRs, is shown in Supplementary Table 2A (CRC vs. NAT), 2B (AD vs. NAT), and 2C (CRC vs. AD).

Figure 1. Differentially methylated regions in WNT pathway gene promoters. Methylation probabilities (β-values) of 100 bp long analyzed regions were calculated with respect to genome wide CpG density dependent Poisson distributions, and are represented on a 0–1 scale. Promoter DMRs are shown A. in colorectal cancer (CRC) vs. normal adjacent tissue (NAT) comparison B. in adenoma (AD) vs. NAT comparison and C. between CRC vs. AD samples. Methylation intensities are illustrated on a color scale: high methylation levels are marked with red, low methylation levels are represented by green. Gene symbols mentioned in the text are listed on the heat maps, the whole gene lists with genomic positions, P- and β-values of DMRs can be seen in Supplementary Table 2.

Promoter methylation alterations have been firstly mentioned, since their effect on gene expression is well described, but most of the detected differentially methylated CpG sites were localized in gene body regions. In CRC samples, 93.75% (150 out of 160) of the analyzed WNT pathway genes had DMRs in their gene body region (all together, 3788 DMRs). In the adenoma group, a similar proportion (96.88%) of WNT pathway genes showed aberrant gene body methylation, with approximately 1.5-fold higher number of DMRs (total of 5737 DMRs).

In silico DNA methylation analysis of WNT pathway genes

CpG sites belonging to the significant DMRs in our MethylCap-sey study were searched using the Illumina BeadChip 450K array data sets from Luo et al.²⁸ (GEO accession number: GSE48684; Table 3) and P- and Δβ-values were determined in CRC vs. NAT, AD vs. NAT, and CRC vs. AD comparisons. Approximately one third of significant promoter DMRs detected by MethylCap-seq were represented in the Illumina BeadChip 450K data set by at least one cg ID (one CpG site; Table 4). Within the represented regions, 58% of the promoter methylation changes could be validated between CRC and NAT samples, 83% between AD and NAT samples, and 30% in CRC vs. AD comparison. In the case of gene body alterations, only approximately 3% of significant DMRs obtained by MethylCap-seq were represented in the BeadChip 450K array, and 25–51% of the methylation alterations could be confirmed between the diagnostic groups (Table 4).

Table 4. Comparison of BeadChip 450K²⁸ data and MethylCap-seq results for WNT pathway genes.

Comparison	Number of DMRs by MethylCap-seq	Number of DMRs represented in BeadChip 450K	Percentage (rate) of DMRs represented on BeadChip 450K	Number of correspondingly changed DMRs (represented on BeadChip 450K )	Percentage (rate) of correspondingly changed DMRs
PROMOTER
CRC vs. NAT	58	19	33% (19/58)	11	58% (11/19)
AD vs. NAT	129	47	36% (47/129)	39	83% (39/47)
CRC vs. AD	77	20	26% (20/77)	6	30% (6/20)
GENE BODY
CRC vs. NAT	3788	115	3.04% (115/3788)	30	26% (30/115)
AD vs. NAT	5737	167	2.91% (167/5737)	85	51% (85/167)
CRC vs. AD	4286	106	2.47% (106/4286)	26	25% (26/106)

Represented DMRs indicate the 100-bp-long differentially methylated regions detected by our MethylCap-seq study which are represented by at least one Illumina BeadChip 450K cg ID (CpG site). Correspondingly changed DMRs indicate differentially methylated regions represented on both analyses that were found to be concordant (similarly hyper- or hypo-methylated in the given comparison).

Gene expression analysis of WNT signaling pathway genes

Gene expression data of colonic tissue samples from the Gene Expression Omnibus (GEO) database were involved in in silico mRNA expression analysis of 160 WNT pathway genes (GEO accession numbers: GSE37364,³⁵ GSE8671,³⁶ GSE18105,³⁷ GSE32323,³⁸ GSE22242,³⁹ GSE9348)⁴⁰ (Table 3). The log Fold Change (logFC) data of significantly differentially expressed WNT pathway genes, together with the methylation data, are shown in Supplementary Table 2 for AD vs. NAT, CRC vs. NAT, and CRC vs. AD comparisons. Compared to promoter methylation alterations, the expression of 10 WNT signaling genes, including CACYBP, CTNNB1, MYC, PRICKLE1, PSEN1, and SFRP1, changed oppositely in CRC vs. NAT comparison (Table 5). The expression of 14 WNT pathway genes, including APC, CHP1, DAAM2, PRICKLE1, RUVBL1, SFRP1, SFRP2, SFRP4, and TLE3, were inversely correlated with promoter methylation status in AD vs. NAT samples (Table 6). In CRC vs. AD tissues, the expression of 8 genes, including CTNNB1, CXXC4, PRICKLE1, VANGL1, and SFRP2, was opposed to promoter methylation differences (Table 7). Inverse relation between promoter methylation and gene expression was found in CACYBP, CHP1, CTNNB1, CXXC4, PRICKLE1, SFRP1, SFRP2, and TLE3 genes, at least in 2 of the 3 comparisons (CRC vs. NAT, AD vs. NAT, and CRC vs. AD).

Table 5. WNT pathway genes showing inverse relation between DNA methylation and mRNA expression (CRC vs. NAT).

					Gene expression
Gene data			DNA methylation		LogFC
Gene symbol	Chr	Start position	P-value	Δβ (CRC vs. NAT)	Affymetrix probe set ID	GSE37364	GSE18105/ GSE32323	GSE22242	GSE9348
CACYBP	1	174967901	0.0289	−0.44	201381_x_at	0.41	0.52	0.34	0.72
					201382_at	−0.23
					210691_s_at	0.79	0.97	0.36	0.86
					211761_s_at	0.35	0.52	0.65	1.01
CHP1	15	41524201	0.0348	0.39	207993_s_at	−1.02	−1.67	−1.98	−2.36
					214665_s_at	−1.01	−1.44		−1.65
CTNNB1	3	41239901	0.0230	−0.38	1554411_at	0.94
		41240001	0.0080	−0.30	201533_at	0.54	−0.45		0.38
		41240101	0.0013	−0.45	223679_at	0.67
CSNK1E	22	38795901	0.0361	−0.27	202332_at	0.82	0.88		1.53
					222015_at	0.44			0.26
					225756_at	0.96	0.87		1.38
					226858_at	1.51	0.53		0.78
DKK4	8	42234701	0.0104	−0.30	206619_at	0.99
		42234801	0.0202	−0.27
FZD6	8	104311901	0.0402	−0.25	203987_at	0.85	2.02		1.07
MYC	8	128746501	0.0239	−0.21	202431_s_at	2.15	3.52	1.52	3.13
		128748901	0.0034	−0.25
PRICKLE1	12	42982701	0.0133	0.45	226065_at	−0.75
		42983001	0.0201	0.43	226069_at	−0.40			−0.33
					230708_at	−0.38
					232811_x_at	−0.13
PSEN1	14	73603901	0.0270	0.34	203460_s_at	−0.74	−0.76		−1.20
					207782_s_at	−0.74	−0.68	−1.06	−1.57
					226577_at	−1.03	−0.89		−1.55
					238816_at	−1.15			−2.07
SFRP1	8	41166001	0.0167	0.44	202035_s_at	−0.36	−1.31	−1.30	−0.14
		41166101	0.0001	0.56	202036_s_at	−2.15	−1.83	−2.41
		41166201	0.0034	0.49	202037_s_at	−1.80	−3.02	−1.21
		41166301	0.0138	0.27
		41166401	0.0013	0.53
		41166501	0.0271	0.44
		41166601	0.0345	0.39
		41166901	0.0420	0.41

Log FC values are indicated only in cases of P < 0.05.

Table 6. WNT pathway genes showing inverse relation between DNA methylation and mRNA expression (adenoma vs. NAT).

					Gene expression
Gene data			DNA methylation		LogFC
Gene symbol	Chr	Start position	P-value	Δβ (AD vs. NAT)	Affymetrix probe set ID	GSE37364	GSE8671	GSE22242
APC	5	112042901	0.0344	0.29	203525_s_at	−0.83	−1.06	−1.00
		112074301	0.0072	0.39	203526_s_at	−1.49	−1.29	−1.72
					203527_s_at	−0.63	−0.54
					215310_at	−0.23	−0.11
					216933_x_at	−0.12	−0.57	−0.91
CACYBP	1	174967301	0.0447	−0.23	201381_x_at	0.23	0.54
		174967501	0.0116	−0.25	210691_s_at	0.31	0.74
		174967801	0.0009	−0.35	211761_s_at		0.49
		174967901	0.0003	−0.50
		174968001	0.0034	−0.44
CHP1	15	41524401	0.0155	0.28	207993_s_at	−0.78		−1.65
					214665_s_at	−0.64	−0.46
CXXC4	4	105415501	0.0080	−0.26	220277_at	0.88	0.62
		105415601	0.0265	−0.24	229774_at	0.48	0.32
DAAM2	6	39760001	0.0024	0.58	212793_at	−1.95	−1.80	−1.28
		39760301	0.0281	0.41
FZD1	7	90893001	0.0154	0.34	204451_at	−1.24	−1.20	−0.74
					204452_s_at	−0.50	−0.29
NFAT5	16	69598901	0.0388	−0.24	215092_s_at			2.09
					224984_at	0.46
PRICKLE1*	12	42982701	0.0005	0.53	226065_at	−1.22	−0.59
		42982801	0.0337	0.36	226069_at	−0.55	−0.41	−0.99
		42983001	0.0046	0.46	230708_at	−0.56	−0.11
		42983101	0.0054	0.54
		42983201	0.0241	0.37
PRKCB	16	23847601	0.0015	0.59	207957_s_at	−1.48	−1.33	−1.24
		23847801	0.0292	0.41	209685_s_at	−1.52	−1.55	−0.86
					227817_at	−0.71	−0.39
					230437_s_at	0.01
RUVBL1	3	127843201	0.0334	−0.35	201614_s_at	1.35	1.75	0.79
SFRP1*	8	41165901	0.0023	0.59	202035_s_at	−0.29	−0.19	−1.30
		41166001	0.0016	0.59	202036_s_at	−2.03	−1.00	−2.35
		41166101	0.0004	0.56	202037_s_at	−1.77	−2.23	−1.52
		41166201	0.0001	0.60	228413_s_at
		41166301	0.0065	0.27
SFRP2*	4	154709401	0.0122	0.36	223121_s_at			−5.87
		154709501	0.0009	0.59	223122_s_at		−0.26	−5.86
		154709601	0.0046	0.55
		154709701	0.0001	0.63
		154709801	0.0012	0.57
SFRP4	7	37955801	0.0359	0.35	204051_s_at	0.01		−3.71
		37955901	0.0067	0.49	204052_s_at	0.02		−2.65
		37956101	0.0219	0.44
TLE3	15	70390701	0.0148	0,29	206472_s_at	−0.62
					212769_at	−0.56	−0.28
					212770_at	−0.59	−0.64
					228340_at	−0.71	−0,72

* + 4 (PRICKLE1); + 9 (SFRP1); +9 (SFRP2) hypermethylated DMRs see in Supplementary Table 3B. Log FC values are indicated only in cases of P<0.05.

Table 7. WNT pathway genes showing inverse relation between DNA methylation and mRNA expression (CRC vs. adenoma).

					Gene expression
Gene data			DNA methylation			LogFC
Gene symbol	Chr	Start position	P-value	Δβ (CRC vs. AD)	Affymterix probe set ID	GSE37364	GSE22242
AXIN2	17	63554401	0.0294	0.25	222696_at	−1.14	−1.99
					224176_s_at	−0.97	−1.04
					224498_x_at	−0.86
CTNNB1	3	41240001	0.0142	−0.20	1554411_at	0.97	1.42
		41240101	0.0071	−0.33	201533_at	0.37	1.30
		41240201	0.0330	−0.11	223679_at	0.47	0.93
CXXC4	4	105415501	0.0164	0.20	220277_at	−0.68	−0.93
		105416801	0.0278	0.21	229774_at	−0.47	−0.90
NFATC1	18	77159401	0.0240	−0.40	209664_x_at	−0.02
		77159501	0.0374	−0.35	211105_s_at		1.22
		77159601	0.0390	−0.37
PRICKLE1	12	42983201	0.0242	−0.31	226065_at	0.47	1.22
		42983301	0.0198	−0.27	226069_at		1.58
		42984301	0.0191	−0.29	230708_at	0.18
		42984401	0.0183	−0.25
		42984501	0.0358	−0.20
SFRP2	4	154709901	0.0429	−0.35	223121_s_at	1.39	6.09
		154710001	0.0484	−0.35	223122_s_at	2.12	6.70
		154710301	0.0213	−0.33
		154710801	0.0143	−0.34
		154711501	0.0267	−0.17
		154711601	0.0171	−0.17
TLE3	15	70390601	0.0341	−0.25	206472_s_at	1.02
		70390701	0.0356	−0.22	212769_at	1.04
					212770_at	0.37
VANGL1	1	116185801	0.0122	0.25	219330_at	−0.48	−0.46
					229134_at	−0.42	−0.31
					229492_at	−0.55	−0.37
					229997_at	−0.36	−0.38

Log FC values are indicated only in cases of P < 0.05.

Demethylation analysis on colon cancer cell lines

Gene expression data of 5 5-aza-2′-deoxycytidine-treated colon cancer cell lines were evaluated (GEO accession numbers: GSE29060,⁴¹ GSE32323,³⁸ GSE14526,⁴² and GSE41588⁴³)(Table 3).

Upregulation of mRNA expression was observed—at least in one of the analyzed colon cancer cell lines after 5-aza-2′-deoxycytidine demethylation treatment—for 7 genes (APC, CHP1, PRICKLE1, PRKCB, PSEN1, SFRP1, and TLE3) showing promoter hypermethylation and downregulated mRNA expression in benign and/or malignant colonic tissue samples. In SW480 and HCT116 colorectal adenocarcinoma cell lines, 86–100% of the genes showed opposite expression changes caused by the demethylation treatment (Fig. 2). GSE29060 and GSE41588 data sets were also evaluated, but mRNA expression of the above-mentioned genes was not restored in HT29 cell line by the demethylating treatment.

Figure 2. Effect of 5-aza-2′-deoxycytidine demethylation treatment on colon cancer cell lines. Gene expression data of colon cancer cell lines under 5-aza-2′-deoxycytidine demethylation treatment analyzed on Affymetrix HGU133 Plus 2.0 microarrays are shown (GEO accession numbers of gene expression data sets: GSE32323,³⁸ GSE14526⁴²). Demethylation agent treatment reversed (upregulated) mRNA expression of 7 WNT pathway genes showing promoter hypermethylation and downregulated mRNA expression in benign and/or malignant colonic tissue samples, principally in SW480 and HCT116 colorectal adenocarcinoma cell lines in varying degrees. A. GSE32323 data set B. GSE14526 data set. Samples are shown in columns, selected transcripts are represented in rows. High mRNA expression intensities are marked in red, low expression levels are shown in green. AZA = 5-aza-2′-deoxycytidine.

APC and β-catenin mutation analysis

Mutation analysis of APC and β-catenin genes was performed on the same sample set that was used in the whole methylome profiling. APC mutations were detected in 0% of NAT, 27% of adenoma, and 29% of CRC samples, while none of the samples showed β-catenin mutation. Two truncated protein resulting in frameshift mutations (COSM13113 and COSM296347) and 2 nonsense mutations with G/T substitution (COSM18775 and COSM18759) were detected. A SNP with G/A substitution (SNPrs41115) was found in 53% of adenoma samples and in 71% of CRC samples at the targeted APC sequence region. The list and details of the detected mutations can be seen in Table 8.

Table 8. Detected mutations in analyzed colonic tissue samples.

Gene symbol	Variant	cDNA position	Amino acid alteration	Affected cDNA sequence	Annotation	Affected samples
APC_1	96:C/T	3871	Q1291stop	ATGTAATCAGACG	COSM19072	—
APC_1	140:A/-	3915	Frameshift from aa 1306 Truncated protein (length 1366 aa instead of 2843 aa)	ATAGCAGAAATA	—	—
APC_1	141:G/T	3916	E1306stop	ATAGCAGAAATA	COSM18760	—
APC_2	31:G/T	3925	E1309stop	ATAAAAGAAAAGA	COSM18775	MCS_20100810_0008_01_220 (AD)
APC_2	33–37:AAAGA/—–	3927	Frameshift from aa 1309 Truncated protein (length 1311 aa instead of 2843 aa)	ATAAAAGAAAAGATTGG	COSM13113	MCS_20100923_0009_01_220 (AD) MCS_20110203_0011_02_312 (CRC) MCS_20110000_0023_02_311 (CRC)
APC_2	70:G/T	3964	E1322stop	TGTGAGCGAAGTTCC	COSM18702	—
APC_2	139:G/T	4033	E1345stop	ATCTTCAGAATCA	COSM18759	MCS_20100427_0005_01_210 (AD)
APC_3	125:G/T	4135	E1379stop	GTTCAGGAGACCCC	COSM18834	—
APC_4	111.5:-/T	4241	Frameshift from aa 1415 Truncated protein (length 1421 aa instead of 2843 aa)	GAATGGTAAGTGGCA	COSM19703	—
APC_5	67:A/-	4314	Frameshift from aa 1439 Truncated protein (length 1421 aa instead of 2843 aa)	TAAAACACCTCCAC	COSM296347	MCS_20110111_0010_01_220 (AD)
APC_5	101:C/T	4348	R1450stop	ACCAAGCGAGAAGTA	COSM13127	—
APC_6	105:G/A	4479	Silent T1493T	GCCACGGAAAGTA	SNP rs41115	MCS_20091104_0002_01_220 (AD) MCS_20100222_0004_01_210 (AD) MCS_20100610_0006_01_220 (AD) MCS_20100810_0008_01_220 (AD) MCS_20100923_0009_01_220 (AD) MCS_20110310_0012_01_210 (AD) MCS_20120109_0014_01_211 (AD) MCS_20120320_0016_01_211 (AD) MCS_20110203_0011_02_312 (CRC) MCS_20100610_0006_02_322 (CRC) MCS_20110111_0010_02_322 (CRC) MCS_20110310_0012_02_312 (CRC) MCS_20110000_0023_02_311 (CRC)
CTNNB1_1	123:C/T	47	P16L	TGGAACCAGACAGA	—	—
CTNNB1_2	87:C/T		S45F	TCCTTCTCTGAGT	COSM5667	—

aa: aminoacid

Discussion

Similarly to other malignancies, colorectal cancer arises as a consequence of accumulating genetic and epigenetic alterations. Beside the conventional mutations, altered DNA methylation has also been revealed as a causative factor in carcinogenesis and tumor progression. Gene silencing mediated by aberrant promoter DNA methylation is one of the key features of carcinogenesis.⁸ Gene body regions are also strongly methylated, which can positively correlate with the gene expression.^8,9 The importance of DNA methylation is further enhanced by the recently initiated trials using low-dose DNA demethylating agents as anticancer drugs.⁹

We focused on DNA methylation analysis of WNT signaling genes/gene promoters in CIMP-negative CRCs and precancerous adenomas. Although since the establishment of CIMP conception the main focus has been on DNA methylation analysis of CIMP-positive proximal colon cancers, notable epigenetic alterations, including promoter DNA hypermethylation, can be detected also in CIMP-negative CRCs representing the majority of sporadic CRCs.^4,44 A reliable plasma methylation marker test with high specificity and sensitivity—such as Epi proColon 2.0 analyzing SEPT9 promoter methylation⁴⁵—should be based on systematic DNA methylation alterations identified and detectable both in CIMP-negative and CIMP-positive CRCs.

In this study, systematic analysis of DNA methylation of WNT signaling pathway genes in relation to CRC development and progression was performed using methyl capture sequencing (MethylCap-seq) and whole genomic expression microarrays. MethylCap-seq can test both promoter and gene body regions in the entire genome. For comparing previous and our DNA methylation results, the in silico analysis of Luo et al.²⁸ data, also including adenomas besides CRC and NAT samples, was performed. Similarly to the work of Luo et al.,²⁸ we compared a group of NAT samples to CRC and adenoma sample groups, instead of doing a pairwise comparison of CRC and NAT samples.

When compared to the methylation array, MethylCap-seq provides a more extensive overview of whole genome methylation, although its resolution and sensitivity is lower than that of BeadChip methylation arrays.³¹ Concerning WNT signaling genes in our experiments, approximately one third of the significant DMRs was represented and could be analyzed on BeadChip array for gene promoters; only ∼3% of the identified DMRs for gene body regions was represented. The coexisting compatible regions in both test systems showed high DNA methylation similarity between the diagnostic groups, especially between adenoma and NAT samples (83%). Concordance between CRC and NAT groups was 58%, but only 30% between CRC and AD. This partial discordance was most probably related to the difference between our patient cohorts concerning CIMP and MSI or mutation status. Regarding resolution, the Illumina BeadChip 450K array detects the methylation status of a given CpG site, while MethylCap-seq measures the summation of methylated CpG sites located within the given 100 bp long analysis window.

In our results, promoter methylation changes affected both canonical (WNT/β-catenin) and noncanonical [planar cell polarity (PCP) and WNT/Ca²⁺] WNT signaling pathways. In CRC samples, significant DMRs were found in 20.63% of the WNT pathway gene promoters, which was even higher in adenoma samples (27.50%). Aberrant promoter methylation was mainly detected in the canonical WNT (22.3%) and WNT/Ca²⁺ pathway genes (22.9%) in CRC, but the PCP pathway was also affected (18.8%). In AD samples, promoter methylation changes were principally observed in WNT/Ca²⁺ pathway genes (40%), while high number of DMRs were also found in canonical and PCP WNT pathways genes (28.1%).

Similarly to our results, promoter hypermethylation of canonical WNT inhibitory genes, such as SFRP1,^11,12 SFRP2,^11,13-15 and SOX17,^11,24 was detected both in CRC and adenoma samples. However, in contrast to the findings of Silva et al.,¹¹ in our study these genes showed elevated promoter methylation in adenoma samples compared to CRC tissue. In accordance with our previous bisulfite sequencing findings,¹⁸ SFRP1 and SFRP2 promoter methylation in CRC and AD tissue samples was also detected by MethylCap-seq at the same regions (Table 9).

Table 9. Methyl capture and bisulfite sequencing¹⁸ data for SFRP1 and SFRP2 promoter regions (overlapping regions).

	METHYL CAPTURE SEQUENCING							BISULFITE SEQUENCING^18
		P-value			Δβ value					DNA methylation difference
Gene symbol	Genomic localization (start) (hg19)	CRC vs. NAT	AD vs. NAT	CRC vs. AD	CRC vs. NAT	AD vs. NAT	CRC vs. AD	Genomic localization (hg19)	Amplicon length (bp)	CRC vs. NAT	AD vs. NAT
SFRP1	WNT_SFRP1_41166101	0.0001	0.0004	—	0.56	0.56	–	chr8:41166121-41166356	235	0.21	0.15
	WNT_SFRP1_41166201	0.0034	0.0001	—	0.49	0.60	–
	WNT_SFRP1_41166301	0.0138	0.0065	—	0.27	0.27	–
SFRP2	WNT_SFRP2_154709501	0.0203	0.0009	—	0.46	0.59	–	chr4:154709513-154709730	217	0.19^*	0.28^*
	WNT_SFRP2_154709601	—	0.0046	—	—	0.55	–
	WNT_SFRP2_154709701	0.0218	0.0001	—	0.47	0.63	–
	WNT_SFRP2_154710201	0.0254	0.0022	—	0.30	0.33	–	chr4:154710290-154710649	359	0.40^*	0.46^*
	WNT_SFRP2_154710301	—	0.0022	0.0213	—	0.58	−0.33

MethylCap-seq analysis of SFRP1 and SFRP2 promoters revealed that all 100 bp-long sequences overlapping with the regions analyzed by bisulfite sequencing showed elevated DNA methylation level in CRC and adenoma samples compared to NAT samples (P < 0.05). In case of MethylCap-seq, only the significant (P < 0.05) methylation changes are represented. The 100-bp analysis windows completely overlapping with bisulfite-sequenced regions are marked in bold.

* Significant methylation differences in bisulfite sequencing analysis¹⁸

Correspondingly to the previous findings,^{11,17, 19-20} we found that another known inhibitor of canonical WNT pathway, DKK2, was significantly hypermethylated in CRC vs. NAT and AD vs. NAT comparisons. However, our promoter estimation considering the sequence in the region between 2000 bp upstream and 1000 bp downstream from the transcription start site (if the positive DNA strand was the coding strand) identified ‘active promoter’ in at least one of the 9 analyzed cell lines, according to the ENCODE ChromHMM results.⁴⁶ The studied DMRs in DKK2 promoter were excluded because they were labeled as ‘poised promoters’ by ENCODE, despite their significance (P-value was between 0.00001 and 0.035 in CRC vs. NAT; P < 0.01 in AD vs. NAT) and high Δβ-values (0.41–0.72 in CRC vs. NAT and 0.49–0.60 in AD vs. NAT). The lack of significant changes in mRNA expression data between groups also supports the tested DKK2 region as an inactive promoter.

Our MethylCap-seq results indicate that promoter DMRs of APC, AXIN1, AXIN2, and CSNK1 genes may participate in the impairment/reduction of the APC-Axin-CSNK1-β-catenin complex leading to the constitutive activation of the canonical WNT pathway. APC hypermethylation has been already described in adenoma and CRC tissue and CRC plasma samples.^19,47,48 In accordance with findings of Judson et al.,⁴⁸ hypermethylation in the active promoter of the APC gene in adenoma samples was detected in our study as well (Table 6). Similar alteration appeared in CRC samples, but in a different region upstream to the APC gene, which was considered as a ‘weak promoter’ or ‘enhancer’. Reduced APC mRNA levels on microarrays indicate silencing by promoter hypermethylation, which impairs APC tumor suppressor function during colorectal carcinogenesis, besides the frequent APC mutations. DNA methylation differences in APC promoter and mutation hot spot regions could not been detected between samples with or without APC mutations. This supports that APC mutations and APC promoter hypermethylation are not exclusive factors, but rather 2 parallel mechanisms leading to the loss of APC function. We found APC mutations in approximately one third of the CRC samples, which is below the mutation rate described by others (75–80%²²). However, we only analyzed a limited region [BLAST: chr5:112,175,118–112,176,026 (Feb. 2009 GRCh37/hg19 Assembly)] with the most prevalent APC mutations. Still, our results were similar to those of Luo et al.²⁸ finding APC mutations in 29.7% of CRC and 27.8% of adenoma tissue samples using larger sample size and testing longer APC regions.

Increased promoter methylation of the AXIN2 gene is thought to be associated with carcinogenesis of MSI+ CRC.⁴⁹ We also detected AXIN2 promoter hypermethylation with MethylCap-seq both in CRC and adenoma tissues compared to NAT samples, but the prevalence of DMRs was higher in gene body regions (25 hypermethylated DMRs, 13 hypomethylated DMRs). In silico analysis also revealed a remarkable upregulation of AXIN2 mRNA in both disease groups, suggesting the tumor promoting function of AXIN2,⁴⁹ especially in MSS CRC cases.

We also found β-catenin hypomethylation in CRC compared to NAT tissues, which was also detected earlier by others using Illumina methylation array and in methylation-sensitive high resolution melting analysis.²⁰ The significantly hypomethylated CpG site (41240163, Illumina ID: cg09678212) in the previous analysis²⁰ is located within WNT_CTNNB1_41240101 100 bp long analysis window in our study. In silico evaluation of Luo et al.²⁸ Illumina methylation array data also resulted in significant (P<0.001) hypomethylation of this CpG site in CRC compared to NAT (Δβ = -0.20) and adenoma samples (Δβ = -0.12). Besides identifying the only β-catenin promoter hypomethylated CpG site found by the Illumina platform, our MethylCap-seq analysis revealed significant hypomethylation in the 300 bp region (WNT_CTNNB1_41239901; WNT_CTNNB1_41240001 and WNT_CTNNB1_41240101) of the CTNNB1 gene promoter. Hypomethylation of this region in CRC compared to NAT and adenoma samples suggests its contribution to aberrantly elevated intracellular β-catenin level in colorectal carcinogenesis. Supporting this, significant elevation of CTNNB1 mRNA levels were observed in CRC compared to normal samples with logFC: 0.38–0.94, along with increased nuclear translocation of β-catenin protein as a sign of the activation of the canonical WNT pathway during malignant transformation. Oncogenic function of canonical WNT pathway activation is known to occur at early stages of carcinogenesis, and can be promoted by DNA methylation alterations from the precancerous adenoma stage.

The PCP pathway, as a noncanonical WNT signaling, can play multiple roles in oncogenesis.⁵⁰ At an early stage, PCP pathway has a tumor suppressive effect as it can antagonize the canonical pathway.⁵⁰ However, during tumor progression, it can be tumorigenic by promoting tumor cell motility and invasion, angiogenesis induction, and metastatic spreading.⁵⁰ DNA methylation changes both in promoters and transcribed regions of PCP pathway genes were detected at the adenoma stage. DMRs were observed in FZD2, FZD7, PRICKLE1, and WNT5A gene promoters, both in adenoma and CRC samples, while promoter hypermethylation on DAAM2 and downregulation of its mRNA level were typical only in adenoma tissue compared to NAT samples. Increased FZD2 and WNT5A levels were measured in cancer cell lines, including CRC and advanced tumors, and were shown to promote epithelial-mesenchymal transition and metastasis.^50,51 FZD7 can enhance migration and invasion of CRC mediated by noncanonical WNT signaling.⁵⁰ PRICKLE1, a disheveled (Dvl) associated protein, may have a tumor suppressor function by antagonizing Dvl recruitment by Frizzled.⁵² It can bind to Dvl3 and facilitate its ubiquitination-mediated degradation, suggesting that PRICKLE1 is a negative regulator of WNT/β-catenin pathway.^52,53 Promoter hypermethylation of PRICKLE1, found in adenoma and CRC samples in our methylation study, also support this idea. DAAM2 and PRICKLE1 seems to be DNA methylation regulated genes as inverse correlation was found between their promoter hypermethylation and reduced mRNA expression.

The noncanonical WNT/Ca²⁺ pathway, which can antagonize β-catenin-dependent signaling and stimulate cell migration,⁵³ was the most affected by promoter methylation changes in adenomas. This involved the WNT5A (ligand), FZD1, 2 and 7 (receptors), and downstream elements, such as PLCB1, PRKCB, and NFAT5 genes in AD samples, whereas promoter hypermethylation of FZD2, FZD7, PLCB1, and PRKCB genes was in silico-validated in adenoma samples.²⁸ Hypermethylation of FZD2, FZD7, PLCB1, and WNT5A gene promoters was also detected in CRC, compared to NAT samples. However, in many cases it was not accompanied by downregulation of mRNA expression. WNT5A has a dual function: it can act as tumor suppressor or as oncogene in a context-dependent manner. Depending on the type and level of co-receptors, it can activate noncanonical or β-catenin-dependent WNT signaling. Although its promoter was found to be hypermethylated both in adenoma and CRC tissues, with only moderate difference in the gene body region, it was probably not a DNA methylation-regulated gene. The remarkable increase on its mRNA expression (AD vs. NAT: logFC: 0.65–1.02; CRC vs. NAT logFC: 1.21–2.58) may be regulated by other mechanisms. FZD2 and FZD7 genes, which were upregulated in CRC (FZD7 is already in adenomas) may promote migration and metastatic spreading.⁵¹ The connection between methylation and mRNA expression of these genes was controversial, since both increased mRNA expression and promoter hypermethylation were measured in CRC, which was further elevated in adenoma tissue samples.

An inverse relation was observed between gene expression and promoter DNA methylation in WNT-related 23 genes. Among others, APC, CHP1, PRICKLE1, PRKCB, PSEN1, and SFRP1 genes were identified as methylation-regulated WNT pathway genes during colorectal normal-adenoma-carcinoma sequence progression in view of 5-aza-2′-deoxycytidine (5-aza) demethylation treatment data of colon cancer cell lines, specifically of HCT116 and SW480. SW480 is a MSS, CIMP-negative colorectal adenocarcinoma cell line⁵⁴ most similar to the clinical samples analyzed in our MethylCap-seq study. Demethylation treatment was capable of inducing re-expression of these genes in CRC cell lines, suggesting a casual relationship between promoter hypermethylation and expression silencing of these genes. Restoration of APC, PRKCB, and SFRP1 expression by 5-aza treatment was also observed earlier in different colon cancer cell lines.^55-57

Despite having detected most of the DMRs in gene body regions, we discuss these only shortly here because their contribution is not well known or controversial. Approximately 95% of the analyzed WNT pathway genes showed aberrant gene body methylation in adenoma and CRC compared to the NAT samples with ∼1.5-fold higher number of DMRs in adenomas. Gene body methylation may positively regulate gene expression^8,9 and can be involved in splicing regulation.^8,58 However, it can also contribute to gene silencing, according to the whole genomic bisulfite and RNA sequencing data of Lou et al.⁵⁹ In the recent study of Wang et al. about the roles of genic methylation in prostate tumorigenesis, 12 groups of genes with collaborative differential methylation patterns, including hyper-methylated promoter/hypermethylated gene body and hypermethylated promoter/hypomethylated gene body, were identified, assuming a gene activating role for gene body methylation.⁶⁰ In line with these findings, we identified strongly hypermethylated gene body associated with upregulated mRNA expression in WNT pathway genes, such as AXIN2, CSNK1E, MYC, NKD1, which expression was downregulated after 5-aza-2′-deoxycytidine treatment in colon cancer cell lines. Also, in other genes, reduced expression was associated with intense gene body hypermethylation reversed by 5-aza demethylation agent (e.g., PPARD) in CRC cell lines. Our study resulted in gene body methylation alterations between CRC and NAT samples in LRP5, TCF7L1, WNT2, and WNT6 WNT signaling pathway genes, on the same, mainly intronic regions, as in a previous study by Farkas et al.²⁰

Further research focusing on disease-specific methylation alterations and gene expression regulatory mechanisms is required to clarify the selective functions of DNA methylation at the gene body and other non-promoter regions.⁵⁹

Conclusions

DNA methylation alterations affect both canonical and noncanonical WNT signaling pathway genes at their promoter and gene body regions. Using MethylCap-seq, which allows extensive discovery of DNA methylation in the entire genome, we identified new aberrant methylation patterns affecting genes including CTNNB1, DAAM2, and PRICKLE1, characterizing the colorectal normal-adenoma-dysplasia sequence progression. Previously described methylation changes in WNT pathway genes during colorectal carcinogenesis, such as SFRP1, SFRP2, SOX17, and APC, could also be successfully verified.

Our results suggest that DNA methylation changes of canonical and noncanonical WNT signaling pathways also contribute to the development and progression of CIMP-negative CRC, besides the frequent APC mutations. Aberrant DNA methylation already appears in adenomas, indicating DNA methylation as an early event of colorectal carcinogenesis. The genomic regions showing epigenetic changes identified in this study might provide the basis for further studies revealing potential early diagnostic and therapeutic targets.

Patients and methods

Clinical samples

A total of 30 colonic tissue samples [15 adenoma, 9 CRC, 6 normal adjacent tissue (NAT)] were included in the MethylCap-seq study (Table 2, Table 3). All adenoma and part of CRC samples were collected during screening colonoscopy; further surgically removed CRC and NAT tissue samples were also analyzed.

After informed consent of untreated patients, colonic biopsy samples were taken during endoscopic intervention and stored in RNALater Stabilization Solution (Ambion, ThermoFisher Scientific, AM7024) at −80°C until use. Biopsy samples from the same site were immediately fixed in buffered formalin for histological evaluation. Histological diagnoses were established by experienced pathologists. Tissue samples from untreated CRC patients were also obtained from surgically removed colon or rectum tumors and from histologically NAT that originated from the area farthest available from the tumor. Tissue samples were snap frozen in liquid nitrogen directly after surgery and were stored at -80°C. Detailed patient specification is described in Table 2. The study was conducted according to the Helsinki declaration and approved by the local ethics committee and government authorities [Regional and Institutional Committee of Science and Research Ethics (TUKEB) Nr.: 69/2008, 202/2009 and 23970/2011 Semmelweis University, Budapest, Hungary].

MSI status was determined according to the results of immunohistochemical staining as described earlier¹⁷ with monoclonal antibodies for MLH1 (1:80), MSH2 (1:80), MSH6 (1:80), PMS2 (1:150, all from Novocastra/Leica Biosystems NCL-L-MLH1, NCL-MSH2, NCL-L-MSH6, NCL-L-PMS2) markers. A tumor sample was scored MSS if all of the 4 above MMR proteins showed expression.

CIMP category was specified using both classic (CDKN2A, MLH1, MINT1, MINT2, and MINT31)^7,32 and Weisenberger (CACNA1G, IGF2, NEUROG1, RUNX3, SOCS1) marker panels,³³ according to the MethylCap-seq data overlapping with regions of CIMP panel methylation-specific PCRs.^7,32,33 A sample was considered as CIMP-positive (CIMP-high) when at least 2 of the 5 markers (in case of the classic panel) and at least 3 of the 5 markers (in case of the new panel) were found to be methylated.

DNA isolation

Tissue samples were homogenized in 2% sodium dodecyl sulfate, and digested with 4 mg/mL proteinase K for 16 hat 56°C. Genomic DNA was isolated using High Pure PCR Template Preparation Kit (Roche Applied Science, 11796828001) according to the manufacturer's instructions.¹⁷ DNA was eluted in 2x with 100 μl RNase- and DNase-free water and stored at -20°C. Isolated DNA samples was quantified by Qubit fluorometer using Qubit dsDNA HS Assay (Invitrogen, ThermoFisher Scientific, Q32854).

Methyl capture sequencing

After fragmentation of 3 μg genomic DNA samples (100–400 bp) using Bioruptor Sonicator (Diagenode, UCD-300), the DNA fragments containing methylated CpGs were isolated in SX-8G IP-Star® Compact Automated System using Auto MethylCap kit (Diagenode, C02020011) using High Elution Buffer according to the manufacturer's instructions.⁶¹

The methylated DNA fraction was purified on QIAquick PCR purification columns (Qiagen, 28104). SYBR Green quantitative RT-PCR was performed using Human meDNA primer pair (TSH2B) (Diagenode, pp-1041-500) for detection of methylated DNA and Human unDNA primer pair (GAPDH) (Diagenode, pp-1044-500) for showing unmethylated DNA. Library preparation was carried out using Illumina TruSeq ChIP Sample Preparation kit (Illumina, IP-202). Cluster generation was done on cBot instrument using TruSeq SR Cluster Kit v3-cBot-HS (Illumina, GD-401-3001) which provides reagents that bind samples to complementary adapter oligos for cluster amplification on single-read flow cells. Next generation sequencing of the methylated DNA fragments was performed on HiScanSQ instrument using TruSeq SBS v3-HS reagents (Illumina, FC-401-3001) according to the manufacturer's instructions. Bowtie2⁶² with default settings was used to map the 100 bp paired and 50 bp unpaired reads to the hg19 human genome reference assembly.⁶³ The generated bam files were sorted and indexed by samtools.⁶⁴ Data were processed by the MEDIPS⁶⁵ bioconductor R package. After quality control, unpaired reads were extended to the length of the average fragment length of the paired samples (250 bp). PCR duplicates were removed, then the coverage data were binned with 100 bp window size. Methylation probabilities (β-values hereafter) were calculated with respect to genome wide CpG density dependent Poisson distributions.

In silico DNA methylation analysis of WNT pathway genes

Illumina BeadChip 450K DNA methylation array data sets from Luo et al.²⁸ were downloaded from NCBI Gene Expression Database and analyzed (GEO accession number: GSE48684). This data set contains the methylation results of 120 colorectal tissue samples including CRC, AD, and NAT samples (Table 3). CpG sites belonging to the significant DMRs in our MethylCap-sey study were searched on the Illumina BeadChip 450K array, and P- and Δβ-values were determined in CRC vs. NAT, AD vs. NAT and CRC vs. AD comparisons.

Gene expression analysis

Gene expression data of colonic tissue samples from the Gene Expression Omnibus (GEO) database were included in the in silico mRNA expression analysis of WNT pathway genes. WNT pathway genes were selected according to the KEGG pathway database. GSE37364 set contains Affymetrix HGU133 Plus2.0 whole transcriptome data of 94 colonic biopsy samples (38 healthy normal, 27 colorectal adenocarcinoma, and 29 tubulovillous/villous adenoma) previously hybridized by our research group.³⁵ Five GEO data sets made by others on the same platform were also analyzed (GSE8671,³⁶ GSE18105,³⁷ GSE32323,³⁸ GSE22242,³⁹ GSE9348⁴⁰). Gene expression data of CRC and healthy normal (N)/ normal adjacent tissue (NAT) samples are available under serial accession numbers GSE18105, GSE32323 and GSE9348. GSE8671 data set contains mRNA expression data of adenoma and healthy normal biopsy samples, GSE22242 dataset includes transcriptome analysis results of pooled normal, adenoma, and CRC biopsy samples. Differentially expressed genes were determined.

Demethylation analysis on colorectal cancer cell lines

Gene expression data of 5-aza-2′-deoxycytidine (5-aza)-treated colorectal cancer cell lines (HT29, COLO320, HCT116, RKO, SW480) were downloaded and in silico analyzed [GEO accession numbers: GSE29060 (previously performed by our research group),⁴¹ GSE32323,³⁸ GSE14526,⁴² and GSE41588⁴³]. GSE29060 and GSE41588 data sets contain the expression data of HT29 cells treated with 5-aza at higher concentration (5–10 μM) for 72 hours and 5 days respectively. In another 2 studies (GSE32323 and GSE14526) 5-aza was applied for 72 hours at lower concentrations (0.5 μM and 3 μM, respectively). Differentially expressed genes after 5-aza demethylation treatment were determined.

Mutation analysis

Mutation hot spot regions of APC and β-catenin genes were amplified using custom made PCR primers (shown in Supplementary Table 3) and amplicons were sequenced with a GS Junior instrument (Roche), as described earlier.⁶⁶ Briefly, Rapid Library Molecular Identifier (RL_MID) adaptors were ligated to the PCR products during library preparation providing a sample specific sequence motif. The quality checking of PCR libraries were made on Agilent Bioanalyzer instrument using High Sensitivity DNA Chip (Agilent, G2939AA, 5067-4626). Emulsion PCR amplification of the amplicon libraries was performed using the Lib-L emPCR Kit (Roche, 05996481001), with 2 DNA molecules per bead ratio following the manufacturer's instructions. Bead enrichment and sequencing were carried out using GS Junior Titanium Sequencing Kit (Roche, 05996554001) according to the method described in the Sequencing Method Manual, GS FLX Titanium Series. Amplicon Variant Analyzer software (Roche) was applied for identification of variants.

Statistical analysis

During statistical evaluation of gene expression and DNA methylation data of 160 WNT signaling pathway genes, false discovery rate (FDR) was applied for Student t-test with the criteria P < 0.05 in all paired comparisons. In case of gene expression for LogFC calculation, the differences between the averages of samples groups were considered (abs ≥ 1 criteria). Methylation alterations between diagnostic groups were characterized by Δβ-values (the differences of the average β-values of sample groups).

Disclosure of potential conflicts of interest

No potential conflicts of interest were disclosed.

Funding

This study was supported by the National Research, Development and Innovation Office (KMR-12-1-2012-0216 grant) and the validation study by the Hungarian Scientific Research Fund (OTKA-K111743 grant).