Association between genome-wide epigenetic and genetic alterations in breast cancer tissue and response to HER2-targeted therapies in HER2-positive breast cancer patients: new findings and a systematic review

Recent evidence suggests that genetic and epigenetic mechanisms might be associated with acquired resistance to cancer therapies. The aim of this study was to assess the association of genome-wide genetic and epigenetic alterations with the response to anti-HER2 agents in HER2-positive breast cancer patients. PubMed was screened for articles published until March 2021 on observational studies investigating the association of genome-wide genetic and epigenetic alterations, measured in breast cancer tissues or blood, with the response to targeted treatment in HER2-positive breast cancer patients. Sixteen studies were included in the review along with ours, in which we compared the genome-wide DNA methylation pattern in breast tumor tissues of patients who acquired resistance to treatment (case group, n = 6) to that of patients who did not develop resistance (control group, n = 6). Among genes identified as differentially methylated between the breast cancer tissue of cases and controls, one of them, PRKACA, was also reported as differentially expressed in two studies included in the review. Although included studies were heterogeneous in terms of methodology and study population, our review suggests that genes of the PI3K pathway may play an important role in developing resistance to anti-HER2 agents in breast cancer patients. Genome-wide genetic and epigenetic alterations measured in breast cancer tissue or blood might be promising markers of resistance to anti-HER2 agents in HER2-positive breast cancer patients. Further studies are needed to confirm these data.


INTRODUCTION
The human epidermal growth factor receptor 2 (HER2) is a transmembrane tyrosine kinase receptor and belongs to the epidermal growth factor receptor (EGFR) family [1] . It comprises an extracellular domain (ECD), a transmembrane segment, and an intracellular region [2] . HER2 gene amplification and receptor overexpression, which occur in approximately 15%-20% of breast cancer patients, are important markers for poor prognosis, including a more aggressive disease and shorter survival [3] . In addition, HER2-positive status is considered a predictive marker of response to HER2-targeted drugs [4] . Detection of receptor overexpression via immunohistochemistry (IHC) and/or HER2 gene amplification using in situ hybridization (ISH) techniques in breast cancer tissue, including fluorescent ISH (FISH), determines patients' eligibility to receive anti-HER2 therapies [5] . Food and Drug Administration (FDA)-approved anti-HER2 agents currently used in clinical settings in combination with chemotherapy comprises recombinant monoclonal antibodies that bind the ECD of HER2, such as trastuzumab (Herceptin®), pertuzumab (Perjeta®) and trastuzumab emtansine (T-DM1, Kadcyla®), and small molecule tyrosine kinase inhibitors, like lapatinib (TYKERB®), that inhibit enzyme function of the intracellular catalytic domain of HER2 and other EGFR members [6,7] .
Although targeted treatment with anti-HER2 agents has significantly improved the disease-free and overall survival rates of metastatic and early-stage HER2-positive breast cancer patients [8][9][10][11] , resistance to anti-HER2 therapy, both primary and acquired, has emerged as a major clinical problem in the treatment of HER2positive breast cancer patients [12][13][14] . Even though several molecular mechanisms of resistance to anti-HER2 agents have been proposed in preclinical models, no clinically applicable strategy to overcome resistance to these targeted treatments has been identified yet [15] . Therefore, there is an urgent need to identify reliable predictive molecular markers of treatment failure with the ultimate goal of developing targeted drugs that can overcome resistance.
The aim of this pilot study was to analyze the association between DNA methylation patterns in breast cancer specimens and response to trastuzumab in a cohort of 12 trastuzumab-treated, non-metastatic HER2-positive breast cancer patients. Additionally, a systematic review that combines these findings with all available published results on the association of genome-wide genetic and epigenetic alterations in breast cancer tissue or blood with response to anti-HER2 treatment in HER2-positive breast cancer patients is also reported.

Study population and data collection
The study population consisted of 12 women (six cases and six controls) selected among 106 trastuzumabtreated patients with non-metastatic, HER2-positive breast cancer diagnosed between July 1, 2005 and December 31, 2010 at the Centre des Maladies du Sein, a specialized breast center in Quebec City, Canada. Information on tumor characteristics and prognostic factors at the time of diagnosis (baseline) and followup information were collected from medical records. The clinical endpoint in this study was disease-free survival (DFS). All breast cancer recurrences (locoregional, contralateral breast, and distant) were considered as events, whereas death (from any cause) before recurrence and loss to follow-up were considered as censoring events.
Over a mean follow-up period of 6.22 years, 22 patients out of 106 experienced recurrence. Eight cases were randomly selected among all patients who developed recurrence during follow-up and six had a sufficient amount of primary breast cancer tissue available for DNA extraction (see below). Of note, baseline characteristics of the six selected cases were comparable to the total population of cases (n = 22) for all characteristics except tumor grade (the proportion of grade III tumors among the selected cases was 67% vs. 41% for the total population of cases) [Supplementary Table 1]. For each case, one control was selected from the 84 patients who had not developed recurrence and were alive at the date of the case's recurrence. Controls were matched to cases for the following factors: age at diagnosis (with 5-year age categories), estrogen receptor (ER) status, year of diagnosis (with 2-year categories), and menopausal status. The number of samples used was determined upon the availability of samples and not evaluated using a statistical sample size calculation. All patients provided written informed consent. Ethical approval of the study was obtained from the Research Ethics Committee of the Centre de Recherche du CHU de Québec (# 2016-2802).

Gene methylation assessment
To ensure that DNA methylation was analyzed to the greatest possible extent in breast cancer tissue and to reduce contamination with other cell types (lymphocytes, adipocytes, fibroblasts), tissue microarray (TMA) blocks containing formalin-fixed, paraffin-embedded (FFPE) breast cancer tissue cores (1 mm in diameter) were constructed for each patient, as previously described [37] . From each TMA block, one section was stained with hematoxylin and eosin (H&E) to verify the cellular composition of the cores. Cores were removed from TMA blocks if they contained abnormal tissue or if epithelial tumor tissue occupied < 70% of the core area before proceeding to DNA extraction. H&E sections were prepared from different levels of the TMA blocks: at the beginning, at regular intervals (every tenth 10-µm-thick serial section), and after the last section. DNA was extracted from tissue cores using GeneJET FFPE DNA Purification kit (ThermoScientific, Ottawa, Canada) with minor modifications to the manufacturer's instructions in which samples were incubated with Digestion buffer for six minutes and incubated with Proteinase K solution for 180 minutes.
DNA samples were sent to Génome Québec Innovation Center (Montreal, Canada). Methylation was measured with the Illumina HumanMethylation450 BeadChip array (Illumina Inc., San Diego, CA, USA) following the manufacturer's instructions for the bisulfite treatment, Infinium FFPE quality control (Illumina FFPE QC kit, Illumina, Inc., CA, USA), and DNA restoration. This BeadChip interrogates 482,421 CpG sites, 3091 non-CpG sites, 65 random SNPs, and covers 21,231 RefSeq genes. It uses two distinct oligonucleotide probes (Infinium I and Infinium II) to assess methylation levels [38] .

Statistical analysis
Raw β-values, defined as the ratio of the methylated probe intensity to the overall intensity (sum of the methylated and unmethylated probe intensities) [39] , were imported into the R statistical programming environment (version 3.2.2). Since M-values [logit transformed β-values, calculated as ] are considered more reliable in the detection rate and true positive rate for both highly methylated and unmethylated CpG sites compared to β-values [39] , M-values were used for statistical analyses.
Quality control was performed with the qcReport function from the minfi package, and none of the 12 samples were excluded due to bad quality control. The Dasen method from the WateRmelon package, also known as data-driven separate normalization, was used to background correct and quantile normalize data based on methylated and unmethylated intensities, separately, by probe types (Infinium I and II) [40] . A probe filtration step was performed to remove CpG sites corresponding to probes that could affect our analysis, including probes with bad detection (detection P-value > 0.01); unique probes having a common single nucleotide polymorphism (SNP) in European individuals at the interrogated CpG loci or the single-base extension according to the list published by Chen et al. [41] ; probes that can hybridize to multiple loci also listed by Chen et al. [41] ; and probes located on X and Y chromosomes. A total of 76,161 unique probes were removed, leaving 406,260 autosomal probes for the analysis. Data were verified for confounding batch effects due to separate chips [42] , and none were observed. All samples passed quality-control tests and were therefore retained in the analysis.
Baseline characteristics between cases and controls were compared using Fisher's exact test for categorical variables, Student's t-test for follow-up time and Wilcoxon-Mann-Whitney test for the other continuous variables. The difference in global methylation levels between median M-values of cases and controls was assessed using a Wilcoxon signed-rank test for paired samples. Differentially methylated probes (DMPs) were identified using LIMMA (robust linear regression method), taking into account the matching factors between cases and controls (i.e., age at diagnosis, ER status, year of diagnosis, and menopausal status). Multiple testing correction was performed using false discovery rate (FDR) estimation (cut-off < 5%). In addition, we used a log 2 -fold change |log2FC| (i.e., the difference between mean M-values measured in breast cancer tissues of resistant patients and controls) > 2.0 as a cut-off to identify probes that were strongly differentially methylated between cases and controls.

Systematic review of published findings
A systematic review was conducted and reported according to the 2020 Preferred Reporting Items for Systematic Reviews and Meta-Analyses (PRISMA) guidelines [43] . Exposure: To be included, a study had to measure response-specific survival using including pathologic complete response (pCR), disease-free survival (DFS), progression-free survival (PFS), or event-free survival (EFS).

Eligibility criteria
Outcome: We considered all assessments of genetic and epigenetic alterations at the genome-wide level in breast cancer tissue or blood, whatever the measurement method.
Types of studies: Any observational or randomized controlled study that assessed the association between genome-wide genetic and epigenetic alterations in breast cancer specimens or blood and the response to anti-HER2 agents in HER2-positive breast cancer patients was eligible. Case reports were excluded.
Only full available articles in English were included.

Information sources
Dragic D searched the PubMed biomedical database from inception to the search date of March 25, 2021, to identify eligible studies.

Search strategy
The search strategy was developed by Dragic D and Furrer D, approved by Diorio C, using controlled vocabulary search terms and free-text words related to HER2-positive breast cancer, genome-wide genetic and epigenetic alterations and treatment outcome [Supplementary Table 2]. No restrictions regarding language were applied.

Selection process
The references identified by the search strategy were selected according to the predefined eligibility criteria in a two-step process: titles and abstracts were screened by one author (Dragic D), and full texts of retained articles were examined by two authors (Dragic D and Furrer D). Disagreements between the two authors were discussed until a consensus was reached, and whenever required, a third review author (Diorio C) was consulted.

Data collection process
We designed a data extraction form for this review. Data were extracted by two authors (Furrer D and Dragic D) for half of the included studies. A third author (Diorio C) was consulted when discrepancies between both authors could not be resolved. For the remaining studies, the extraction was done by one author (Dragic D), and when needed, a second author (Diorio C) was involved. Additionally, the authors of studies of interest that lacked data to evaluate eligibility (n = 1) or other measures needed for this review (n = 3) were contacted to obtain the necessary information.

Data items
For all selected articles, study characteristics (study design and sample size), patient's characteristics (age, stage, ER and PR status, menopausal status, and treatment received), assessment of genetic and epigenetic alterations (tissue processing, DNA, RNA or miRNA extraction method, assessment method, and parameters used), as well as statistical methods and study results, were collected. The study's definition of response to targeted treatment was recorded.

Risk of bias assessment
Studies included in the review were assessed for risk of bias using the "Risk Of Bias in Non-randomized Studies of Interventions (ROBINS-I) tool [44] . The following domains were assessed: selection of participants in the study, exposure measurement, outcome measurement, potential confounding factors, missing data, and selective reporting.
Assessment of the risk of bias was done by two authors (Furrer D and Dragic D) for half of the included studies. Inconsistencies were discussed to reach a consensus. For the remaining studies, the assessment was done by one author (Dragic D), and when required, a second reviewer (Furrer D) was consulted.

Assessment of heterogeneity
Differences between studies, including study design, patient characteristics (age, menopausal status, ethnicity, and treatment received), tumor characteristics (stage, ER status), assessment of genome-wide genetic and epigenetic alterations (tissue processing, extraction method, and measurement method), and different levels of risk of bias were considered for exploring possible sources of heterogeneity.

Synthesis methods
Considering that high heterogeneity between studies was expected, quantitative data synthesis was not considered appropriate. Instead, we adopted a qualitative systematic review approach to investigate the relationship between epigenetic and genetic alterations and response to HER2-targeted therapies in HER2positive breast cancer patients. The selection process was detailed using the PRISMA 2020 flow diagram. Extracted data were first reported in a table gathering summary characteristics of all included studies. Additional information specific to the epigenetic or genetic method used was detailed in several tables. Results and genes identified by several studies were also highlighted in a table. If pathway analysis was not presented, we performed pathway analysis using the list of the differentially expressed genes reported by the study authors and the PANTHER online software (Protein Analysis Through Evolutionary Relationships). P-values < 0.05 were considered significant.

Pilot study of new findings
Genome-wide DNA methylation data in 12 breast cancer specimens were obtained from six trastuzumabtreated HER2-positive breast cancer patients who experienced recurrence during follow-up (cases) and six individually matched patients who had not developed recurrence and were alive at the date of the case's recurrence (controls). Baseline characteristics of cases and controls are summarized in Table 1. Baseline characteristics for both groups were comparable in clinicopathological characteristics (tumor grade, lymph node status, and tumor size) and treatment received. Compared to controls, a higher proportion of cases (50%) had a body mass index > 25 kg/m 2 , although not statistically significant.
Global methylation levels between cases and controls were not statistically different: the median M-values of cases were 0.487, and the median M-values of controls were 0.504 (P-value: 0.844). At probe methylation levels, we identified 2,009 CpGs (1,382 genes) that were differentially methylated between cases and controls: 1,200 DMPs (885 genes) were significantly hypermethylated and 809 DMPs (497 genes) were significantly hypomethylated in tumor tissues of cases compared to those of controls after multiple testing correction (FDR < 0.05).

Study selection
Of the 758 references retrieved by electronic search in PubMed, we reviewed 52 full-text documents, and fifteen met the eligibility criteria [ Figure 1]. This review also included our pilot study that assessed the association between DNA methylation patterns in breast cancer specimens and response to trastuzumab in a cohort of 12 HER2-positive breast cancer patients.
DFS was reported in four studies including ours [46,52,54] , RFS in two studies [49,53] , and pCR in ten studies [45,47,48,50,51,[55][56][57][58][59] . DFS was defined as the time between the diagnosis and recurrence (locoregional recurrence, recurrence in the contralateral breast, and distant breast cancer recurrence) or death; RFS was defined as the time from the start of trastuzumab treatment to the first local, regional or distant recurrence event; and pCR as the absence of invasive breast cancer in the breast and axillary lymph nodes at the time of surgery. Besides ours, the other fifteen included studies were published between 2007 and 2020 and involved between 9 and 849 HER2-positive, anti-HER2 therapy-treated breast cancer patients [ Table 3].
Studies of genome-wide DNA methylation and response to targeted treatment Characteristics of the study (our new findings) that evaluated the association between genome-wide DNA methylation patterns measured in breast cancer tissues and response to targeted treatment in HER2positive, trastuzumab-treated breast cancer patients are presented in Supplementary Table 3. The study was retrospective. The mean age of included HER2-positive breast cancer patients was 51.3 years. The proportion of ER-positive breast cancer patients was 67%, and the proportion of positive lymph node status was 92%. Sixty-seven percent of HER2-positive breast cancer patients were postmenopausal. DNA methylation was assessed using Infinium HumanMethylation450 BeadChip array. Breast cancer patients were treated with adjuvant chemotherapy and trastuzumab. The study reported DFS. HER2-positive breast cancer patients were of Caucasian ethnicity. Tumor cell fraction was ≥ 70%.

Study of genome-wide lincRNA expression profile and response to targeted treatment
Characteristics of the single study that examined the association between genome-wide lincRNA expression profile and response to targeted treatment in HER2-positive trastuzumab-treated breast cancer patients are presented in Supplementary Table 6. The study had a retrospective design. lincRNA expression was measured genome-wide in fresh frozen breast cancer tissues using the Illumina HiSeq 2500 platform in a cohort of 13 HER2-positive breast cancer patients. Patients were treated with trastuzumab and pCR was reported. Tumor cell fraction was not reported.

Study of genome-wide CNA profile and response to targeted treatment
Characteristics of the two studies that examined the association between genome-wide CNA profile and response to HER2-targeted therapies in HER2-positive breast cancer patients are presented in Supplementary Table 7. The studies had a prospective design [48,56] . CNA profile was assessed genome-wide in fresh frozen breast cancer tissues using Affymetrix Genome-wide Human SNP array [48] or Illumina HiSeq [56] in cohorts of 68 and 11 HER2-positive breast cancer patients. Patients were treated with trastuzumab, lapatinib or both in one study [48] and trastuzumab in the other [56] . pCR was reported. Tumor cell fraction was not reported.

Study of genome-wide protein expression profile and response to targeted treatment
Characteristics of the only study that examined the association between genome-wide protein expression profile and response to targeted treatment in HER2-positive trastuzumab-treated breast cancer patients are presented in Supplementary Table 8. The study had a retrospective design [57] . The protein expression profile was assessed genome-widely in blood using TMT-6 plex Isobaric Label Reagent Set in a cohort of 6 HER2positive breast cancer patients. Patients were treated with trastuzumab. pCR was reported.

Study of genome-wide mutation profile and response to targeted treatment
Characteristics of the two studies that examined the association between genome-wide somatic or germline mutation profile and response to targeted treatment in HER2-positive breast cancer patients are presented in Supplementary Table 9. The studies had a prospective design [58,59] . Somatic mutation profile was assessed genome-wide in fresh frozen breast cancer tissues for the two studies. Germline mutation profile was assessed genome-widely in blood using Illumina HiSeq 2000 [59] , and whole genome sequencing (WGS) and whole exome sequencing (WES) + Illumina HiSeq 2000 platform [58] in cohorts of 48 and 203 HER2-positive breast cancer patients. Patients were treated with trastuzumab. pCR was reported. Tumor cell fraction was not reported in one study [58] and was at least 10% in the other study [59] .
Gámez-Pozo et al. and Sorokin et al. identified several pathways associated with response to trastuzumab, including those involved in EGF receptor signaling, PI3K, apoptosis signaling, and p53 [47,52] . Gámez-Pozo et al. observed that the PI3K pathway was the most strongly associated with treatment response [47] [Supplementary  [49,52] . For two [46,48] out of the four included studies that did not report pathway analysis [45,46,48,50] , we performed gene ontology analysis using the list of genes reported as being differentially expressed by the authors using PANTHER. We observed that the Notch signaling pathway was overrepresented in both studies [46,48] , an observation also reported by Sorokin et al. [52] . We observed that Wnt signaling was overexpressed in Khoury et al. [46] , as reported by Sorokin et al. [52] . We did not perform pathway analysis for the two remaining studies [45,53] , as the number of differentially expressed genes (n = 11) was too small to perform the analysis [45] or not reported [53] .
In our study, we observed overlap between the identified differentially expressed genes and the strongly differentially methylated (i.e., |log2FC| > 2.0) genes. PRKACA was hypermethylated in our study (within the TSS region as well as the gene body), upregulated in one study [50] , and downregulated in another study [47] .
We observed that among the 344 genes that were reported to be differentially expressed (at the mRNA level) between breast cancer tissues of cases and controls in at least two studies, 170 were predicted to be targeted by miRNA identified as differentially expressed between breast cancer tissues of cases and controls in Ohzawa et al. [55] or Du et al. [54] [Supplementary Table 14]. We observed that among the 15 genes identified as differentially methylated in our study, five (KCNH7, ADAMTS2, SIX2, DOCK1 and ZNF598) were predicted to be targeted by miRNA identified as differentially expressed in breast cancer tissues of cases compared to controls in the study of Ohzawa et al. [55] or Du et al. [54] [Supplementary Table 14].
Genome-wide lincRNA expression profile in breast cancer tissues and response to targeted treatment in HER2positive breast cancer patients treated with trastuzumab Merry and collaborators observed that 371 lincRNAs were differentially expressed in non-pCR samples compared to pCR samples, where 33 lincRNAs showed decreased expression and 338 increased expression [50] [Supplementary Among the differentially expressed lincRNAs reported by Merry et al. [50] , we observed 44 genes that were predicted to be target genes of miRNA identified as differentially expressed in the study of Ohzawa et al. [55] or Du et al. [54] [Supplementary Table 14].
None of the genes whose lincRNAs were reported to be differentially expressed between patients showing pCR and those showing non-pCR in Merry et al. was reported as differentially methylated with a |log2FC| > 2.0 between cases and controls in our study [50] . However, when we consider the entire list of differentially methylated genes in our study, regardless of log 2 -fold change, we identified six overlapping genes. Among these six genes, five were hypermethylated (GABRA5, ZIC5, GRAMD4, RSPH3, and VCAN), and one was hypomethylated (CSMD1) in breast cancer samples of cases compared to controls.
We observed that among the 557 genes that were reported to show CN alteration in non-pCR samples compared to pCR-samples in the study of Guarneri et al. [48] , 279 were predicted to be targeted by miRNA identified as differentially expressed in breast cancer tissues of cases compared to controls in the study of Du et al. [54] or Ohzawa et al. [55] [Supplementary Table 14].
No overlap was observed between genes showing CN alterations in non-pCR samples compared to pCRsamples in Guarneri et al. and genes identified as differentially methylated with a |log2FC| > 2.0 in our study [48] . However, when we consider the whole list of differentially methylated genes, we identified 27 genes that were differentially methylated in breast cancer tissues of cases compared to controls in our study among the 557 genes showing CNV variations in the study of Guarneri et al. [48] .

Genome-wide protein expression profile in blood of breast cancer cases compared to controls
Out of the 18 genes that were reported as differentially expressed (five downregulated and 13 upregulated) in the blood of breast cancer cases compared to controls by Yang et al. [57] [Supplementary Table 17], three (APOB, SLC3A2, CST3) were differentially expressed in breast cancer tissues of cases compared to control in at least two studies. Three (LDHA, DBF4B, and MASP1) were predicted to be targeted by miRNA identified as differentially expressed in breast cancer tissues of cases compared to controls in the study of Du and collaborators [54] or Ohzawa and collaborators [55] [Supplementary Table 14].
None of the differentially expressed genes in the study of Yang et al. were differentially methylated in our study [57] .

Genome-wide somatic and germline mutations profile in breast cancer tissues of cases compared to controls
Whereas Shi and collaborators observed that higher somatic mutation frequency in the PIK3CA gene in the breast tissues of cases was associated with trastuzumab resistance [59] , Lesurf et al. reported that no somatic or germline mutations were associated with response to trastuzumab in breast cancer tissues of cases compared to controls [58] [Supplementary Table 18].
The guanine nucleotide exchange factor DOCK1 (Dedicator of cytokinesis protein 1) is involved in cytoskeletal rearrangements required for phagocytosis of apoptotic cells and cell mobility [70] . A recent study reported that DOCK1 inhibition leads to suppressed migration of the triple-negative breast cancer cell lines MDA-MB-157 and MDA-MB-231 [68,71] . ADMATS2 (ADAM metallopeptidase with thrombospondin type 1 motif 2) belongs to the ADAM metallopeptidase with thrombospondin type 1 motif and processes collagen precursors into mature collagen molecules [72] . It has been proposed that ADAMTS2 exerts an anti-tumor effect by inhibiting intratumoral vascularization [73] . The PLEC gene encodes the pectin protein, which plays a role in maintaining tissue integrity [74] . A recent study suggests that a PLEC gene polymorphism (rs138924815) might increase the risk for familial testicular cancer [75] . USP4 (Ubiquitin carboxyl-terminal hydrolase 4) is a deubiquitinating enzyme that removes conjugated ubiquitin from target proteins [76] . It has been reported that USP4 expression was decreased in breast cancer tissue samples compared to paired normal breast tissues [63] . Moreover, USP4 expression was associated with decreased proliferation in two HER2-negative breast cancer cell lines (MCF7 and BT549) [63] . The PRKACA gene encodes for a protein kinase that plays a role in controlling cellular processes such as glucose metabolism and cellular division [72] . Of note, one recent study suggests that PRKACA expression might be associated with the development of trastuzumab resistance in HER2-positive breast cancer patients [64] . The authors observed that in a subgroup of HER2-positive breast cancer patients who developed trastuzumab resistance (three out of five patients), PRKACA expression was highly increased in the breast cancer sample obtained after the onset of trastuzumab resistance compared to the pre-treatment sample. Considering that in our study, the PRKACA gene was hypermethylated within the gene body and that hypermethylation within this gene region often promotes gene elongation and, therefore, gene expression, we can postulate that our results might be concordant with those reported by Moody et al. [64] . Although in our study, DNA methylation was exclusively measured in pre-treatment samples and Moody et al. observed increased PRKACA expression only in breast cancer samples obtained after the onset of recurrence [64] .
Among all genes that we identified as strongly differentially methylated between trastuzumab-resistant and trastuzumab-sensitive HER2-positive breast cancer patients, one of them, PRKACA, was reported as being higher expressed in breast cancer tissues of cases compared to controls in the study conducted by Merry and collaborators [50] and as being lower expressed in breast cancer tissues of cases compared to controls in the study of Gámez-Pozo et al. [47] . In our study, the PRKACA gene was hypermethylated within the TSS region and the gene body. Our observation could be partly concordant with the results reported by Merry et al. [50] (but not with those of Gámez-Pozo et al. [47] ), as hypermethylation in the gene body (but not within the TSS region) is usually associated with increased gene expression [34] . Interestingly, the results of Merry and collaborators [50] (but not those of Gámez-Pozo et al. [47] ) are concordant with the study above [64] , where the authors observed that PRKACA expression was increased in breast cancer samples of HER2positive trastuzumab-resistant breast cancer patients.
In one of the studies retained in our systematic review, the authors created a predictive model to differentiate HER2-positive trastuzumab-treated breast cancer patients with a higher risk of relapse from those with a lower risk in a cohort of 53 patients [49] . The validity of this predictive model was then confirmed in an independent and bigger data set. The authors observed that differentially expressed genes in the breast cancer tissues of patients identified as at low risk for relapse in the independent data set using this model were associated with the immune system. When we performed gene enrichment analysis of genes showing CN alteration in the study of Guarneri and collaborators, we observed that pathways associated with the immune response (inflammation mediated by chemokine and cytokine signaling pathway) were overrepresented. The immune system's involvement in response to anti-HER2 agents in HER2-positive breast cancer patients has also been reported in other studies [77][78][79][80][81][82] .
To our knowledge, we conducted the first systematic review on the association of epigenetic and genetic alterations in breast cancer tissues or blood with the response to anti-HER2 agents in HER2-positive breast cancer patients. Sixteen studies were included in this review, and very few overlaps between studies were found. The most consistent results were the higher expression of GOLGA2 and PHF21A genes and the higher expression and CN gain of MIA3, WDR26 and C1orf133 genes observed in three different studies. Among these five genes, WDR26 has been shown to promote breast cancer growth and metastasis via the PI3K/AKT signaling pathway [83] . Gámez-Pozo and collaborators also observed that genes identified as differentially expressed in breast tumor tissues of cases compared to controls were overrepresented in the PI3K pathway [47] . Interestingly, the study by Shi and collaborators [59] also reported that PI3K mutations were associated with response to anti-HER2 agents in HER2-positive breast cancer patients and other previously published studies [84][85][86] . Taken together, these results suggest that genes of the PI3K pathway might play a relevant role in the development of resistance to anti-HER2 agents in breast cancer patients. Exploring the ramifications of these and other findings in larger cohorts or datasets like TCGA should be considered in future studies.
Several factors might explain why we only observed a few overlaps in our systematic review. Tumor cell content varied from > 40% to > 80% between studies, and in several publications (n = 4), this information was not provided. Contamination with cell types other than breast cancer cells can modify the observed pattern of genetic and epigenetic markers, as DNA methylation and other epigenetic or genetic markers widely vary across tissues and cellular types [87] . Moreover, the type of outcome evaluated in the study might also play a role, as the assessment of pCR in the neoadjuvant setting might mainly identify patients who did not primarily respond to targeted treatment (primary resistance), whereas the evaluation of DFS in the adjuvant setting might allow identifying patients who initially respond to targeted treatment but who develop resistance over time (acquired resistance). Moreover, epigenetic and genetic markers can be influenced by clinicopathological data, including age, stage, ER status, menopausal status, and ethnicity [88][89][90][91][92][93][94][95][96] .
Unfortunately, it was difficult to evaluate this aspect, as patients' clinicopathological data were not extensively reported in the majority of the included studies. The risk of bias in most studies was due to confounding.

CONCLUSION
In conclusion, although the sample size of the present pilot study was small and despite the lack of validation cohort, using a high-throughput analysis, we identified genes that were differentially methylated in breast cancer tissues of HER2-positive trastuzumab-treated breast cancer patients who developed resistance toward this drug compared to those who responded to targeted therapy. One of the most differentially methylated genes, PRKACA, has been reported to be differentially expressed in breast cancer tissues of trastuzumab-resistant compared to trastuzumab-sensitive HER2-positive breast cancer patients in two studies included in our systematic review. Although we identified very few genes that overlap between studies, our review suggests that some of the genes acting in the PI3K pathway, such as PRKACA [97] , might play an important role in developing resistance to anti-HER2 agents in breast cancer patients. Although the associations between PI3KCA mutations and PI3K dysfunctions and anti-HER2 treatment resistance are well documented in the literature [86,[98][99][100] , further studies on this topic are needed, which may help to unveil carcinogenic mechanisms involved in this pathway.
Although the observations reported in the retained studies were only marginally concordant, our work and the studies presented in this article suggest that knowledge gathered from these high-throughput studies could be useful for the identification of novel biomarkers of trastuzumab resistance. This might promote the development of new targeted drugs that could be administered to trastuzumab-resistant HER2-positive breast cancer patients.