Combined HER3-EGFR score in triple-negative breast cancer provides prognostic and predictive significance superior to individual biomarkers

Epidermal growth factor receptor (EGFR) and human epidermal growth factor receptor 3 (HER3) have been investigated as triple-negative breast cancer (TNBC) biomarkers. Reduced EGFR levels can be compensated by increases in HER3; thus, assaying EGFR and HER3 together may improve prognostic value. In a multi-institutional cohort of 510 TNBC patients, we analyzed the impact of HER3, EGFR, or combined HER3-EGFR protein expression in pre-treatment samples on breast cancer-specific and distant metastasis-free survival (BCSS and DMFS, respectively). A subset of 60 TNBC samples were RNA-sequenced using massive parallel sequencing. The combined HER3-EGFR score outperformed individual HER3 and EGFR scores, with high HER3-EGFR score independently predicting worse BCSS (Hazard Ratio [HR] = 2.30, p = 0.006) and DMFS (HR = 1.78, p = 0.041, respectively). TNBCs with high HER3-EGFR scores exhibited significantly suppressed ATM signaling and differential expression of a network predicted to be controlled by low TXN activity, resulting in activation of EGFR, PARP1, and caspases and inhibition of p53 and NFκB. Nuclear PARP1 protein levels were higher in HER3-EGFR-high TNBCs based on immunohistochemistry (p = 0.036). Assessing HER3 and EGFR protein expression in combination may identify which adjuvant chemotherapy-treated TNBC patients have a higher risk of treatment resistance and may benefit from a dual HER3-EGFR inhibitor and a PARP1 inhibitor.

risk-predictive biomarker in TNBC [2][3][4][5] . HER3 is the only receptor in the family that is catalytically inactive and requires dimerization with other members in order to be activated 6 . HER3 overexpression has been reported in approximately 20-30% of invasive breast carcinomas 7 . However, the value of HER3 and EGFR as prognostic biomarkers in TNBC remains uncertain. One study found that HER3-positivity was an independent predictor of poor disease-free survival in breast cancer after adjusting for breast cancer subtype and other potential confounders, whereas EGFR-positivity did not independently predict disease-free survival 8 . Reduced expression of EGFR members can be compensated by increased expression of HER3, which promotes resistance to EGFR inhibitors 8,9 . Since EGFR expression varies with that of HER3, these proteins may need to be considered jointly to be useful biomarkers. Thus, the aims of the present study were to (i) determine whether combined HER3-EGFR protein expression in pre-treatment resection specimens could independently predict survival in TNBC patients, (ii) compare the prognostic value of HER3-EGFR to individual HER3 and EGFR protein expression, and (iii) identify potential targets for therapy based on stratification by combined HER3-EGFR protein expression in RNA-sequenced pre-treatment TNBCs.

Materials and Methods
Ethics approval and consent to participate. Nottingham cohort: This study was approved by the Nottingham Research Ethics Committee 2 under the title "Development of a molecular genetic classification of breast cancer, " as well as the Institutional Review Boards (IRB) of the respective institutions. All samples from Nottingham used in this study were pseudo-anonymized and collected prior to 2006 and therefore under the Human Tissue Act informed patient consent was not needed. Release of data was also pseudo-anonymized as per Human Tissue Act regulations. Stavanger cohort: The use of the Stavanger material was approved by the Regional Ethical Committee (REC) (2010/1241). REC also approved that informed patient consent was not needed. All insights in a patient's journal were monitored electronically, and all, except the treating physician, are required to state the reason why they read a patient's journal. This log is always available for all patients. All data and biological material provided were pseudo-anonymized. Emory cohort All samples from the Emory cohort were archival tissue for which informed consent was not required as approved by the Emory IRB. All data and biological material provided were pseudo-anonymized.
Datasets, specimens, immunohistochemistry, and scoring. TNBCs were identified as lacking ER and PR expression and lacking HER2 over-expression and/or gene amplification (i.e., ER−/PR−/HER2−). Datasets were compiled from Nottingham City Hospital, Nottingham UK (n = 302); Stavanger University Hospital, Stavanger, Norway (n = 104); and Emory University Hospital, Atlanta, GA, US (n = 104), all consisting of TNBC patients with primary operable invasive disease. Biomarkers were evaluated prior to any treatment. For the Nottingham and Emory cohorts, specimens comprised tissue microarrays, whereas for the Stavanger cohort, specimens comprised full-face tissue sections. All slides were centrally reviewed.
Immunohistochemistry and scoring. Nottingham specimens were from patients belonging to the well-characterized Nottingham Tenovus Primary Breast Carcinoma Series 10 , for which HER3 and EGFR scores from tissue microarrays had been previously determined as described by Abd El-Rehim et al. 11 . Tissue microarrays (Emory cohort) or full-face formalin-fixed paraffin-embedded tissue slides (Stavanger cohort) were deparaffinized and then rehydrated in serial ethanol baths (100%, 90%, 75% and 50%). For HER3, antigen retrieval was performed using the Dako Target Retrieval Solution at pH 9 in a 90 °C water bath for 35 min. Then, slides were incubated with monoclonal mouse anti-human HER3 antibody (DAK-H3-IC) at 1:25 for 30 min. For EGFR, antigen retrieval was performed using Proteinase K (Agilent, S3020) for 6 min on a hot rack. Then, samples were incubated with mouse anti-human EGFR antibody (Life Technologies, S3020) at 1:20 for 30 min. Finally, the Dako EnVision ™ + Dual Link System-HRP system (K4065) was used for both anti-HER3 and anti-EGFR-labeled slides according to the manufacturer's instructions, followed by counterstaining with hematoxylin. Stavanger and Emory specimens were scored as in the Nottingham series (H-scores) by the same experienced breast pathologist (for HER3) or independently by two experienced breast pathologists (for EGFR), with the average score used in downstream analysis. Combined HER3-EGFR protein expression was computed as the sum of the individual HER3 and EGFR H-scores. Descriptive statistics can be found in Table S1.
Exploring differences in biomarker scores by cohort. Mean HER3, EGFR, and HER3-EGFR protein expression varied by cohort (ANOVA p < 10 −6 for all); specifically, according to Tamhane's T2 post-hoc test, the means were higher in the Nottingham cohort than the other cohorts (p < 0.05 for all); mean HER3 was higher in the Stavanger than Emory cohort (p < 10 −6 ); mean EGFR was higher in the Emory than Stavanger cohort (p = 1 × 10 −6 ); and mean HER3-EGFR did not differ between Stavanger and Emory cohorts (p = 0.40) (Tables S2  and S3). Variation in HER3, EGFR, and HER3-EGFR was explored via categorical regression with optimal scaling using the CATREG procedure. Specifically, HER3, EGFR, HER3-EGFR, age at diagnosis were scaled numerically and discretized by multiplying; Nottingham grade and AJCC stage were scaled ordinally and discretized by ranking; adjuvant chemotherapy (yes/no) was scaled nominally and discretized by grouping (uniform). The initial configuration was set as multiple systemic starts. All other settings were default. Variation in expression of HER3, EGFR, and HER3-EGFR was primarily explained by the cohort rather than other patient or clinicopathologic variables (Table S4).
RNA sequencing. RNA-sequencing was performed on 60 formalin-fixed paraffin-embedded TNBC pre-treatment samples for which HER3 and EGFR scores were available (full-face slides), which were taken from patients in the Nottingham cohort who were eventually treated with adjuvant CMF chemotherapy. Cases were selected based on which had sufficient material remaining for RNA-seq. RNA-seq data were deposited in ArrayExpress (accession: E-MTAB-6729), where detailed methods can be found. In brief, samples were processed by the Emory Integrated Genomics Core using the Mag-Bind XP FFPE RNA isolation kit (Omega), KingFisher Scientific RepoRtS | (2020) 10:3009 | https://doi.org/10.1038/s41598-020-59514-1 www.nature.com/scientificreports www.nature.com/scientificreports/ Flex magnetic particle separator (ThermoFisher), TruSeq RNA Access library kit (Illumina), Agilent 2200 TapeStation (Agilent Technologies), the QuantiFluor dsDNA System (Promega), and the Agilent High Sensitivity D1000 ScreenTape on the Agilent 2200 Tapestation (Agilent Technologies), and the HiSeq2500 (Illumina, Inc.) per the manufacturers' instructions.  www.nature.com/scientificreports www.nature.com/scientificreports/ Analysis of RNA-seq data. The Salmon index command was used for transcriptome index was building on GRCh38.P10, after which alignment-free transcript abundance was quantified 12 . Gene-level abundance was estimated using tximport 13 . Batch effects were removed using the SVA package 14 . The DESeq2 approach 15 was used to determine differential expression with and without adjusting for age at diagnosis and AJCC stage (Tables S5 and  S6). Ingenuity Pathway Analysis (IPA) was used to identify differentially regulated canonical pathways and causal networks based on 1,378 transcripts (out of 35,590) differentially expressed at the FDR q < 0.05 level in age-and stage-adjusted differential expression analysis.
IPA causal network analysis. Causal network analysis was performed in IPA with the settings adjusted to include only genes, RNAs, and proteins (e.g., rather than drugs or functions). The expression log 2 ratio used to calculate directionality (Z-score). The list of predicted causal networks was filtered to include only hits with significant z-scores (Z-score > 2) without apparent bias. This is to say, we excluded regulators with |µ| < 0.25, where "bias" or  www.nature.com/scientificreports www.nature.com/scientificreports/ testing in subgroups based on adjuvant chemotherapy status (none vs. treatment received) using univariate and multivariate (age and stage-adjusted) Cox regression. Survival analyses were stratified by cohort (i.e., different baseline survival functions were computed for each stratum) to adjust for the potentially confounding influence of cohort differences (due to correlation with protein expression levels and independent prognostic value) without estimating their effects on survival. Results of survival analyses were considered significant if p < 0.05. Analyses were conducted using IBM SPSS Statistics v. 21.

Analysis of PARP1 by IHC. Nuclear non-cleaved PARP1 H-scores in pre-treatment FFPE resection samples
were available from TNBC patients in the Nottingham cohort who ultimately received adjuvant chemotherapy, the staining and scoring of which have previously been described by Green et al. 16 . PARP1 H-scores were compared between low/high HER3-EGFR groups by a two-tailed Mann Whitney U test (n = 38 and n = 93, respectively) using SPSS.

Results
In Kaplan-Meier analysis, high combined HER3-EGFR protein expression in pre-treatment resection specimens from TNBC patients conferred worse 10-year BCSS (83.1% vs 69.2%, log-rank p = 0.017) and DMFS (80.8% vs 70.4%, log-rank p = 0.05) after adjuvant chemotherapy. Moreover, patients with low HER-EGFR has a better prognosis compared to other groups (Fig. 1). In univariate Cox models, high HER3-EGFR was associated with 2.50-fold increased risk of dying from breast cancer after adjuvant chemotherapy (p = 0.003) but not when adjuvant chemotherapy was not administered (Table 1). By contrast, HER3 and EGFR individually had no effect on BCSS regardless of adjuvant chemotherapy status; similar results were obtained in analyses of DMFS. High The significance of each pathway (bar graphs; left y-axis), the log 2 fold-change (line graph; right y-axis), and z-score/activity pattern (bar color). (B) Stacked bar chart depicting the percentage (bar graph; left y-axis) of pathway components up-, down-, or not differentially expressed (bar color), with the number of pathway components given over each respective bar, and the significance of the differential regulation (line graph; right y-axis). Only transcripts with q < 0.05 were entered in the analysis; all pathways depicted exhibit p < 0.05.
HER3-EGFR was associated with 1.95-fold increased risk of distant metastasis after adjuvant chemotherapy (p = 0.022) but not when adjuvant chemotherapy was not administered (Table 1). Neither HER3 nor EGFR significantly impacted DMFS regardless of adjuvant chemotherapy status.
Next, we analyzed differences in pathways, upstream regulators, and networks between HER3-EGFR groups in TNBC patients who received adjuvant chemotherapy to reveal potentially actionable biology in the HER3-EGFR-high, poor-prognosis group. RNA-seq data for a subset of 60 TNBCs in the Nottingham cohort were analyzed (21 with low HER3-EGFR, 39 with high HER3-EGFR). Neither HER3 nor EGFR was differentially regulated based on HER3-EGFR (q = 0.97 and q = 0.95, respectively; Table S5), reflecting a disconnect between transcript levels and expression of these proteins at the plasma membrane. After adjusting for age and stage in differential expression analysis (Table S6), the top differentially regulated canonical pathway was ATM signaling, which was lower in the HER3-EGFR-high group (z = −1.27, p = 0.004; Fig. 2). Aligned with this finding, DNA damage-induced 14-3-3σ signaling was significantly deregulated by HER3-EGFR group, with about half of the pathway genes downregulated and half upregulated (p = 0.017; Fig. 2). The predicted top master regulator whose activity could explain observed transcriptional differences was the enzyme thioredoxin (TXN) ( Table 3). Intriguingly, this causal network includes EGFR, DNA damage sensor PARP1, and caspases 3 and 8 (all predicted to be activated), and p53 and NFκB (predicted to be inhibited) ( Fig. 3 and Table 3). Consistent with these results, the top molecular/cellular function was Cell Death and Survival, with significant activation of genes involved in death/apoptosis of cardiac and brain cells (activation z scores = 0.24 and 0.50, p = 0.014 and 0.016, respectively). Thus, TNBCs exhibiting high HER3-EGFR protein expression may be susceptible to inhibitors of the DNA damage response. Of note, the second-top hit in causal network analysis was reticulon-1 (RTN1), predicted to be activated and a regulator of TP53, TERT, and various MAPKs.  Table 3. Predicted causal networks explaining observed gene expression differences based on combined HER3-EGFR protein expression. www.nature.com/scientificreports www.nature.com/scientificreports/ Finally, we compared uncleaved nuclear PARP1 H-scores between HER3-EGFR groups and found that PARP1 was significantly higher in HER3-EGFR TNBCs from patients who ultimately received adjuvant chemotherapy (p = 0.036; mean rank = 40.24 vs. 53.39 in HER3-EGFR-low and high cases, respectively), consistent with our findings from IPA causal network analysis.

Discussion
Our finding that combined HER3-EGFR protein expression, but not individual HER3 or EGFR protein expression, independently predicts worse BCSS and DMFS following adjuvant chemotherapy for TNBC suggests that HER3 and EGFR should be considered jointly since their expression has been empirically demonstrated to be linked, with decreased expression in one protein inducing compensatory upregulation of the other and thereby promoting resistance to targeted therapies 8,9 . We discovered that a transcriptional network in the HER3-EGFR-high group may be controlled by master regulator TXN (predicted to be inhibited), which may cause activation of PARP1, EGFR, and caspases 3 and 8, while also inhibiting p53, NFκB, AP1, SMAD3, and HIF1A, proteins that are all involved in regulating apoptosis, which may be abnormal in TNBCs with high HER3-EGFR protein expression. Intriguingly, TXN-negative breast cancers are more responsive to docetaxel 17 . Thus, regimens including docetaxel may be effective for HER3-EGFR-high TNBCs. Furthermore, we also identified the potential master transcriptional regulator RTN1, a protein that promotes endoplasmic reticulum-mediated apoptosis and is implicated in neurodegenerative diseases and cancer 18 . Thus, RTN1 may also represent a novel therapeutic target for TNBCs with high HER3-EGFR protein expression.
Limitations of the present study include its retrospective nature and differences in tissue fixation and immunohistochemistry between cohorts. However, strengths include the relatively large sample size, multivariate analysis that adjusted for possible confounders including cohort, and analysis of RNA-sequencing data. This work lays the foundation for in vitro and in vivo studies of combinatorial regimens including dual HER3/EGFR inhibitors, PARP1 inhibitors, and docetaxel-based chemotherapy in TNBCs exhibiting high combined HER3-EGFR protein expression. Importantly, PARP inhibitors have thus far not shown promise in unselected TNBC patients 19 . Therefore, a lack of technically feasible and cost-effective biomarkers to guide selection of TNBC patients for anti-PARP therapy is a critical barrier to progress in the field, which this study may help to address. Our results justify a retrospective analysis of HER3-EGFR in clinical trials or could be the basis for translational sub-projects in upcoming studies for patients with TNBC. Altogether, this work highlights the clinical value of assessing protein expression of HER3 and EGFR in combination which may potentially guide the selection of targeted drugs (dual HER3-EGFR and PARP1 inhibitors) and cytotoxic agents for TNBC patients with poor prognosis after adjuvant chemotherapy.

Data availability
The RNAseq and clinical data are freely available on ArrayExpress (accession: E-MTAB-6729).