Case report: 16-yr life history and genomic evolution of an ER+ HER2- breast cancer.

Metastatic breast cancer is one of the leading causes of cancer-related death in women. Limited studies have been done on the genomic evolution between primary and metastatic breast cancer. We reconstructed the genomic evolution through the 16-yr history of an ER+ HER2- breast cancer patient to investigate molecular mechanisms of disease relapse and treatment resistance after long-term exposure to hormonal therapy. Genomic and transcriptome profiling was performed on primary breast tumor (2002), initial recurrence (2012), and liver metastasis (2015) samples. Cell-free DNA analysis was performed at 11 time points (2015-2017). Mutational analysis revealed a low mutational burden in the primary tumor that doubled at the time of progression, with driver mutations in PI3K-Akt and RAS-RAF signaling pathways. Phylogenetic analysis showed an early branching off between primary tumor and metastasis. Liquid biopsies, although initially negative, started to detect an ESR1 E380Q mutation in 2016 with increasing allele frequency until the end of 2017. Transcriptome analysis revealed 721 (193 up, 528 down) genes to be differentially expressed between primary tumor and first relapse. The most significantly down-regulated genes were TFF1 and PGR, indicating resistance to aromatase inhibitor (AI) therapy. The most up-regulated genes included PTHLH, S100P, and SOX2, promoting tumor growth and metastasis. This phylogenetic reconstruction of the life history of a single patient's cancer as well as monitoring tumor progression through liquid biopsies allowed for uncovering the molecular mechanisms leading to initial relapse, metastatic spread, and treatment resistance.

INTRODUCTION

Breast cancer is the most common malignancy affecting women in the United States, causing the most cancer-related deaths (Howlader et al. n.d.; Breast Cancer Statistics and Resource 2017). However, patients normally do not die from the primary disease, but as a consequence of metastasis (Fouad et al. 2015). Hormone receptor–positive disease accounts for >60% of all breast cancer (American Cancer Society n.d.). The 20-yr recurrence risk of distant metastasis in estrogen receptor (ER)-positive breast cancer is between 10% and 40% and depends on the size of the primary tumor and number of positive lymph nodes (Pan et al. 2017). More than one-half of the recurrence happens after 5 yr of initial diagnosis (Dowling et al. 2019). Little is known about the triggers of metastasis or how the metastatic disease is different from the primary tumor. Developing sequencing technology made it possible to identify the cancer drivers and closely track disease evolution on the genomic and transcriptomic level.

Here we present a patient who was diagnosed with stage III ER-positive breast cancer and developed distant recurrences after 10 years of remission. The patient survived late-stage disease for another 6 years and passed away as a result of extensive disease burden and cancer-related complications. Tissue samples were collected from the patient's breast tumor at initial diagnosis, first recurrence in supraclavicular lymph node, and distant metastasis in the liver. Eleven liquid biopsies were collected after distant metastasis. DNA and RNA from the samples were extracted, sequenced, and analyzed. From the genomic and transcriptomic analysis, we were able to identify potential disease drivers and resistance mechanisms at each stage of the disease, reconstruct the cancer evolutionary tree, and investigate the relationship between the metastases and primary tumor.

RESULTS

Clinical Presentation

A 46-yr-old premenopausal Caucasian female presented with stage III invasive ductal carcinoma of the breast in January 2002. The patient‘s treatment history is summarized in Figure 1A. Solid and liquid samples collected from 2002 to 2017 are gathered in Table 1.

Figure 1.

Patient disease history with treatment and tumor molecular profiles. (A) Disease and treatment history, with red explosions indicating disease progression and an × mark on the time axis indicating treatment on hold because of adverse events. T1, M1, and M2 on the time line indicate when solid tumor samples were collected. Treatment duration is plotted according to the time axis on top. Treatment strip color is assigned per drug class: blue for chemotherapies; red for estrogen receptor antagonists; pink for aromatase inhibitors; purple for mTOR inhibitor; yellow for tyrosine kinase inhibitors; lilac for MEK inhibitor; and green for immunotherapy. (B) Driver mutations from solid tissue samples and liquid biopsies (cfDNA). The table is shaded with a stream graph representing the mutation load during disease progression. Driver mutations are shown with colored bands with the band widths representing mutation allele frequency.

Table 1.

Sample information

Tumor sample						Test	Average coverage
Sample	Specimen site	Collection year	ER	PR	HER2
T1	Right breast	2002	+	+	–	WES	190×
						RNA-seq	46M
M1	Right scalene lymph node	2012	+	–	–	WES	181×
						RNA-seq	75M
M2	Liver	2015	+	+	–	Targeted Sequencing (FoundationOne 315 gene)	462×
Liquid biopsies	Blood	2015–2016	NA	NA	NA	Targeted Gene Sequencing (Guardant 70 gene)	Min 15,000×
Liquid biopsies	Blood	2017	NA	NA	NA	Targeted Gene Sequencing (FoundationAct 67 gene)	10,000×

Solid and liquid samples collected from 2002 to 2017.

At initial diagnosis, the patient's tumor was characterized as ER-positive, progesterone receptor (PR)-positive, and human epidermal growth factor receptor 2 (HER2)-negative by immunohistochemistry (IHC). Right breast mastectomy and axillary lymph node dissection were performed. The patient was treated adjuvantly with doxorubicin and cyclophosphamide followed by paclitaxel (AC-T). She was then maintained on tamoxifen for 2 years, followed by 8 years of aromatase inhibitor therapy (anastrozole).

In May 2012, the patient returned to the clinic with recurrent disease at the right supraclavicular lymph node and possible metastatic nodules in the right upper lobe of the lung and anterior mediastinum. The lymph node metastasis was ER-positive, PR-negative, and HER2-negative. She started on single-agent exemestane. Everolimus was added in July 2012. Radiation therapy to the regional lymphatics and right chest wall was completed in December 2012. In January 2013, restaging positron emission tomography (PET) demonstrated a complete metabolic response. The disease was considered stable for 24 mo under this regimen.

In May 2014, progressive liver lesions were shown on computerized tomography (CT) scan. In response to disease progression, the patient was moved to the next line of treatment with fulvestrant and exemestane.

In September 2014, evidence of progression in the liver, pleura, and bones was found through CT scan. The treatment was switched to single-agent capecitabine. The patient was stable for 6 mo.

In March 2015, new lesions were found in the patient's liver. The previous disease in her mediastinal lymph nodes increased in size. The liver lesion tested ER-positive, PR-positive, and HER2-negative by IHC. In April 2015, treatment was changed to paclitaxel plus everolimus and sorafenib.

To guide treatment and track tumor evolution, genomic testing was performed first through targeted DNA sequencing of a liver biopsy (through FoundationOne test) and subsequently liquid biopsies (through Guardant360 and FoundationAct tests) from August 2015 to November 2017. Mutational findings are summarized in Figure 1B.

In August 2015, all treatments were held for 2 wk because of severe neutropenia. The regimen was then changed to everolimus plus vinorelbine and crizotinib.

In February 2016, because of progression in the liver, the patient started a new line of treatment with trametinib, eribulin, and everolimus. A CT scan in April 2016 suggested a 20% reduction in liver tumor volume. In July 2016, crizotinib was added to the regimen as a result of progression in the liver and the right pleura.

In early September 2016, the patient finished 11 cycles of eribulin. Her CT scan showed slight progression on the chest wall. She was switched to carboplatin plus exemestane and cabozantinib. Everolimus was added to the regimen at the end of December 2016.

In late February 2017, new lesions were detected on the patient's bone scan. The patient started a new regimen of atezolizumab, fulvestrant, everolimus, and cabozantinib. Protein-bound paclitaxel was added in early May. In June 2017, everolimus was held because of thrombocytopenia and neutropenia. Later in June, slight progression was seen on a CT scan, and trametinib was added.

In October 2017, all therapy was held because of her implant that had eroded through the chest wall, which was then surgically removed. The patient continued to have an open area at the right breast but wished to continue with aggressive treatment after surgery. Treatment was restarted at the beginning of November 2017. At the end of November, CT showed progression in the lungs and retroperitoneal lymph nodes. In mid-December 2017, the patient's therapy was changed to protein-bound paclitaxel, cabozantinib, atezolizumab, abemaciclib, and fulvestrant.

At the end of January 2018, the patient struggled with wound healing and a large pneumothorax developed on the right with communication at the wound site. There were no good options to treat the pneumothorax based on her status. Because of the significant tumor burden in the liver, the patient was not capable of tolerating further treatment. The patient was referred to hospice and passed away a week later.

The patient was enrolled into the multicenter clinical trial “Identifying Molecular Drivers of Cancer (CCD)” (NCT02470715). Therapies were selected by incorporating recommendations of a molecular tumor board consisting of oncologists, pharmacists, nurses, genetic counselors, bioinformaticians, patient advocates, and molecular biologists. The therapies ultimately given were furthermore based on the treating physician's consideration of patient preferences, drug toxicities, and availability (i.e., insurance coverage). For administration of drug combinations, we routinely followed Nikanjam et al. (2016), or other data where available. Patients generally did not receive treatment without at least safety data being available. In addition, patients were monitored closely, and adverse event management was planned on the basis of theoretical drug metabolism, with the result that no treatment-related mortality occurred.

Germline Variants

The patient was identified as having the single-nucleotide polymorphism (SNP) FGFR4 p.G388R (rs351855-G/A). No germline alteration was detected among genes commonly included in breast cancer–related or hereditary cancer panels.

Tumor Mutation Burden

The tumor mutation burden at diagnosis in 2002 (T1) was 4 mutations/Mb, which increased to 16 mutations/Mb at first recurrence in 2012 (M1) and 17 mutations/Mb in the liver biopsy in 2015 (M2).

Somatic Mutational Signatures

Somatic mutational signatures are unique combinations of nucleotide mutations generated by certain mutational processes (Alexandrov et al. 2013). Among the 30 COSMIC mutational signatures (Alexandrov et al. 2013; Cosmic 2020), APOBEC signature was detected in all three tumor samples (signature 13 in T1, signature 2 and signature 13 in M1 and M2). Signature 1 (common in cancer) and Signature 10 (POLE) were in all three samples. Signature 7 (mutation possibly due to ultraviolet light exposure) is found in T1 and M1 (Supplemental Fig. 1).

Somatic Mutations

The mutational profile from T1 and M1 samples were from whole-exome sequencing, whereas the M2 liver sample was tested through FoundationOne assay (324 genes) (Fig. 1B; Supplemental Files 1 and 2). We limited our analysis to driver mutations in genes from the FoundationOne assay. Multiple driver mutations affecting the PI3K–Akt signaling pathway were detected (PIK3CA, PTEN, PIK3R2, PIK3CG, MTOR). The phylogenetic tree (Fig. 2) represents the shared history of known and predicted driver mutations among T1, M1, and M2. The two known drivers, PIK3CA Q546R and TP53 E180K, are present in all three solid tumor samples, with TP53 E180K being subclonal (detected in a fraction of the tumor cells in the sample) in T1 but clonal (detected in all tumor cells of the sample) in M1 and M2. The tree shows early branching off between T1 and M1/M2, with T1 containing four private predicted drivers in ALK, APC, CDKN1B, and ERBB2. As the disease progressed, the total number of driver mutations raised from n = 8 in T1 to 10 and 14 in M1 and M2, and the number of private driver mutations were 4, 6, and 12, respectively, signifying evolution in the tumor at each location and progressive disease evolution with time.

Figure 2.

Phylogenetic reconstruction of the tumor history of the patient. Text in bold represents known driver variants as determined by Cancer Genome Interpreter.

The proportion of clonal mutations among the drivers increased at disease recurrence and decreased again when disease further progressed with more driver mutations. At initial diagnosis, only one out of nine driver mutations was clonal and all the others were subclonal. The tumor sample from the first recurrence (M1) shared four of the driver mutations from T1. The clonal status of the driver mutations are presented in Supplemental File 3.

Liquid biopsies (circulating tumor DNA [ctDNA], n = 11) were performed from August 2015 to November 2017 (Fig. 1B; Supplemental File 1; Supplemental Fig. 2). The first two liquid biopsies were negative with initial ctDNA mutations being detected in November 2015. The number of mutations with known cancer-related activity increased to the 20s over time. An ESR1 E380Q mutation, which was not previously detected, appeared in the ctDNA by the end of 2016 after exemestane treatment was restarted and had rapidly increased in allele frequency by the end of 2017 (from 0.9% to 28.3%) (Supplemental Fig. 2).

Copy-Number Variations

No copy-number variations were detected in the analyzed samples.

Differential Expression Analysis

RNA-seq analysis of T1 and M1 revealed gene expression changes between initial diagnosis and disease recurrence after 10 years. Seven hundred and twenty-one genes were significantly differentially expressed between T1 and M1 (P < 0.05), among which 298 genes are differentially expressed with log2 fold change ≥3 with 184 genes up-regulated and 114 genes down-regulated in M1 (Supplemental File 4). Table 2 lists the top 15 differentially expressed genes that are reported to be breast cancer–related or intersected with the Cosmic Cancer Gene Census (https://cancer.sanger.ac.uk/census). At disease recurrence, several tumor-suppressor genes were down-regulated, including RGS7, PTPRT, and TFF1. Down-regulation of estrogen-dependent genes (CCND1, TFF1, PGR, GREB1) might be the result of prolonged use of aromatase inhibitors (Mackay et al. 2007). Among the up-regulated genes in M1, SOX2, MUC16, PTHLH, S100P, and CEACAM6 have been reported to promote tumor growth and metastasis in various cancer types (Arumugam et al. 2005; Dowen et al. 2005; Fuentes et al. 2007; Lewis-Wambi et al. 2008; Leis et al. 2012; Urosevic et al. 2014; Liu et al. 2017; Wuebben and Rizzino 2017; Kanwal et al. 2018).

Table 2.

Top differentially regulated genes between T1 and M1

Up-regulated genes				Down-regulated genes
Gene Symbol	EntrezID	Log2FC	FDR	Gene Symbol	EntrezID	Log2FC	FDR
SOX2	6657	7.60	4.34 × 10−25	TFF1	7031	−6.37	2.08 × 10−27
PTHLH	5744	7.33	6.96 × 10−24	PGR	5241	−6.05	6.48 × 10−25
MUC16	94,025	4.55	3.03 × 10−07	GREB1	9687	−5.72	1.53 × 10−22
S100P	6286	4.08	2.16 × 10−05	SCGB2A2	4250	−5.04	2.67 × 10−17
SSX1	6756	3.67	4.49 × 10−04	CYP24A1	1591	−4.89	1.52 × 10−16
CEACAM6	4680	3.63	5.45 × 10−04	DACH1	1602	−4.51	3.56 × 10−14
DCAF12L2	340,578	3.37	3.80 × 10−03	AREG	374	−4.50	3.48 × 10−14
TNFSF10	8743	3.30	5.39 × 10−03	RGS7	6000	−4.00	2.45 × 10−11
MLLT3	4300	3.21	9.23 × 10−03	PTPRT	11,122	−3.94	4.19 × 10−11
CEACAM1	634	3.15	1.42 × 10−02	NTRK3	4916	−3.83	1.57 × 10−10
MUC4	4585	3.05	2.45 × 10−02	SCUBE2	57,758	−3.73	5.18 × 10−10
HES6	55,502	3.04	3.04 × 10−02	KRT5	3852	−3.51	6.90 × 10−09
HUNK	30,811	3.01	3.29 × 10−02	IGF1R	3480	−3.19	1.42 × 10−07
TMPRSS2	7113	2.94	4.81 × 10−02	COX6C	1345	−3.07	4.81 × 10−07
KLF5	688	2.93	4.85 × 10−02	CCND1	595	−2.93	1.73 × 10−06

All genes listed in this table are in Cancer Gene Census (CGS) by Catalog of Somatic Mutations in Cancer (COSMIC) or have been reported to be related to breast cancer.

(FC) Fold change, (FDR) false discovery rate.

Pathway analysis between T1 and M1 showed enrichment (P < 0.05) for the hallmark gene sets “Estrogen Response early and late” and “Epithelial mesenchymal transition.” A KEGG pathway analysis suggested that the “Olfactory transduction pathway” was up-regulated in M1. Down-regulated pathways include “Protein digestion and absorption,” “Neuroactive ligand–receptor interaction,” “ECM–receptor interaction,” “PI3K–Akt signaling pathway,” and “Estrogen signaling pathway.” Intersections of the differentially regulated pathways are presented in Supplemental File 5.

DISCUSSION

Breast cancer is the leading cause of cancer death in women around the world, the vast majority of those deaths are caused by metastatic disease (Breast Cancer - Metastatic - Statistics 2017). Through phylogenetic reconstruction of the presented patient's cancer, an increased mutational rate over time, from the primary tumor at diagnosis to samples taken at relapse, is revealed. In addition, a high level of tumor heterogeneity was found with metastatic sites sharing only two driver mutations with the primary tumor but having acquired additional de novo driver mutations. We also demonstrate how transcriptome profiling aids in uncovering the mechanisms that lead to the patient's initial relapse, indicating expression profile changes that promote tumor progression and metastasis as well as reduced sensitivity and resistance to therapy.

Patient Case Tumor Evolution

It is known that tumors are a heterogenous group of cells in regards to their genomic profile (Dagogo-Jack and Shaw 2018). Mutations present only in a subgroup of the tumor cells are defined as subclonal mutations, whereas mutations present in all of the tumor cells in the sample are defined as clonal mutations (Lawrence et al. 2013). A clonality analysis of the tumor sample at diagnosis supported this fact with PIK3CA Q546R identified as the only clonal mutation with an additional eight subclonal driver mutations.

As disease progressed, the patient's cancer evolved over time with changes to the mutational profile. According to the branching evolution theory, clones diverge from common ancestor and expand simultaneously while accumulating new driver mutations. Clones that carry treatment resistant mutations gain advantages and expand (Davis et al. 2017). As seen in this case, the TP53 E180K mutation, which was initially subclonal in 2002 (T1), presented clonal in the relapse sample in 2012 (M1) as well as later samples. TP53 encodes the tumor protein p53 (p53). Hormonal therapy in ER-positive breast cancer depends on p53-mediated cell cycle arrest (el-Deiry et al. 1993; Varma and Conrad 2000). Oncogenic mutations of TP53, such as the one observed here, can decrease the effect of hormonal therapy (Love 1989). In this case, tumor cells that harbor the TP53 mutation may have gained survival advantages and expanded while the patient was maintained on hormonal therapy in the 10 years prior to relapse. We also observed an accumulation of de novo mutations. Although the number of driver mutations increased with disease progression, the majority of driver mutations in metastatic samples were de novo. Half of the driver mutations (n = 4) in the 2002 biopsy were not detected in later samples. By 2015, PIK3CA Q546R and TP53 E180K were the only two driver mutations shared with previous samples (Fig. 1B; 2). A similar trend of emerging de novo mutations as well as competing and evolving subclonal mutations was reported previously within a larger patient cohort (Miller et al. 2016). The evolving mutational landscape can be the result of natural disease progression or selective pressure from treatment. No matter the cause, considering tumor heterogeneity and the dynamic changes in tumor genetic alterations, genomic profiling from single tumor tissue biopsies may not capture all alterations within the tumor. Repeat liquid biopsies may overcome this limitation (Esposito et al. 2016).

A sample-based phylogenetic reconstruction of the patient's tumor evolution (Fig. 2) revealed an early branching off between the initial tumor and metastatic samples, suggesting a change of mutation profile between primary tumor and metastasis. Further changes are observed between M1 and M2 as well. The increase in tumor mutational burden, numbers of total driver and private driver mutations with the progression of disease, show that despite continuous therapy, the disease continued to evolve and become heterogenous, before and after recurrence. Although multiregional sequencing of each tumor biopsy could provide more insight and resolution into such heterogeneity, the tests performed here represent what is routinely available to a clinical provider. Tracking the presence and clonality of actionable variants such as PIK3CA Q546R and TP53 E180K can be very beneficial in the context of precision oncology.

Mutational Signatures

Apart from the driver mutations discussed above, most mutations identified were classified as passenger mutations, which are generally considered to not be involved in tumorigenesis (Pon and Marra 2015). However, the combinations of mutation types arising from specific mutagenesis processes can be characterized as somatic mutational signatures. Inference of somatic mutational signatures based on COSMIC mutational signatures on the three solid tissue samples (T1, M1, M2), revealing presence of APOBEC-related mutational signatures. Apolipoprotein B mRNA editing enzyme (APOBEC) is a family of cytidine deaminases. The APOBEC protein family protects humans from viral infections through causing lethal viral cDNA hypermutations (Harris and Liddament 2004) but can be a driving force of somatic hypermutation in cancer (Roberts and Gordenin 2014), resulting in tumor evolution, higher heterogeneity and rising treatment resistance, and poor prognosis (Swanton et al. 2015). The identification of somatic mutational signatures holds its value in supporting new therapeutic options and revealing possible resistance. In this patient case, the detection of APOBEC activity via mutational signature analysis could have presented the option of using PARP inhibitors or ATR inhibitors (Kanu et al. 2016; Buisson et al. 2017; Green et al. 2017; Nikkilä et al. 2017; Ma et al. 2018) and utilizing aggressive treatment plans (Swanton et al. 2015). In addition, tamoxifen administration is suggested to be used with close monitoring because of potential APOBEC3B-mediated tamoxifen resistance (Law et al. 2016).

Germline Analysis

Germline analysis identified the SNP FGFR4 p.G388R (rs351855-G/A), which is a common SNP with global minor allele frequency of 0.29952 (A), changing glycine (G) at FGFR4 codon 388 to arginine (R) (Ulaganathan et al. 2015). FGFR4, as a member of the fibroblast growth factor receptor family, is a tyrosine kinase receptor for fibroblast growth factors (FGFs). When activated, FGFR4 fires up the downstream signaling, including the PI3K–Akt and RAS–RAF pathways. It can also phosphorylate signal transducer and activator of transcription 3 (STAT3) directly (Tang et al. 2018). With FGFR4 Arg388 allele, as seen in this patient, the receptor exposes an additional membrane-proximal cytoplasmic STAT3 binding site. The additional STAT3 binding site recruits STAT3 to the cell membrane and increases STAT3 phosphorylation (Ulaganathan et al. 2015). When overactivated, STAT3 promotes tumorigenesis through up-regulating tumorigenesis genes, releasing inflammatory cytokines and contributing to epithelial–mesenchymal transformation (EMT) (Real et al. 2002; Leslie et al. 2006; Chen and Han 2008; Ogura et al. 2008; Thiery et al. 2009).

The overactivation of STAT3 due to the germline variation at FGFR4 may have acted as a resistance mechanism to everolimus treatment in this patient (Fig. 3). In this patient, by blocking mTORC1 with everolimus, one would expect down-regulated transcriptional activity of STAT3, which is one of the downstream effectors of mTORC1. However, with the FGFR4 G388R variant, STAT3 was phosphorylated independent of mTORC1. Blocking mTORC1 with everolimus in this patient might have not been as efficient as for patients with wild-type FGFR4. Because of sample limitation, we were unable to examine the actual phosphorylation level of STAT3. The effect of FGFR4 G388R in cancer recurrence in this patient remains a hypothesis.

Figure 3.

Pathways affected by mutations. DNA-based analysis from solid tumors revealed multiple somatic mutations in the PI3K–Akt signaling pathway (PIK3CA, PIK3R1, PTEN, MTOR) across all samples. Other mutations were detected in both solid and liquid samples, affecting RAS–RAF, cell cycle, and other critical pathways regulating apoptosis and cell proliferation. Germline variant FGFR4 p.G388R was detected. It creates an additional binding site for STAT3 at the cell membrane so it can be phosphorylated without upstream signals from preclinical studies. If the above conclusion holds true in this patient, the activated STAT3 might have established a resistance mechanism for treatment targeting the PI3K-mTOR pathway.

Transcriptome Profiling

In addition to genomic alterations identified through DNA profiling, transcriptome profiling can reveal changes in gene expression levels, which may affect protein expression and signal transduction. For this case report, we were especially interested in how gene expression is regulated, as the patient did not present with DNA copy-number changes, as well as which genes were differentially expressed between the metastatic and primary tumors. Comparing the transcriptome profile from the first recurrence in 2012 (M1) with the sample from the initial diagnosis in 2002 (T1), we identified multiple oncogenes to be up-regulated as well as tumor suppressors to be down-regulated.

Down-Regulation of PGR in M1 Was in Agreement with Results from IHC Staining

TFF1 and CCND1, both estrogen-dependent genes, were down-regulated in M1, possibly because of aromatase inhibitors used for 10 years (Mackay et al. 2007). However, the anticancer effect of down-regulated CCND1 could have been diminished as a result of a RB1 Q217* nonsense mutation appearing in 2012. The RB1 gene encodes retinoblastoma protein (Rb). Rb functions as a tumor suppressor through binding to E2F, a transcription factor, to prevent cell division caused by E2F-dependent transcription (Dick and Rubin 2013). CCND1 phosphorylates Rb, which dissociates the Rb–E2F complex. The released E2F would then start the subsequent transcription and promote cell division (Ortega et al. 2002). Down-regulation of CCND1 with hormonal therapy is supposed to stop this process by reducing Rb phosphorylation and limiting the free E2F level, but RB1 Q217* causes an early termination of the peptide leading to a missing E2F binding domain (aa380–787) (Dick and Rubin 2013). Rb thus loses the function of forming the Rb–E2F complex and the cellular free E2F level is not affected by the phosphorylation status of Rb or CCND1 expression level.

In the presented case, we hypothesize that although the CCND1 transcription level was successfully decreased by hormonal therapy, the RB1 mutation, which was not present at diagnosis, was one of the causes for treatment resistance and disease recurrence in 2012.

Among the differentially regulated genes, down-regulation of GREB1 and up-regulation of SOX2 may have presented an additional resistance mechanism to the patient's hormonal therapy and treatments targeting the PI3K–AKT pathway. GREB1 is a tumor promoter. In ER+ breast cancer, ESR1, when activated by E2, binds to estrogen response elements upstream of the GREB1 promotor and up-regulates GREB1 expression (Deschênes et al. 2007; Sun et al. 2007; Yang et al. 2014). Higher E2 level is correlated with GREB1 overexpression (Sun et al. 2007). When up-regulated, GREB1 promotes cell proliferation and migration. GREB1 knockdown does not affect migration but can decrease cell proliferation rate and slow down disease progression (Hodgkinson et al. 2018). The GREB1 down-regulation observed here in 2012 was probably achieved by the prolonged hormonal therapy the patient was on. In 2002 GREB1 expression was elevated as a result of ER-positive status. From 2002 to 2012, the patient was on tamoxifen and anastrozole. Decreased E2 production and ESR1 activity lead to down-regulation of GREB1 expression. Down-regulation of GREB1 was supposed to control disease progression; however, a recent study suggested that the growth suppression caused by GREB1 knockdown can be reversed by constitutive activation of the PI3K–Akt pathway (Haines et al. n.d.), which was the case in this patient.

SOX2 was the topmost up-regulated gene in 2012 with a log2 fold change of 7.6 compared to the 2002 sample. SOX2 overexpression promotes cell invasion through activating PI3K–Akt pathway (Yang et al. 2014). In this patient, the growth suppression caused by GREB1 down-regulation through hormonal therapy may have been reversed by SOX2 up-regulation, presenting a resistance mechanism for hormonal therapy. Because of sample limitation, we were not able to assess SOX2 expression status in later stages of the tumor evolution, but if SOX2 up-regulation remained after 2012, it might have also caused resistance toward everolimus, which was targeting the PI3K–Akt pathway through blocking mTOR (Strimpakos et al. 2009).

In addition, we performed pathway enrichment analysis, which revealed EMT-related genes enriched in the 2012 sample. EMT is one of the important processes in cancer development when cells acquire mesenchymal properties and lose the original shape and reduce intercellular adhesion. The loosely bonded, small, and flexible cells are more capable of invasion and spreading the disease to distance (Thiery et al. 2009). One of the critical events in EMT is the up-regulation of PI3K–Akt pathway (Larue and Bellacosa 2005), which was present in this patient as discussed above.

Copy-number variations (CNVs) are commonly found in advanced stage breast cancer (Zhang et al. 2018). However, no CNV was detected in this patient. It is an unusual finding, but because analysis in the advanced stage samples was limited to either targeted sequencing or liquid biopsies, CNV changes outside the targeted regions might have gone undetected.

Clinical Implications

Unraveling the molecular makeup of a patient's tumor has become an important factor to deliver treatment plans that best suit the patient (Schmidt et al. 2016; Sicklick et al. 2019). One of the best clinical responses received for this patient was achieved in 2016 with 20% tumor reduction as determined by imaging. At the time, the patient was treated with a combination of trametinib, everolimus, and eribulin. According to the molecular profile of the tumor, during this period, the PI3K–mTOR and RAS–RAF pathways were two of the major pathways affected by mutations in PIK3CA, PTEN, mTOR, and NF1. We hypothesize that at that time this regimen precisely targeted the tumor cells, with the MEK inhibitor trametinib blocking the RAS–RAF pathway, everolimus acting on mTORC1 to down-regulate the PI3K–mTOR pathway, and eribulin disrupting the microtubule network and decreasing cell division, which lead to a reduction in tumor size.

As stated, the tumor genome changed with disease progression. To achieve the best clinical response, one hypothesis would be to constantly change treatment strategy in response to the changing molecular profile of tumor cells (Xu et al. 2017; Kosovec et al. 2018; Xu et al. 2019; Zhang et al. 2019). However, taking repeated tissue biopsies is often not feasible. Liquid biopsy, which tests for ctDNA in the plasma, provides a noninvasive, less expensive alternative to solid tissue biopsy since being introduced in the mid-2010s (Davis et al. 2017; Sheridan 2019). We monitored this patient with repeat liquid biopsies after distant metastasis was detected in mid-2015. Initial ctDNA profiling was negative. When the patient presented with active disease in the liver, mutations started to be picked up in the plasma samples 2 months later, with an overall trend of increasing mutation numbers, as well as increasing allele frequencies, which we considered as an indication of increased tumor burden. A positive association between plasma ctDNA and tumor burden has been reported in other studies as well (Bettegowda et al. 2014; Lee et al. 2018).

Treatment resistance is the leading cause of therapeutic failure (Diaz et al. 2012; Aparicio and Caldas 2013). New driver mutations are an important source of treatment resistance. Detecting rising resistance mechanisms and adjusting treatment plans accordingly in a timely manner is crucial to control the disease (Housman et al. 2014). Through liquid biopsy, we observed an emerging ESR1 E380Q mutation that could be a possible acquired resistance mechanism to aromatase inhibitor therapy in this patient (Li et al. 2013; Jeselsohn et al. 2015). The patient was on aromatase inhibitors (anastrozole and exemestane) for 10 years before distant metastasis appeared in 2014. Her treatment was switched to chemotherapy and kinase inhibitors for 2 years (09.2014–09.2016) before exemestane was added back to the regimen in September 2016 when the patient progressed in the liver and right pleura. The ESR1 E380Q mutation, which had not been detected in previous samples, in neither solid (2002–2015) nor liquid biopsies (08.2015–07.2016), appeared in the first ctDNA test (11.2016) after exemestane was added in September 2016 (Supplemental Fig. 2). ESR1 encodes for estrogen receptor alpha (ERɑ). ERɑ, when combined with estrogen, travels to the nucleus and promotes transcription of genes involved in cell proliferation and survival through binding to the promoter sequence of the genes (Ip et al. 1979; Katzenellenbogen et al. 1987). This process is considered one of the key components of tumorigenesis in ER-positive breast cancer (Rose et al. 1980; Sunderland and Osborne 1991). Aromatase inhibitors are commonly used to stop this process by reducing estrogen production (Buzdar 2003). However, alterations in the ligand binding domain of ESR1, such as observed here, would lead to estrogen-independent constitutive transactivation (Toy et al. 2013), in which case estrogen is not needed in ERɑ activation, and lower estrogen levels through aromatase inhibitor will not decrease the activity of ERɑ to promote tumor growth. ESR1 ligand binding domain mutations have been observed as acquired resistance mutations to aromatase inhibitors (Jeselsohn et al. 2015). A study of 171 patients with advanced breast cancer found ESR1 mutations exclusively in ER-positive breast cancer patients previously exposed to aromatase inhibitors (Schiavon et al. 2015). Another study observed that even with confirmed ESR1 E380Q mutation in tumor tissue, the variant went undetected in plasma DNA (Takeshita et al. 2017). The authors believed that this was because the ESR1 E380Q mutation in the tumor was only in a subgroup of cancer cells. In our case, we observed the emergence of ESR1 E380Q mutation in the first liquid biopsy after exemestane was restarted. The allele frequency of the mutation started low at 0.9 but quickly increased in the following tests (Supplemental Fig. 2). We cannot rule out the possibility that the ESR1 mutation was present in the primary tumor cells, went undetected because of tumor heterogeneity or low concentration of ctDNA, and expanded after exemestane was restarted. Based on previous research (Miller et al. 2016), ESR1 E380Q is more likely an acquired resistant mutation after prolonged exposure to aromatase inhibitors. Eventually, exemestane was switched to fulvestrant, an ER antagonist, to control ER overactivation and overcome resistance to aromatase inhibitors.

Another potential treatment resistance mechanism we observed in this patient was the occurrence of PTEN mutations, which were firstly detected in liver metastatic tissue in 2015 and later observed in ctDNA from 2016 on. It has been suggested that PTEN loss of function could lead to resistance to the PI3K inhibitor alpelisib (BYL719) (Juric et al. 2015). Given the PIK3CA driver mutations, this patient would have been a candidate for PI3K inhibition if the drug was available during her treatment. However, she may not have gotten the best clinical outcome from PI3K inhibition because of her PTEN mutations.

Genomic and transcriptome profiling makes it possible to generate patient-specific cancer molecular profiles. In addition to “traditional” pathology (IHC, FISH)-derived biomarkers, those genomics-derived molecular profiles provide additional insights on oncogenesis and disease progression. The relevance between molecular aberrations, disease progression, and potential treatment can be direct, when proven by clinical trials, or indirect, if the evidence is still being investigated in early-stage trials or preclinical studies only. As an example, the germline variant FGFR4 G388R as identified in this patient would be considered as indirect evidence (preclinical study). Molecular data that has not been demonstrated to have significant clinical utility can still be informative to a clinician, regarding which pathways are activated and potentially driving tumor growth and progression. These biomarkers would not be utilized directly to make treatment decisions, but could be used to inform the clinician on potential resistance mechanisms and possibility of a poor response.

Conclusion

In summary, we demonstrate that molecular profiling can provide valuable information to guide decision making for a patient's cancer treatment. Although genomic profiling has been utilized in clinical practice more in recent years, transcriptome analysis has rarely been considered. As of today, no transcriptome profiling test has been approved for clinical application. In this case report, with detailed analysis of information acquired from germline DNA alterations, somatic mutations, and gene expression profiling, we hope to emphasize the importance of closely monitoring cancer genomic profile and the significance of including transcriptome profiling while treating cancer patients.

METHODS

Sample Collection

Tumor tissues were collected from the primary breast tumor (T1) in 2002, initial recurrence at right supraclavicular lymph node (M1) in 2012, and further metastasis in the liver (M2) in 2015. Peripheral blood was collected for cell-free DNA analysis and matched germline sequencing. Sample names, collection time, and pathology characteristics, together with genomic tests performed and average coverage, are presented in Table 1.

Tissue Preparation and Sequencing

DNA/RNA extraction and sequencing from T1, M1, and normal blood (buffy coat) was performed at Avera Institute of Human Genetics. In short, we utilized one rapid run flow cell on the Illumina HiSeq2500 platform (Illumina) in 2 × 125 fashion. Lane 1 housed whole-exome libraries (Clinical Research Exome v1; Agilent), whereas lane 2 retained the RNA libraries (TruSeq RNA Access; Illumina). Appropriate adapter trimming was performed and a spike-in of 1% PhiX control for each well was used.

Sample M2 was fixed in formaldehyde and embedded in paraffin wax. The formalin-fixed, paraffin-embedded (FFPE) sample was then sent to Foundation Medicine for sequencing using FoundationOne assay, which covers 3769 exons of 236 cancer-related genes and 47 introns of 19 genes with frequent rearrangements, at an average depth of 250× (Frampton et al. 2013).

ctDNA was isolated and sequenced by the Guardant360 assay from August 2015 to November 2016. The remaining samples from February 2017 to November 2017 were sequenced through the FoundationAct assay. Guardant360 is a 73-gene ctDNA sequencing (average coverage = 8000×) test with analytic validation (Lanman et al. 2015; Odegaard et al. 2018). FoundationAct tests for ctDNA with a coverage of 62 genes to the depth of 5000× (Clark et al. 2018).

Variant/Copy-Number Calling and Annotation

Raw exome sequencing data from the two tumor (T1, M1) and matched normal samples were processed through multiple tools: BBDuk from BBTools suite 36.32 (BBDuk Guide n.d.; SourceForge n.d.) to trim adapters, FastQC 0.11.5 (Babraham Bioinformatics n.d.) and Qualimap 2.2 (García-Alcalde et al. 2012) for quality control, BWA-MEM 0.7.15 (Li 2013) to generate mapped alignments to the human genome reference (hs37d5), sambamba 0.6.3 (Tarasov et al. 2015) to mark duplicate reads, ABRA 0.97 (Mose et al. 2014) to realign indel regions. VarDict 1.5.1 (Lai et al. 2016) was used to call somatic and germline variants, whereas CNVkit 0.9.5 (Talevich et al. 2016) was used to call copy-number variants, and ANNOVAR (Wang et al. 2010) was used to annotate all variants.

A previously published method of checking for cross-sample contamination was used to verify tumor-normal matching (Stephens et al. 2018). The number of nonsynonymous variants per 1 Mb of DNA adjusted to the exome panel size was calculated as the tumor mutation burden of a sample. Tumor purity was computationally estimated by Sequenza (Favero et al. 2015). Palimpsest (Shinde et al. 2018) was used to calculate cancer cell fraction (CCF) and determine clonality status for the somatic variants per sample. A phylogenetic tree for the patient was constructed by LICHeE (Popic et al. 2015) based on CCF data for variants classified as drivers by CancerGenomeInterpreter (Tamborero et al. 2018). The R package deconstructSigs (Rosenthal et al. 2016) was used to infer Catalog of Somatic Mutations in Cancer (COSMIC) mutational signatures v2 (Cosmic 2020).

Differential Gene Expression

iDEG (individualized differentially expressed genes) (Li et al. 2019), a method to compare two conditions without replicates, was used to identify differentially expressed genes (DEGs) between T1 and M1.

Pathway Analysis

For downstream pathway analysis, we selected DEGs with a false discovery rate of ≤0.05 and absolute fold change greater than three. Gene sets from KEGG, Reactome, Biocarta, MSigDB Hallmark, Progeny, and Wikipathways were downloaded in gmt format. DEGs and their normalized gene expression values at the time points T1 and M1 were tested using MixEnrich (Li et al. 2017) to identify bidirectional and concordantly dysregulated pathways. DEGs were tested for their enrichment in KEGG pathways using gProfileR (Reimand et al. 2007) package in R.

Drug–Gene Interaction

Up-regulated DEGs were used to check for their interaction with drugs in DGIdb (Drug Gene Interaction database) (Griffith et al. 2013; Wagner et al. 2016; Cotto et al. 2018).

ADDITIONAL INFORMATION

Data Deposition and Access

Raw sequencing data has been uploaded to European Genome-phenome Archive (https://www.ebi.ac.uk/ega/) and is available under accession number EGAS00001004624.

Ethics Statement

The patient provided written and informed consent. This study was approved by the WIRB (protocol #20140659).

Acknowledgments

The authors are grateful to the patient and her family for donating medical data and samples to this research. We thank Karla Flak for her assistance. The work of the authors was supported by the Avera McKennan Hospital & University Health Center and Avera McKennan Foundation.

Author Contributions

Study concept and design was done by B.X., C.B.W., and T.M. Sample acquisition, preparation, and sequencing were performed by C.D., S.T., S.V., and J.P. The data were analyzed by B.X., A.A., P.S., and T.M. Administrative, technical, or material support was provided by G.D., E.E., and C.B.W. The study was supervised by E.E., C.B.W., and T.M. The drafting and revision of the manuscript was done by all authors. All authors have read and agreed to the published version of the manuscript.

Competing Interest Statement

The authors have declared no competing interest.