Genome methylation patterns in male breast cancer - Identification of an epitype with hypermethylation of polycomb target genes.

Male breast cancer (MBC) is a rare disease that shares both similarities and differences with female breast cancer (FBC). The aim of this study was to assess genome-wide DNA methylation profiles in MBC and compare them with the previously identified transcriptional subgroups of MBC, luminal M1 and M2, as well as the intrinsic subtypes of FBC. Illumina's 450K Infinium arrays were applied to 47 MBC and 188 FBC tumors. Unsupervised clustering of the most variable CpGs among MBC tumors revealed two stable epitypes, designated ME1 and ME2. The methylation patterns differed significantly between the groups and were closely associated with the transcriptional subgroups luminal M1 and M2. Tumors in the ME1 group were more proliferative and aggressive than ME2 tumors, and showed a tendency toward inferior survival. ME1 tumors also displayed hypermethylation of PRC2 target genes and high expression of EZH2, one of the core components of PRC2. Upon combined analysis of MBC and FBC tumors, ME1 MBCs clustered among luminal B FBC tumors and ME2 MBCs clustered within the predominantly luminal A FBC cluster. The majority of the MBC tumors remained grouped together within the clusters rather than being interspersed among the FBC tumors. Differences in the genomic location of methylated CpGs, as well as in the regulation of central canonical pathways may explain the separation between MBC and FBC tumors in the respective clusters. These findings further suggest that MBC is not readily defined using conventional criteria applied to FBC.

Introduction

Male breast cancer (MBC) is a rare disease, accounting for <1% of all breast cancer cases (Cancerfonden, 2013; Siegel et al., 2013). The majority of the tumors are estrogen receptor (ER) and progesterone receptor (PR) positive and the most common histologic type is invasive ductal carcinoma (Ruddy and Winer, 2013). MBC is similar to female breast cancer (FBC) in many ways, e.g. histologically, but there are also many distinct differences, including the distribution of age and grade (Nilsson et al., 2011; Giordano, 2005). Further, on a molecular level, it has been shown that MBC, like FBC, is a heterogeneous disease and differences between FBC and MBC on transcriptional (Johansson et al., 2012; Callari et al., 2010), copy number (Johansson et al., 2011; Tommasi et al., 2010) and microRNA levels (Fassan et al., 2009; Lehmann et al., 2010) have been reported. In a recent gene expression profiling study of MBC we described two new subgroups, luminal M1 and M2, that did not resemble any of the intrinsic subtypes reported in FBC, and as such may be specific to breast cancer in men. Luminal M1 tumors (70% of the MBC tumors) seemed to be more aggressive and were associated with worse prognosis and also appeared to have a less activated ER pathway, while luminal M2 tumors (30% of the MBC tumors) displayed an up‐regulated immune response and a more activated ER pathway (Johansson et al., 2012).

Epigenetic changes like promoter methylation and histone modifications play crucial roles in tumor progression (Hansen et al., 2011; Jaenisch and Bird, 2003). Epigenetic changes are stable and heritable during cell division even though they do not involve DNA mutations. As such, epigenetic alterations constitute an additional layer of regulation of the genome (Jaenisch and Bird, 2003). Histone modification by polycomb group (PcG) proteins can mediate gene silencing and maintenance of cellular identity. This can be carried out by the polycomb repressive complex 2 (PRC2) that catalyzes trimethylation of lysine 27 on histone H3 (H3K27me3), a mark of transcriptionally silent chromatin (Kondo et al., 2008). The core components of the PRC2 complex include EZH2, SUZ12, and EED (Kuzmichev et al., 2004) and target genes for the PRC2 complex in embryonic stem cells are involved in cell fate decisions, development and differentiation (Bracken et al., 2006). Poised promoters are associated with polycomb repressed regions and are marked with both H3K27me3 and trimethylation of lysine 4 on histone H3 (H3K4me3), thereby forming a bivalent domain which is associated with silencing of developmental genes in stem cells, while at the same time keeping them poised for activation (Bernstein et al., 2006; Ernst et al., 2011). It has been shown that genes that are marked with H3K27me3 in stem cells exhibit increased frequency of de novo methylation upon differentiation (Mohn et al., 2008). Moreover, methylation may influence the expression of cancer related genes and is hence a potential druggable target (Jovanovic et al., 2010). CpG islands, genomic regions that contain high frequencies of CpG sites, often become hypermethylated during tumorigenesis which may result in silencing of important tumor suppressor genes, while CpGs in other contexts tend to be hypomethylated genome‐wide, thereby potentially affecting chromosomal stability (Jones, 2002). DNA methylation studies have previously been limited to a small number of CpG islands in individual genes, but new microarray and sequencing techniques have recently made it possible to study global genome methylation patterns in CpGs at many different locations across the genome. Epigenetic changes have not previously been studied on a global level in MBC. The aim of the present study was to assess and interrogate genome‐wide methylation patterns in MBC and compare them with the transcriptionally derived subgroups of MBC, luminal M1 and M2, as well as the intrinsic subgroups of FBC.

Materials and methods

Patient and tumor material

All cases of MBC diagnosed between 1983 and 2008 at three Swedish hospitals (Skåne University Hospital, Uppsala University Hospital, and Örebro Hospital), with sufficient tumor material available were collected. DNA was extracted as previously described and applied to an array comparative genomic hybridization study (Johansson et al., 2011). The 47 tumors with sufficient remaining DNA were used in the present study. Patient and tumor characteristics are summarized in Table 1 (see (Nilsson et al., 2013) for further details).

Table 1

Clinico‐pathological characteristics in the whole MBC cohort, and in the ME1 and ME2 epitypes, respectively.

	All MBCs (N = 47)		ME1 (N = 23)		ME2 (N = 24)
	N	%	N	%	N	%
Age at diagnosis
Mean	68		70		67
Range	42–92		42–91		49–92
ER status
Positive	38	81	21	91	17	71
Negative	2	4	0	0	2	8
N/A	7	15	2	9	5	21
PR status
Positive	33	70	17	74	16	67
Negative	7	15	4	17	3	13
N/A	7	15	2	9	5	21
HER2 status
Positive	2	4	2	4	0	0
Negative	22	47	11	23	11	46
N/A	23	49	10	21	13	54
BRCA2 mutation status
Positive	3	6	2	9	1	4
Negative	5	11	4	17	1	4
N/A	39	83	17	74	22	92
Histological grade
I	1	2	0	0	1	4
II	17	36	6	26	11	46
III	12	26	8	35	4	17
N/A	17	36	9	39	8	33
Lymph node status
N0	15	32	7	30	8	33
N+	26	55	13	57	13	54
N/A	6	13	3	13	3	13
Adjuvant endocrine therapy
Yes	26	55	13	57	13	54
No	16	34	7	30	9	38
N/A	5	11	3	13	2	8
Adjuvant chemotherapy
Yes	4	9	1	4	3	13
No	38	80	19	83	19	79
N/A	5	11	3	13	2	8
Post‐operative radiotherapy
Yes	21	45	11	48	10	42
No	22	47	10	43	12	50
N/A	4	8	2	9	2	8

Fresh frozen tumor tissue from FBC was obtained from the South Sweden Breast Cancer Group tissue bank at the Department of Oncology, Skåne University Hospital (Lund, Sweden). DNA was extracted as described previously (Holm et al., 2010; Holm et al., 2015). The 11 normal tissue samples were comprised of two non‐malignant breast tissue samples from two female patients with breast cancer and nine non‐malignant cell lines analyzed in duplicate including adult dermal fibroblasts, adult epidermal keratinocytes, dark, light and medium epidermal melanocytes, mammary epithelial cells, mammary fibroblasts, mammary endothelial cells and bone marrow derived mesenchymal stem cells (ScienCell Research Laboratories, Carlsbad, CA).

The study was approved by the regional Ethics Committee in Lund (2012/89).

Illumina 450K methylation arrays

DNA quality and concentration were assessed with a NanoDrop ND1000 (Thermo Scientific, Waltham, MA) and DNA was bisulfite‐converted using the EZ DNA Methylation Kit (Zymo Research, Orange, CA) following the manufacturer's protocol. Samples were run on Illumina's Infinium HumanMethylation450K BeadChips (Illumina, CA) at the SCIBLU Microarray Facility at Lund University as described elsewhere (Sandoval et al., 2011), providing the methylation status of more than 480,000 cytosines distributed over the whole genome (Dedeurwaerder et al., 2011; Sandoval et al., 2011), In short, MBC, FBC and normal tissue samples were whole‐genome amplified, enzymatically digested and hybridized to the arrays and single nucleotide extension was performed. Finally, the chips were scanned using a two‐color channel Illumina HiScan SQ scanner and the array output was read with Illumina GenomeStudio (2011.1) and background normalization was conducted. The data were exported and read into R (R Development Core Team, 2013). Samples with >90,000 (MBC) or >10,000 (FBC) missing probes with a detection p‐value >0.05 were removed, resulting in 47 MBC, 188 FBC and 11 normal tissue samples. The methylation score for each probe was represented with a β‐value calculated from the raw intensity of methylated (M) and unmethylated (U) signal:

The β‐value may thus be any value between 0 (completely unmethylated) and 1 (completely methylated). β‐values that could not be calculated were imputed by using the R package impute. The Illumina 450K array uses two types of Infinium assay designs (Infinium type I and II), which were normalized individually using a peak‐based correction, similar to that described by Dedeurwaerder et al. (2011). Briefly, the β‐values were smoothed using Epanechnikov smoothing kernel for both assay design types to estimate unmethylated and methylated peaks and a linear scaling was used to move the unmethylated peak to 0 and the methylated peak to 1, with β‐values in‐between hence stretched. β‐values above 1 were set to 1 and β‐values below 0 were set to 0. The FBC and MBC samples were handled in three batches during bisulfite treatment and hybridization and samples were randomized across batches such that each batch represented all types of male and female breast cancer, making it possible to adjust for the batches. A principal component analysis (PCA) was run and associations between principal components and technical and biological annotations were evaluated for detection of potential technical confounding factors (Lauss et al., 2013), (Supplementary Figure S1A). The batches were identified as a technical confounding factor and this was corrected for using a supervised empirical Bayes method (ComBat, Supplementary Figure S1B) (Johnson et al., 2006). The β‐values above 1 were again reset. When MBC samples were analyzed separately the dataset will from here on be called ‘MBC dataset’ and the dataset containing the combined MBC and FBC tumors will be called ‘MBC and FBC dataset’. Probes in the MBC dataset containing SNPs at the target CpG were removed, as were probes cross‐hybridizing between autosomes and sex chromosomes, resulting in 425,240 probes (Price et al., 2013). For subsequent analyses, β‐values for probes with a p‐value >0.05 were defined as missing values. For the combined MBC and FBC dataset, the sex chromosome probes, probes that contained SNPs at the target CpG and probes that cross‐hybridized between autosomes and sex chromosomes were removed, resulting in 415,080 probes (Price et al., 2013).

Finally, β‐values were divided into three categories: hypomethylated = 0 (β‐values of 0 to ≤0.3), heterogeneously methylated = 0.5 (β‐values of >0.3 to <0.7) and hypermethylated = 1 (β‐values of ≥0.7–1.0). Hyper‐ and hypomethylation frequencies for each sample were calculated as the fraction of CpGs with the value 1 or 0, respectively.

Gene expression and aCGH datasets

Gene expression (GEX) and array comparative genomic hybridization (aCGH) data were available for the MBC tumors, and these datasets were normalized as previously described (Johansson et al., 2012, 2011). Low varying features in both the gene expression and methylation datasets were removed and all CpG probes and gene expression probes were matched on gene symbols, resulting in 4374 unique genes and 34,916 unique CpG probes. The Pearson correlation between gene expression and methylation for each unique gene was calculated and permuted to obtain data‐driven false discovery rates.

Data analysis

The 2.5% most variable CpG probes (n = 10,507) were selected from the MBC dataset and unsupervised clustering was performed using a recursively partitioned mixture model (RPMM) from the RPMM R‐packages (Houseman et al., 2008) to identify subgroups of MBC. RPMM is a model based unsupervised clustering method specifically developed for beta‐distributed DNA methylation data (Houseman et al., 2008) and has been successfully used to identify stable epitype subgroups in diverse tumor types (Avissar‐Whiting et al., 2011, 2009, 2010, 2012, 2009). For comparison, unsupervised hierarchical clustering (HCL) using Euclidean distance and Ward's algorithm was also performed. The ‘MBC and FBC dataset’ was clustered on the 10,000 CpGs that varied the most across MBC and FBC tumors.

To further understand the biology of the MBC epitypes, a rank‐based module activity score was calculated for each MBC tumor for eight gene expression modules reflecting biologically relevant transcriptional programs found in FBC as previously described (Fredlund et al., 2008). Scores for seven gene expression modules associated with key biological processes in FBC were also calculated to discover biologically meaningful differences between MBC epitype subgroups as described previously (Desmedt et al., 2008). Annotations for genomic locations of CpG probes were derived from Illumina (Illumina, CA). A set of 654 PRC2 target genes in embryonic stem cells identified by Lee et al. using chromatin immunoprecipitation arrays as bound by SUZ12, EED and H3K27me3, was also used (Lee et al., 2006). We also annotated all probes for embryonic stem cell chromatin states (poised promoter, active promoter, enhancer, transcribed, insulator and heterochromatin, respectively) (Ernst et al., 2011).

We performed BioFunctional analyses using AMIGO (http://amigo.geneontology.org/amigo) GO annotations and Molecular signature database (MSigDB) enrichment p‐values were calculated using Fisher's exact test. The transcription factor binding site analysis program Systematic Motif Analysis Retrieval Tool (SMART) was used as previously described to search for transcription factor binding site (TFBS) motifs among differentially methylated genes. Promoter regions were defined as the genomic interval from −1500 to +500 bp relative to the putative transcription start sites (Veerla et al., 2010).

All the statistical calculations and figures were made in R (R Development Core Team, 2013). All P‐values are two‐sided. The survival analysis was performed using the R‐package survcomp with distant disease free survival (DMFS) as end‐point.

Immunohistochemistry of EZH2

A tissue microarray (TMA) with two 1 mm cores from each of 220 MBC tumors was constructed as described (Nilsson et al., 2011). Sections of 3–4 μm were cut, transferred to SuperFrost Plus slides, dried at room temperature and then baked for 2 h at 60 °C. The DAKO Envision horseradish peroxidase rabbit/mouse kit (DAKO, Glostrup, Denmark) and a Dakocytomation Autostainer (DAKO) were used for the staining procedure. The TMAs were stained with a purified mouse anti‐EZH2 monoclonal antibody (clone 11, BD Transduction Laboratories, Franklin Lakes, NJ) at a 1:25 dilution after antigen retrieval at high pH as described elsewhere (Holm et al., 2012). EZH2 staining was scored by one reader (IJ) in a blinded manner. The percentage of positively stained tumor cells was evaluated and scored as: 0 (0%), 1 (1–10%), 2 (11–25%), 3 (26–50%), 4 (51–75%) or 5 (>75%). Tumors were considered positive for EZH2 if >50% of the cancer cells showed nuclear staining (Supplementary Figure S2).

Results

Unsupervised clustering revealed two stable epitypes of male breast cancer

RPMM clustering of the ∼10,000 most variable CpG probes among the MBC tumors revealed two stable and equally sized epitype groups, from here on called ME1 (23 tumors) and ME2 (24 tumors) (Figure 1). Varying the number of most variable probes among the MBC tumors between 5000 and 13,000 resulted in identical groups (data not shown). HCL was also performed using the same CpG probes and the groups identified were compared with the RPMM groups to investigate the stability of the groups. This comparison resulted in significant overlap, with 96% (45/47) of the tumors showing agreement in cluster assignment. Furthermore, the two epitypes were associated with the previously defined gene expression subgroups luminal M1 and M2 (Fisher's exact test, P = 0.011) (Johansson et al., 2012). Although the methylation derived and transcriptional subgroups were associated and the ME1 group contained the majority of the luminal M1 tumors (87%), the ME2 group contained a mixture of MBC tumors from the luminal M1 and M2 transcriptional subgroups. The normal tissues showed a homogenous methylation pattern of the most varying probes across MBC tumors (Figure 1), and ME2 tumors displayed more similarities with normal cells across these CpG probes than ME1 tumors.

Figure 1

RPMM clustering of MBC tumors revealed two stable epitypes. The large heatmap to the left displays the ∼10,000 most variable CpG probes in the 47 MBC tumors, and the heatmap to the right displays the same CpGs in 11 normal tissues (two non‐malignant breast tissue samples from two females and nine cell lines: adult dermal fibroblasts, adult epidermal keratinocytes, dark, light and medium epidermal melanocytes, mammary epithelial cells, mammary fibroblasts, mammary endothelial cells and mesenchymal bone marrow stem cells). β‐values in the heatmaps range from unmethylated (blue) to methylated (yellow). The bars on the left hand represent CpG islands (black = island; dark grey = shore; light grey = shelf; white = open sea) and embryonic stem cell chromatin states (purple = poised promoter; red = active promoter; yellow = enhancer; green = transcribed; blue = insulator; white = heterochromatin), respectively.

Clinical and biological context of the MBC epitypes

A Kaplan–Meier survival analysis suggested inferior survival for the ME1 epitype, although this was not statistically significant (Figure 2A; P = 0.22). The ME1 group also displayed a significantly higher frequency of tumor cells positive for the proliferation marker cyclin A (P = 0.012, Wilcoxon test), a higher fraction of genome altered (FGA; P = 0.0045, Wilcoxon test) and a higher S‐phase fraction (P = 0.035, Wilcoxon test) (Figure 2B–D). Taken together, these findings indicate that ME1 tumors are more proliferative and aggressive, in line with previous findings of inferior outcome within the transcriptionally defined luminal M1 subgroup compared to the luminal M2 subgroup (Johansson et al., 2012). However, no difference in age at diagnosis, tumor size, Nottingham Histological grade (NHG), ER, PR or number of lymph node metastases was found between the epitypes (Table 1).

Figure 2

Differences between the two MBC epitypes. (A) Kaplan–Meier survival analysis of distant metastasis free survival for the MBC patients stratified by epitype. (B) Fraction of cyclin A positive cells in the two epitypes. (C) Fraction of genome altered (FGA) in the two epitypes. (D) S‐phase fraction in the two epitypes.

To further investigate the biology of the epitypes, we examined gene modules previously defined as biologically relevant in breast cancer (Desmedt et al., 2008; Fredlund et al., 2012). From the Fredlund study, only the two modules representing proliferation – mitotic checkpoint and mitotic progression – were significantly different between the epitypes, with the ME1 tumors displaying the highest activity (Figure 3A–B; P = 2.2e‐06 and P = 1.5e‐05, respectively, Wilcoxon test). The activity of the proliferation module from Desmedt et al. was also higher in ME1 tumors than in ME2 tumors (Figure 3C; P = 2.8e‐7, Wilcoxon test), further strengthening the picture of ME1 tumors constituting a more proliferative group of MBC. In addition, the ER score was lower among ME1 tumors (Figure 3D; P = 0.048, Wilcoxon test), again indicating that the epitypes capture the heterogeneity among MBCs in a manner similar to gene expression profiles (Johansson et al., 2012). A SAM analysis was performed between the epitypes, resulting in 4674 CpGs significantly more methylated in ME1 tumors compared to ME2 tumors (FDR < 0.05). The most significant functional annotations associated with these genes involved the epigenetic gene silencing mark H3K27me3 (P = 4.4e‐153) and transcriptional regulation by HOX/homeobox genes (P = 1.6e‐22). Further functional pathway analyses revealed significant associations with the canonical pathways cell adhesion (P = 5.6e‐5), WNT signaling (P = 2.8e‐4), TGFβ signaling (P < 0.001), focal adhesion (P < 0.005), MAPK signaling (P < 0.005), FGFR ligand binding and activation (P < 0.007) and hedgehog signaling (P < 0.007). A significant enrichment of genes targeted by the polycomb‐group (PcG) family (including EED, SUZ12, PRC2) involved in maintaining transcriptional repression was observed among genes hypermethylated in ME1 compared to ME2 MBCs, suggesting differences in transcriptional regulation mechanisms between the epitypes. A search for TFBS motifs among CpG probes hypermethylated in ME1 vs. ME2 tumors converged on transcription factors within the MAP/ERK and protein kinase families, including among the top ten hits PRKCα, MAPK1, PKCα and PKCβ, and the histone demethylase KDM5B (also known as JARID1B) (all P‐values<5.0e‐7; Fisher's exact test). A complete list of transcription factors and their co‐factors enriched among differentially methylated genes is provided in Supplementary Table S3.

Figure 3

Functional module activity scores for the two MBC epitypes. Module activity scores for mitotic checkpoint (A) and mitotic progression (B) from Fredlund et al. (2012), and proliferation (C) and ER signaling (D) from Desmedt et al. (2008), respectively, in the two epitypes.

Male breast cancer epitypes and polycomb regulated genes

MBC tumors of the ME2 epitype displayed significantly more hypomethylated CpGs than ME1 epitype tumors (P = 0.0038, Wilcoxon test), while ME1 tumors showed a weak trend toward more hypermethylated CpGs (P = 0.091, Wilcoxon test). ME1 tumors specifically showed hypermethylation on poised promoters whereas ME2 tumors and normal tissue samples did not (Figure 1). Interestingly, there was a significant enrichment of poised promoters among the CpG island probes (regions rich in CpG sites close to the promoter regions of the genes) that varied the most (P = 2.4e‐5, Fisher's exact test). To explore this further, the genomic distribution of the differentially methylated probes was plotted, revealing significant differences between the two MBC epitypes. Hypermethylation within transcription start sites/islands was more common in tumors of the ME1 epitype, while CpGs in shelves and open seas were significantly less hypermethylated in these tumors (Figure 4; P = 5e‐33, and P = 8e‐260, respectively, Fisher's exact test). No CpGs were more methylated in ME2 tumors compared to ME1 tumors (FDR<0.05). Furthermore, when analyzing the methylation level of the 13,937 probes from the 654 PRC2 target genes identified by Lee et al. (Lee et al., 2006), the ME1 group had a significantly higher average β‐value (P = 1.3e‐7, Wilcoxon test) and correspondingly the gene expression level of the PRC2 target genes displayed the opposite pattern (Figure 5A–B; P = 7.2e‐5, Wilcoxon test). Furthermore, the expression of the core member of PRC2, EZH2, was also significantly increased in the ME1 group (Figure 5C; P = 3.3e‐7, Fisher's exact test). Based on the findings that high mRNA levels of EZH2 correlated to the ME1 epitype, we investigated the protein expression of EZH2 in 220 MBC tumors assembled in a TMA; 79/220 (36%) tumors were EZH2 positive and 122/220 (55%) tumors were EZH2 negative (data were missing for 19 (9%) tumors). As has previously been shown in breast cancer as well as other tumor types, EZH2 positive MBC tumors had significantly higher mRNA levels of EZH2 than EZH2 negative tumors (P = 0.0011, Wilcoxon test) and a difference in the protein expression of EZH2 was suggested between the subgroups, with more ME1 tumors being EZH2 positive compared to ME2 tumors (P = 0.056, Fisher's exact test). No difference in DMFS was seen between EZH2 positive and negative tumors (data not shown).

Figure 4

Genomic distribution of methylated CpGs in MBC tumors of the ME1 epitype. (A) Probes annotated based on gene components, promoters (TSS1500, TSS200), first exon, 5′UTR, body and 3′UTR or intergenic (P = 5e‐33). (B) Probes annotated based on the categorical regions islands, shores, shelves and open sea (P = 8e‐260). Blue, hypermethylatyed in ME1 compared to ME2; grey, SAM background.

Figure 5

PRC2 targets in the two MBC epitypes. Average relative methylation (A) and average relative expression (B), respectively, of polycomb probes in the two epitypes of MBC; and relative expression of EZH2 (C) in the two epitypes of male breast cancer.

DNA methylation and gene expression

Next, the correlation between methylation and gene expression was studied globally in the MBC tumors. Among the 34,916 unique CpG probes with matched gene expression data, 4914 negatively correlated CpGs (Pearson correlation ≤ −0.38, FDR <0.05) and 1531 positively correlated CpGs were identified (Pearson correlation ≥0.38, FDR <0.05). As expected, the majority of the CpGs that correlated (76%) were hence negatively correlated to gene expression levels and significantly more of the negatively correlated CpGs were island probes compared to the positively correlated CpGs (P = 1.1e‐7, Fisher's exact test). Among these, many HOX genes, including HOXB2, were negatively correlated, with Pearson correlations < −0.5.

Overall methylation patterns and co‐clustering of male and female breast cancers

Next, male and female breast cancers were compared to establish the degree of similarity, and whether MBCs resemble the previously described methylation subgroups in breast cancer (Holm et al., 2010). Genome‐wide DNA methylation data were obtained from 47 MBCs, 188 FBCs and 11 normal tissue samples. No significant overall difference in methylation of the most variable CpG probes between MBC and FBC was observed, as 41.6% of probes were hypermethylated and 39.7% hypomethylated in MBC tumors, and 41.6% of probes were hypermethylated and 38.9% hypomethylated in FBC tumors. All breast cancers displayed comparable methylation patterns across the genome, with low methylation levels near transcription start sites and high methylation levels in gene bodies, 3′UTRs and intergenic regions (Supplementary Figure S4).

Co‐clustering of the MBC and FBC tumors based on the 10,000 most differentially methylated probes across MBC and FBC resulted in the MBC tumors being split among the two clusters of ER positive (luminal) FBCs. 21/47 (45%) MBCs clustered together with mainly luminal A FBCs (middle cluster in Figure 6), and 26/47 (55%) MBCs clustered among the predominantly luminal B FBCs (right cluster in Figure 6), while none of the MBCs clustered within the basal FBC cluster (left cluster in Figure 6). In addition, the majority of the MBCs in the luminal B (right) cluster were of the more proliferative/aggressive ME1 epitype, while MBCs of the ME2 epitype were closely correlated with the predominantly luminal A (middle) cluster, confirming the stability of the MBC epitypes and the established correlations between breast cancer subtypes and disease aggressiveness.

Figure 6

Hierarchical clustering of MBC and FBC tumors. The heatmap displays the 47 MBC and 188 FBC tumors clustered on 10,000 CpG probes. β‐values in the heatmap range from unmethylated (blue) to methylated (yellow).

Interestingly, the majority of the MBC tumors were grouped together within the respective clusters rather than being interspersed among the FBC tumors, suggesting that they nevertheless differed from the FBC tumors in the respective clusters to some extent (Figure 6). In an attempt to investigate whether differences in ER status may account for this separation, we constructed a predictor of ER status in the FBCs based on genes with reported methylation differences between ER positive and negative tumors (Rønneberg et al., 2011) using promoter CpGs (up to 1500 upstream of TSS) and applied it to the MBCs. Using leave‐one‐out cross validation (LOOCV), 77% of the FBCs were correctly classified, while 11% of the ER positive cases were incorrectly assigned to the ER negative class, and 12% of the ER negative cases were incorrectly assigned to the ER positive class (data not shown). Among 40 MBCs with available ER status (38 of which were ER positive and 2 ER negative by IHC) 32 were correctly classified; 6 ER positive cases were misclassified as ER negative and the 2 ER negative cases were incorrectly classified as ER positive (Supplementary Figure S5). All of the misclassified cases were of the ME2 epitype. The overall misclassification rate of 20% was hence similar to the performance of the predictor in FBC.

Supervised analyses of male vs. female breast cancer epitypes

To explore the potential differences between MBCs and FBCs in the two clusters, SAM analyses were performed to identify differentially methylated probes. Elucidation of the locations of differentially methylated probes revealed less hypermethylated probes in islands (16.0% vs. 26.6%) and more hypermethylated probes in open seas (58.4% vs. 47.0%) in MBCs compared to FBCs in the middle (luminal A/ME2) cluster (Figure 7A–B; P = 8e‐9, Fisher's exact test). Conversely, more hypermethylated islands (28.7% vs. 21.1%) and less hypermethylated open seas (41.9% vs. 48.6%) were observed in MBCs compared to FBCs in the right (luminal B/ME1) cluster (Figure 7C–D; P = 3e‐4, Fisher's exact test), suggesting differential methylation patterns not only between the MBC epitypes, but also between MBCs and FBCs that cluster together.

Figure 7

Genomic distribution of methylated CpGs in FBC vs. MBC tumors in the middle cluster (A–B), and right cluster (C–D) fromFigure 6, respectively. Probes were annotated based on gene components, promoters (TSS1500, TSS200), first exon, 5′UTR, body and 3′UTR or intergenic (A, P = 0.05; C, P = 3e–12) or based on the categorical regions islands, shores, shelves and open sea (B, P = 8e–9; D, P = 3e–4). Purple, hypermethylated in FBC compared to MBC; green, hypermethylated in MBC compared to FBC.

Finally, to explore whether functional or molecular pathways were differentially affected by this diversity, pathway analyses comparing MBC and FBC tumors in the middle cluster and in the right cluster in Figure 6, respectively, were performed. Hypermethylation of CpG sites in genes enriched in the NOTCH (P = 8e‐4), ERα (P = 0.0011), HDAC (P = 0.0052), C‐MYB (P = 0.017), MTOR (P = 0.022), ETS (P = 0.022), TGFβ (P = 0.03), MAPK (P = 0.035), IL1 (P = 0.036) and MYC (P = 0.048) canonical pathways were associated with FBCs compared to MBCs in the middle cluster in Figure 6 containing mainly luminal A FBCs and ME2 MBCs, as well as genes associated with basal vs. luminal breast cancer (P = 1.1e‐10), H3K27me3 (P = 6.4e‐7), apocrine vs. luminal breast cancer (1.7e‐6) and endocrine/tamoxifen therapy resistance (P = 1.0e‐5). Gene ontology analyses comparing MBC and FBC tumors in the right (predominantly luminal B/ME1 MBCs) cluster in Figure 6 revealed hypermethylation of genes involved in the same canonical pathways as in the middle cluster, as well EED (P = 4.0e‐5) and SUZ12 (P = 6.7e‐5) targets and genes related to response to androgen (P = 0.0012) and TP53 targets (P = 0.0081) in FBCs compared to MBCs.

Discussion

MBC is a disease that is not well studied on the genomic level and moreover, no global studies on the epigenetic level have been performed to date. Epigenetic changes are stable and heritable during cell division (Jovanovic et al., 2010), and cytosine methylation is one of the most commonly studied epigenetic changes in cancer (Esteller, 2007). Using Illumina 450K Infinium arrays to study global methylation patterns in MBC we identified two stable and equally sized epitypes among the MBC tumors through unsupervised clustering, ME1 and ME2. The epitypes were associated with the previously described transcriptional subgroups luminal M1 and M2. We previously showed that luminal M1 tumors appeared more aggressive and were associated with worse prognosis and also appeared to have a less activated ER pathway, while luminal M2 tumors displayed an up‐regulated immune response and a more activated ER pathway, despite >90% of the MBCs being ER positive by immunohistochemistry (Johansson et al., 2012). The ME1 tumors, which were associated with the luminal M1 subgroup, also had a low score for the ER module (Desmedt et al., 2008). There was no significant difference in survival between the epitypes, although ME1 tumors showed a trend toward worse prognosis compared to ME2 tumors. ME1 tumors were also found to be more proliferative, as reflected in significantly higher levels of the proliferation marker cyclin A and higher S‐phase fraction. Higher mitotic checkpoint activity, known to correlate with poor prognosis in luminal FBC (Fredlund et al., 2012), was also observed among ME1 tumors. Furthermore, the higher fraction of genome altered observed in ME1 tumors is in line with chromosomal instability in breast cancer cells being associated with increased expression of mitotic checkpoint genes (Yuan et al., 2006). All these findings indicate that, among MBCs, ME1 tumors are more proliferative and genomically unstable and may be more aggressive than ME2 tumors. Moreover, ME2 tumors displayed methylation patterns more similar to normal cells than ME1 tumors. Overall, ME1 tumors displayed a higher frequency of hypermethylation and lower frequency of hypomethylation than ME2 tumors. Interestingly, hypermethylation of CpGs in promoter islands was significantly more common in ME1 tumors, potentially also suggesting different mechanisms of gene regulation among some MBCs. In line with these data, a recent study showed that high promoter methylation frequencies of 25 tumor suppressor genes were associated with an aggressive phenotype and poor survival in a cohort of approximately 100 MBCs (Kornegoor et al., 2012).

The most significant GO terms associated with genes found to be more methylated in ME1 tumors in relation to ME2 tumors included homeobox and HOX genes. HOX genes have been implicated in a variety of cellular processes, and have been shown to be hypermethylated in aggressive tumors, including breast cancer (Fackler et al., 2011), urothelial cancer (Lauss et al., 2012), ovarian cancer (Fiegl et al., 2008) and leukemia (Strathdee et al., 2007). Of interest, genes with a strong negative correlation between methylation and expression included HOX genes, further supporting that HOX genes may be silenced by methylation in MBC tumors. Similar findings have been reported in urothelial cancer, where a strong negative correlation was found for HOXB2 (Lauss et al., 2012). Furthermore, HOXB2 was found to be a negative regulator of tumor growth in both basal and luminal breast cancer cell lines in a functional in vivo screen (Boimel et al., 2011). Taken together, these findings imply that down‐regulation of HOXB2 by hypermethylation may be one mechanism whereby certain MBCs regulate tumor growth.

Trimethylation of H3K27 is catalyzed by PRC2 and a set of 654 PRC2 target genes displayed significantly higher average relative DNA methylation levels in the ME1 tumors. Additionally, ME1 tumors also displayed specific hypermethylation of poised promoters, a set of promoters that are bivalently marked by H3K27 and H3K4 trimethylation in embryonic stem cells. The ME1 tumors also showed significantly lower levels of expression of PRC2 target genes, indicating that they could be silenced by de novo methylation. It has been shown previously that genes hypermethylated in tumors compared to normal cells are mainly polycomb target genes (Wolff et al., 2010; Kron et al., 2013; Avissar‐Whiting et al., 2011), and also that polycomb target methylation occurs in more aggressive forms of various cancer types (Kron et al., 2013; Lauss et al., 2012). ME1 tumors in the present study displayed significantly higher average mRNA and protein levels of EZH2, suggesting that ME1 tumors may have more PRC2, and therefore most likely a higher degree of H3K27me3. This finding is supported by the functional pathway enrichment analyses, demonstrating a significant enrichment of H3K27me3 among CpGs more methylated in ME1 tumors compared to ME2 tumors. These findings indicate that developmental processes in the ME1 tumors may be repressed by H3K27me3 through the PRC2 complex and then further repressed by de novo DNA methylation of PRC2 target genes. Methylation, and thereby decreased expression of polycomb target genes in conjunction with high mRNA and protein levels of EZH2 may be one mechanism whereby more undifferentiated and aggressive tumors develop. An inverse correlation between EZH2 expression and global H3K27 trimethylation has been demonstrated in breast cancer, and high expression of EZH2/low H3K27me3 has been associated with inferior survival (Holm et al., 2012). Functional pathway analyses also identified differential activation of the polycomb family group members, supporting the importance of their role in MBC tumorigenesis. Moreover, significant differences in the activation of fundamental signaling pathways in cancer, including the WNT, TGFβ, MAPK and hedgehog signaling pathways as well as cellular adhesion, involved in tumor aggressiveness and metastatic dissemination were identified between MBC tumors of the two epitypes, potentially suggesting significant differences in tumor development and supporting a more aggressive phenotype in tumors of the ME1 epitype. In line with this, TFBS motif analyses identified the transcription factor JARID1B, a luminal lineage driving oncogene associated with poor prognosis in breast cancer (Yamamoto et al., 2014), in ME1 tumors. Interestingly, JARID1B has also been implicated in regulation of genome stability through its role in promoting double‐strand break signaling (Li et al., 2014), and ME1 tumors displayed a more unstable genome, with higher FGA levels, than ME2 tumors. Finally, JARID1B, also known as the H3K4 demethylase KDM5B, has been shown to control epithelial–mesenchymal transition (EMT) in cancer cells (Enkhbaatar et al., 2014), a key event in tumorigenesis. These findings together suggest that JARID1B may be a novel molecular target for therapy in e.g. MBC.

Overexpression of EZH2 has been shown in many tumor types including breast, urothelial and prostate cancer and is associated with more aggressive disease (Chang and Hung, 2011; Holm et al., 2010; Lauss et al., 2012; Varambally et al., 2002). The highest levels of EZH2 in FBC have been found in ER negative tumors (Holm et al., 2010; Collett et al., 2006; Puppe et al., 2009), although when ER positive luminal A and B tumors were compared, EZH2 expression was significantly higher in luminal B tumors (Holm et al., 2010). EZH2 is not only associated with cancer stem cell formation and expansion of an aggressive cancer stem cell population promoting cancer progression and metastasis (Chang et al., 2011), it also regulates several target genes with a variety of functions including promoting EMT and tumor angiogenesis (Cao et al., 2008; Chang and Hung, 2011; Lu et al., 2010). Specifically blocking EZH2 expression or activity in tumor cells may thus represent a promising strategy for anti‐cancer treatment (Chang and Hung, 2011). In fact, McCabe et al. showed in a recent study that GSK126, a small molecule inhibitor of EZH2 methyltransferase activity, decreased global H3K27me3 levels and reactivated silenced PRC2 target genes in lymphoma (McCabe et al., 2012).

Since the MBC and FBC tumors were analyzed in a single experiment, a direct comparison and analysis of male and female breast cancer could be performed. HCL clustering of the male and female breast cancers together based on the 10,000 most varying CpG probes revealed that MBC tumors clustered among ER positive, luminal FBCs. Three main clusters were identified: one with mainly basal FBC tumors, one with the majority of the luminal B FBC tumors and ME1 MBC tumors, and one with the majority of the luminal An FBC tumors and ME2 MBC tumors. This is in line with previous reports of DNA methylation in FBC; three groups associated with the luminal A, luminal B and basal FBC tumors were identified, while HER2 and normal‐like FBC tumors were distributed among the three groups (Holm et al., 2010). Considering that the vast majority of the MBC tumors were ER positive, it is not surprising that they clustered among ER positive (luminal) FBC tumors. Nevertheless, while luminal A and B FBC tumors are known to exhibit active ER signaling, ME1 MBC tumors displayed a low ER module activity score, in line with the observation in the transcriptionally derived subgroup luminal M1, which correlated with ME1 herein (Johansson et al., 2012). A low ER score is characteristic of basal FBC tumors and although the ME1 tumors had a low ER score, it was not as low as for basal FBC tumors and the ME1 tumors did not cluster with basal FBC tumors on the methylation level. MBCs of the ME2 epitype had a high score for the ER module, similar to luminal A and B FBC tumors, which they also mainly clustered together with on the methylation level. Using a previously published predictor of ER status based on CpGs associated with ER status in FBC (Rønneberg et al., 2011), the vast majority of MBCs in the present study were assigned to the ER positive class, in line with the IHC based assessment, indicating a strong correlation between ER expression and methylation patterns shared across both male and female breast cancer. Nonetheless, the MBC tumors tended to cluster together rather than being inter‐mixed with the FBC tumors in the respective clusters, and significant differences in key signaling pathways, including the ERα, NOTCH, HDAC, MTOR, ETS, TGFβ and MYC canonical pathways, as well as gene signatures distinguishing between basal and luminal breast cancer and endocrine therapy response between MBCs and FBCs clustering together suggest that MBCs in fact display features representative of the whole spectrum of breast cancer in women, and may explain the phenotypic and clinical heterogeneity of MBC. Differences in the genomic distribution of hypermethylated CpGs between MBCs and FBCs also support the notion that DNA methylation in subsets of male and female breast cancers may be associated with different chromatin states, with ME1 MBCs displaying more frequent hypermethylation of CpGs located in promoter islands (and less in open seas) compared to the predominantly luminal B FBCs with which they clustered, and inversely, ME2 MBCs displaying less hypermethylation of CpGs located in promoter islands (and more in open seas) compared to the luminal A FBCs with which they clustered.

Taken together, these findings suggest that MBC and FBC are not only different on the genomic and transcriptional levels, but also harbor differences on an epigenetic level. Although the MBCs were divided between the predominantly luminal A and B FBC clusters, respectively, significant differences in genomic distribution of methylated CpGs and affected functional pathways were apparent between breast cancers in men and women, likely reflecting fundamental underlying differences.

Conclusions

In conclusion, we have identified two epitypes among MBC tumors, ME1 and ME2. The ME1 epitype was potentially more aggressive, with specific hypermethylation of polycomb target genes and high expression of EZH2. MBC tumors clustered together with luminal FBC tumors; however, despite clustering together, the MBC tumors were grouped together within the clusters rather than being interspersed among the FBC tumors. These findings underscore the heterogeneity of MBCs and that they display features observed across the spectrum of FBCs; they are therefore not readily defined using conventional criteria applied to FBC.

The following are the supplementary data related to this article:

Supplementary Table S3 Transcription factor enrichment analysis based on genes hypermethylated in ME1 vs. ME2 MBC tumors.

Click here for additional data file.

Supplementary Figure S1 A principal component analysis (PCA) was performed and associations between principal components and technical and biological annotations were evaluated, whereupon a platform specific bias was detected in the main principal component (A). After adjustment using ComBat, no technical variation was found among the main principal components (B).

Click here for additional data file.

Supplementary Figure S2 Representative IHC staining of EZH2 protein in MBC illustrating negative (A) and positive (B) staining, respectively.

Click here for additional data file.

Supplementary Figure S4 Mean DNA methylation levels across the genome in MBC (A) and FBC (B) tumors, respectively. Mean β‐values for different types of probes are shown.

Click here for additional data file.

Supplementary Figure S5 Classification of MBC tumors using a predictor for ER status derived from FBC tumors (Rønneberg et al. 2011). Euclidean distances for ER positive and ER negative predictions, respectively, are plotted.

Click here for additional data file.