Abstract Genome-wide methylation profiling is used in breast cancer (BC) studies, because DNA methylation is a crucial epigenetic regulator of gene expression, involved in many diseases including BC. We investigated genome-wide methylation profiles in both canine mammary tumor (CMT) tissues and peripheral blood mononuclear cells (PBMCs) using reduced representation bisulfite sequencing (RRBS) and found unique CMT-enriched methylation signatures. A total of 2.2–4.2 million cytosine–phosphate–guanine (CpG) sites were analyzed in both CMT tissues and PBMCs, which included 40,000 and 28,000 differentially methylated regions (DMRs) associated with 341 and 247 promoters of differentially methylated genes (DMGs) in CMT tissues and PBMCs, respectively. Genes related to apoptosis and ion transmembrane transport were hypermethylated, but cell proliferation and oncogene were hypomethylated in tumor tissues. Gene ontology analysis using DMGs in PBMCs revealed significant methylation changes in the subset of immune cells and host defense system-related genes, especially chemokine signaling pathway-related genes. Moreover, a number of CMT tissue-enriched DMRs were identified from the promoter regions of various microRNAs (miRNAs), including cfa-mir-96 and cfa-mir-149, which were reported as cancer-associated miRNAs in humans. We also identified novel miRNAs associated with CMT which can be candidates for new miRNAs associated with human BC. This study may provide new insight for a better understanding of aberrant methylation associated with both human BC and CMT, as well as possible targets for methylation-based BC diagnostic markers. Keywords: DNA methylation, canine mammary tumor, miRNA in cancer 1. Introduction Breast cancer (BC) is one of the most frequently diagnosed cancers in women, comprising 25% of total cancer in women [[30]1]. Although many studies, projects, and clinical trials demonstrated the association between BC and some mutations on cancer driver genes such as BRCA1/2 and PIK3CA, as well as familial history, only 5–10% of BC can be explained by genetic aberrations while the rest remains unclear [[31]2]. Cancer studies using laboratory animal models, specifically rodents, were developed in diverse ways, such as the use of transgenic and targeted knockout mice [[32]3]. However, it is difficult to mimic human environmental conditions using a laboratory animal model system. Companion animals such as dogs and cats, thus, became great model animals because of their closely shared environment with their owners [[33]4]. Furthermore, canine mammary tumor (CMT) is a well-known animal model for BC because of its similarities in genetic, biological, and clinical features to BC [[34]5]. Benign and malignant CMTs occur at a similar frequency, but classification was subdivided into malignant tumor subtypes by the World Health Organization (WHO). The 2011 classification subdivided CMTs into seven benign and 23 malignant subtypes, according to their prognostic values. Complex carcinoma is the most frequently diagnosed type of CMT (~20%), and it was subdivided into five subtypes (carcinoma complex type, carcinoma mixed type, carcinoma and malignant myoepithelioma, adenosquamous carcinoma, and intraductal papillary carcinoma) by the 2011 classification [[35]6]. It is well known that the complex type of CMT, which includes a myoepithelial cell population and comprises ~20% of total CMT, is a great model for human myoepithelial cell-type BC, since myoepithelial cell-type BC is so rare in humans (~0.1%) [[36]7]. Over the last decade, tumor immunity became a focus point due to the ability of immune systems to detect and destroy various cancer cells. Although the immune system mainly has a protective role, it does not always suppress tumor growth and dissemination successfully [[37]8]. On the contrary, it was reported that tumor cells can co-opt immune cells such as macrophages into the tumor stroma, resulting in the promotion of endothelial cell proliferation in the tumor stroma and angiogenesis, which can enhance tumor cell mobility and invasiveness [[38]9,[39]10]. Moreover, platelets, tumor-associated macrophages, and regulatory T cells are known to increase the survival of circulating tumor cells, which promotes the establishment of metastatic foci [[40]11]. Thus, understanding complex immune responses between tumors and the immune system is crucial to developing new diagnostic biomarkers and anti-cancer immunotherapy. Diverse environmental factors including nutrition, air pollution, and toxicants are involved in carcinogenesis through varying epigenetic alterations of DNA such as DNA methylation, histone deacetylation, and non-coding RNA regulation [[41]12]. In particular, DNA methylation is highlighted due to its nature of early response to cancer [[42]13]. The methylation profiles of the promoter regions of tumor suppressor genes such as PTEN and BRCA2 were examined in different types of BC, and many current studies suggest that a few specific cytosine–phosphate–guanine (CpG) sites might be sufficient to provide prognostic information on BC with high accuracy and specificity [[43]14,[44]15]. Moreover, the Encyclopedia of DNA element (ENCODE) project revealed that the regulatory regions of genes such as promoters and enhancers are more conserved between human and dog models than human and rodent models. Indeed, a deeper understanding of methylation aberrations on gene regulatory regions is essential. MicroRNA (miRNA), of which the length is short (20–24 nucleotides long), is a well-studied non-coding RNA. Functionally, miRNAs mediate gene silencing via recognition of the 3′ untranslated region (UTR) of mRNAs, followed by guiding Argonaute proteins to these sites [[45]16]. MicroRNAs can repress multiple mRNAs which harbor complementary sites through two post-transcriptional mechanisms, mRNA cleavage and translational repression. Almost 2600 mature human miRNAs are currently identified and many of them are conserved in vertebrates and invertebrates [[46]17]. Moreover, numerous studies illustrated the function of oncogenic or tumor suppressive miRNA [[47]18]. For example, highly expressed miR-34a promotes apoptosis [[48]19], whereas miR-27a inhibits the expression of ZERB10/RINZF, which induces Sp factor expression followed by an increase in antipoetic and angiogenetic molecules [[49]20]. Methylation is known as an important epigenetic regulator of miRNA expression in cancer [[50]21], but miRNA studies across species are still veiled. The current pilot study profiled methylation landscapes using reduced representative bisulfite sequencing (RRBS) in the mammary tumor and adjacent normal tissues, as well as peripheral blood mononuclear cells (PBMCs), of dogs. These comparative methylome data provide a biological link to CMT-specific differential methylation regions (DMRs) in both tissues and PBMCs. We suggest that this approach using the datasets of CMTs and matching normal tissues with PBMCs will provide new candidates for CMT biomarkers and help form a better understanding of human BC. 2. Results 2.1. RRBS Demonstrated CMT-Related Methylation Profiles in Tissues and Blood Cells The sets of specimens, CMTs, adjacent normal tissues, and PBMCs, were collected from two individual patient dogs, and normal blood samples were drawn from dogs of matching breed in healthy condition. RRBS was performed on two sets of specimens obtained from dogs diagnosed as a carcinoma, complex type of CMT, for which degree of malignancy is low (grade I), to discover the cancer-specific methylation signature, and the result was validated by bisulfite conversion PCR in extended samples ([51]Figure 1 and [52]Table 1). Since a very limited number of samples were analyzed in the present study, we only focused on the carcinoma, complex type subtype and the Schnauzer dog breed in order to minimize variations. Figure 1. [53]Figure 1 [54]Open in a new tab Scheme of research flow and specimens. DMR: differentially methylated region, PBMC: peripheral blood mononuclear cell, RRBS: reduced representation bisulfite sequencing. Table 1. Patient information. Breed Sex Age Diagnosis PBMC Tissue Schnauzer FS 13 Carcinoma complex type O O Schnauzer FS 9 Carcinoma complex type O O Schnauzer MN 10 Healthy O O Schnauzer FS 10 Healthy O O Maltese F 12 Carcinoma complex type - O Dachshund FS 14 Carcinoma complex type - O Dachshund F 3 Carcinoma complex type - O [55]Open in a new tab Clinical characteristics of complex carcinoma and normal canine samples used in this study. Tissues are composed of canine mammary tumors (CMTs) and adjacent normal tissues. PBMC: peripheral blood mononuclear cell; FS: female spayed; MN: male neutralized; F: female. An averaged total of 2.9 gigabytes of data were sequenced, which covered 6.6–12.3% of canine genomic regions. Sequence statistics are summarized in [56]Table S1 and the reads’ quality check is visualized in [57]Figure S1. Assayed CpGs were mostly annotated in the intergenic and intron regions of the genome (97.9%), while ~2.2% were annotated in promoter and exon regions. The percentage methylation histogram was bimodal at both ends and the CpG base histogram showed typical shape of good quality in RRBS without PCR duplication ([58]Figure S1A,B). A CMT methylation signature was defined from a CMT-related DMR as a genomic region with a ≥15% change in methylation relative to matched normal with false discovery rate (FDR)-adjusted q-value ≤0.05. CMT methylation signatures were profiled in the comparison between CMTs and adjacent normal tissues and between PBMCs obtained from CMT patient and healthy dogs. A total of 41,044 and 27,667 DMRs were identified from the comparison of tissues and PBMCs, respectively. As known, genomic DNA in cancer cells when compared to normal cells is hypomethylated in human; in concordance, hypomethylation (26,562 in tissues and 17,631 in PBMCs) was found almost twice as often hypermethylation (14,482 in tissues and 10,036 in PBMCs) in CMT-related DMRs. Since the two healthy dogs consisted of both male and female sex, we computed CMT-related DMRs excluding DMRs on the X chromosome ([59]Table S2). Hierarchical clustering and heatmap analysis using the CMT-related DMRs successfully differentiated PBMCs in healthy condition from those in CMT conditions, as well as CMT tissues from adjacent normal tissues ([60]Figure 2A,B). Figure 2. [61]Figure 2 [62]Open in a new tab Genome-wide DNA methylation changes in CMT tissue and PBMCs. (A,B) Hierarchical clustering and heatmap of DMRs defined by difference ≥15%, false discovery rate (FDR)-adjusted q-value ≤ 0.05, and read coverage ≥10 in CMT tissues. The color indicates the percentage methylation: hyper (in red) and hypo (in blue). N: normal, C: CMT. Heatmaps for tissues (A) and PBMCs (B) were generated by complete linkage clustering using a Pearson correlation. (C,D) Volcano plots of all DMRs separated according to their genomic features. The x-axis represents methylation percentage and the y-axis represents −log2 (q-value). Methylation changes in promoter, exon, intron, and intergenic regions were examined; DMRs are depicted as blue (hypo) or red (hyper) dots in CMT tissues (C) and PBMCs (D). CMT: canine mammary tumor; PBMC: peripheral blood mononuclear cells; DMR: differentially methylated region. The list of CMT-related DMRs for which methylation level was significantly changed in CMTs only is summarized in [63]Table S2 (Supplementary Materials). It will be useful to investigate epigenetic alterations noticed from cancer liquid biopsies using cell-free DNA (cfDNA). To better understand the methylation landscape of CMT, all DMRs were separated according to genomic distribution in both gene and intergenic regions. We identified that DMRs were 35–38% hypermethylated and 62–65% hypomethylated in tissues, while 27–38% were hypermethylated and 62–73% hypomethylated in PBMCs. The genomic distribution of DMRs is depicted by volcano plot in [64]Figure 2C,D, and the total number of DMRs is summarized in [65]Table S3. The largest number of DMRs was observed in intergenic regions, followed by intron and promoter regions. These data were consistent with a previous study [[66]22]. The number of DMRs found in tumor tissues was higher than in PBMCs. However, the percentage distribution of hyper- and hypomethylated DMRs was similar between tissues and PBMCs (35% and 36% were hypermethylated in CMT tissues and PBMCs, respectively). Unexpectedly, more than half of the assayed region was obtained from the intergenic region after restriction enzyme digestion enrichment for CpG island (CGI)regions. This result could be explained by reports of canine intergenic CGI density being higher than that in the human genome [[67]23]. We also compared CG content in the intergenic DMR regions determined by chromosomal average. CG content of assayed intergenic DMRs was ~60%, which is higher than the chromosomal average (~40%), and means that CGI enrichment was successful ([68]Figure S2). On the other hand, more than half of the promoter DMRs were primarily hypomethylated in CMTs when compared to matched normal. Additionally, intron regions were also significantly enriched in DMRs. This may suggest that methylation dynamically regulated gene expressions at not only the promoter region but also the intron region. 2.2. Differentially Regulated DNA Methylations in CMT Tissues and PBMCs DMGs are defined as genes existing within 2 kb up- and downstream of a DMR. A total of 341 DMGs in tissues (202 hypo- and 139 hypermethylated) and PBMCs (175 hypo- and 72 hypermethylated) were subjected to gene ontology (GO) and Kyoto encyclopedia of genes and genomes (KEGG) pathway analysis ([69]Figure 3, [70]Table S4). Figure 3. [71]Figure 3 [72]Open in a new tab Enrichment analysis of gene ontology (GO) and KEGG pathway using DMGs in CMT. Gene ontology (GO) and KEGG pathway analysis in tissue and PBMC DMGs are depicted with −log (p-value), where the y-axis represents terms and IDs. From tissues, 139 hyper- and 202 hypomethylated DMGs, and, from PBMCs, 72 hyper- and 175 hypomethylated DMGs were applied. The top three significant GO_BP terms are depicted. Colors represent hyper (red) and hypo (blue) methylations. GO analysis in biological processes in tissues (A) and in PBMCs (B). KEGG pathway analysis in tissues (C) and PBMCs (D). KEGG: Kyoto encyclopedia of genes and genomes; DMG: differentially methylated gene; BP: biological process; CMT: canine mammary tumor; differentially methylated genes. Cancer-associated terms in GO enrichment analysis are illustrated in [73]Figure 3A,B. In tissues, DMGs were significantly enriched in biological processes (BP) of extrinsic apoptotic signaling pathway via death domain receptors (GO:0008625), positive regulation of protein kinase B signaling (GO:0051897), and positive regulation of cell proliferation (GO:0008284) ([74]Figure 3A). This is noticeable because a different study reported that protein kinase B signaling, including IGF2, TNF, and ZP3, has activity that inhibits apoptosis in human BC via acting growth factor transducer [[75]24]. The terms of hypomethylated DMGs enriched in positive regulation of cell proliferation and KRT and MMP families that involved tumor metastasis promoting through extra cellular matrix (ECM) remodeling [[76]25] were also hypomethylated ([77]Table S4, Supplementary Materials). Hypermethylated genes were enriched in positive regulation of cytosolic Ca^2+ concentration (GO:0007204), negative regulation of K^+ transmembrane transport (GO:1901380), and positive regulation of apoptotic process (GO:0043065) ([78]Figure 3A). Aberrant methylation of K^+ transmembrane transport-related genes, a frequently observed feature in human cancer, was also hypermethylated in CMTs [[79]26]. CAV1/3 and KCNH2 were enriched in this term; the CAV1 gene is thought to function as a tumor suppressor or modifier gene in mammary epithelia [[80]27]. In PBMCs, hypomethylated genes were highly clustered in positive regulation of cytosolic Ca^2+ concentration (GO:0007204), cell–cell signaling (GO:0007267), and macrophage chemotaxis (GO:0048246). Hypermethylated genes were enriched in negative regulation of apoptotic process (GO:0043066) and cell growth (GO:0016049) ([81]Figure 3B). It is an interesting feature that positive regulation of cytosolic Ca^2+ concentration GO was hypomethylated in cancer PBMCs while it was hypermethylated in cancer tissues. Moreover, through KEGG pathway analyses, we determined that the immune response-related pathways were hypomethylated in both CMT tissues and PBMCs ([82]Figure 3C,D). Interestingly, cytokine–cytokine receptor interaction was hypomethylated and was likely to be enriched in CMT tissues, indicating that cytokine–cytokine receptor interaction genes are demanded for crosstalk in immune cells. On the other hand, chemokine signaling was also enriched in PBMCs ([83]Figure 3C). Aberrant regulation of chemokine signaling via DNA methylation dictated that many signaling pathways in immune cells are affected in the tumor environment. Although many studies were conducted on cancer immune systems, there is limited knowledge on the methylation status of specific gene sets directly associated with cancer cells. 2.3. Cancer-Associated miRNA Expressions Regulated by Methylation in CMT One interesting finding in the KEGG pathway enrichment analysis was that the term of microRNAs in cancer (KEEG:05206) was highly enriched in both hypo- and hypermethylated DMGs in CMT tissues, and hypermethylated DMGs in PBMCs ([84]Figure 3C,D). Importantly, we found that the lists of miRNAs enriched in each hypo- or hypermethylated DMG were clearly separated according to their roles in cancer, such as oncogenic and tumor suppressive miRNAs. Many of the identified miRNAs, including miR-141, miR-142, and miR-10B, which are all known as cancer-associated miRNAs, were enriched in hypomethylated DMGs. On the other hand, miR-22, miR-149, and miR-124-2, which are all reported with cancer suppressive activity, were found in hypermethylated DMGs [[85]28,[86]29,[87]30,[88]31]. These results suggest that cancer-enriched methylation profiles may directly link to cancer development via regulating cancer-associated miRNA expression ([89]Table 2). Table 2. Potential functions of different microRNAs (miRNAs) in cancer progression. Name Methylation Status in Canine Function in Human References