Intron 1–Mediated Regulation of EGFR Expression in EGFR-Dependent Malignancies Is Mediated by AP-1 and BET Proteins

This article is featured in Highlights of This Issue, p. 2143

Mammalian cells contain thousands of transcriptional control elements known as enhancers responsible for the regulation of gene expression (1). Enhancers contain open chromatin (2) and are marked by specific chromatin modifications, including monomethylation of histone H3 at lysine 4 (H3K4me1) and acetylation of histone 3 at lysine 27 (H3K27Ac; ref. 3). Active enhancers are typically identified by a specific enrichment of H3K27Ac (4) and contain high levels of enhancer-associated transcription factors (TF; ref. 5) that are often cancer-specific (6). The activator protein-1 (AP-1) family of oncogenic TFs activates transcription of different genes through homodimers and heterodimers of Jun, Fos, and ATF family members (7) and has been shown to be critical for maintenance of GBM transcriptional heterogeneity (8). These factors bind to and modulate enhancer activity in combination with other chromatin-associated proteins including BRD4 (9, 10), the YAP/TAZ/TEAD complex (11, 12), and the SWI/SNF (BAF) complex (13).

The epidermal growth factor receptor (EGFR) is a transmembrane tyrosine kinase whose downstream signaling pathways modulate a wide range of cellular activities, including growth, migration, and survival (14). EGFR is frequently overexpressed in a variety of cancer types, including cancers of the head and neck (HNSCC) and glioblastoma (GBM; ref. 15). Overexpression of EGFR is detectable in as much as 84% of HNSCC tumors, with mutation and/or amplification occurring in approximately 31% of these tumors (16), suggesting high EGFR protein levels are driven by transcriptional control mechanisms in HNSCC. GBM possesses a stronger correlation between EGFR copy number and expression, with mutation and/or amplification occurring in approximately 46% of GBM, over 90% of which overexpress the protein (17). Recent large-scale analysis of cancer epigenomes identified a significant relationship between somatic copy-number alterations (SCNA) and enhancer expression, with the most significant increases in enhancer expression occurring in tumors that have high aneuploidy and high mutation load (18). HNSCC and GBM have high SCNA frequency (19) and a high frequency of EGFR gene alterations (20), indicating epigenome hyperactivity may play a role in overexpression of EGFR.

In spite of the prevalence of EGFR dependency in solid tumors, few studies have attempted to elucidate the mechanisms of transcriptional control of the gene. Early studies identified cis-acting elements that may regulate EGFR, including CA dinucleotide repeats (21), intron 1 DNase I hypersensitive sites (22), and cooperative promoter-upstream and intron 1 enhancers (23). Lack of utilization of next-generation sequencing (NGS) techniques limits the scope of these studies and argued for a larger breadth of analysis. Recently, EGFR superenhancers were identified in various cancer types, including cervical, glioma, and lung (8, 10, 24–26), however systematic mapping of these super enhancers was lacking. The increased interest in enhancer-mediated expression of EGFR strongly argues for a clearer understanding of the regions critical for enhancer function at EGFR and the TFs that mediate these effects.

In this study, we interrogated the transcriptional regulation control of the EGFR locus. We utilized cell line models of solid tumors which commonly overexpress EGFR to identify the cis- and trans-acting factors involved in the transcriptional control of the EGFR gene and identified two critical constituent enhancers located within intron 1 of EGFR. We characterized these enhancers by mapping their domains with luciferase reporter assays, chromatin interaction assays, and CRISPR/Cas9-mediated genetic perturbation. Finally, genome-wide motif analysis implicated AP-1 TF and BET (bromodomain and extraterminal domain) protein family members whose presence and activity in EGFR intron 1 were further functionally validated. The results revealed that two constituent enhancers regulate EGFR transcriptional control, proposing an epigenetic regulation of the EGFR gene in cancer.

GBM cell lines U87 (RRID:CVCL_0022), T98G (RRID:CVCL_0556), and LN229 (RRID:CVCL_0393) were purchased from ATCC. LNZ308 (RRID:CVCL_0394) was provided by Dr. Erwin Van Meir (Emory University, Atlanta, GA). SF767 (RRID:CVCL_6950) was provided by Dr. Mitch Berger (UCSF Brain Tumor Center, San Francisco, CA). Head and neck cell lines HN12 (RRID:CVCL_5518), Cal27 (RRID:CVCL_1107), and Detroit562 (RRID:CVCL_1171; provided by Dr. Silvio Gutkind, UCSD Moores Cancer Center, San Diego, CA) were maintained in DMEM (HyClone, #SH30022.01) supplemented with 10% fetal bovine serum (Atlanta Biologicals, #S12450) and 1% penicillin–streptomycin (Gibco, #15140-122) and grown as adherent cultures. GBM neurosphere cell lines GSC23 (RRID:CVCL_DR59, Dr. Fred Lang, University of Texas MD Anderson Cancer Center, Houston, TX) and TS576 (Dr. Cameron Brennan, Memorial Sloan Kettering Cancer Center, New York, NY) were maintained in DMEM/F12 (Gibco, #11320-033) supplemented with B27 supplement (Gibco, #12504-044) and 1% penicillin–streptomycin and grown in suspension. Mycoplasma testing was performed with the PlasmoTest kit (InvivoGen, #rep-pt1) and found to be negative. All experiments are performed within 10 passages of the original frozen stock or post-manipulation.

DNA fragments tested in the luciferase reporter assay were cloned from human genomic DNA (Promega, #G3041). PCR products were cloned downstream of Firefly luciferase in the pGL4.24 minimal promoter vector (Promega, #E8421) using the SalI (NEB, #R3138S) site. Constructs were sequence confirmed by Sanger sequencing using the pGL4.24-R primer. pMIEG3-JunDN (RRID: Addgene_40350) and 3xAP1pGL3 (RRID: Addgene_40342) were gifts from Alexander Dent. pMIEG3-Empty was created by removing the JunDN sequence by EcoRI digestion. For each transfection reaction, 100 ng control plasmid expressing Renilla luciferase (Promega, #E2241) and 1 μg Firefly luciferase construct were cotransfected with Lipofectamine 2000 (Thermo Fisher, #11668030) into 2 × 10⁵ cells in a 12-well plate. After 24 hours, cells were collected in 1× PLB. Luciferase activity was measured by the Dual-Luciferase Reporter Assay System (Promega, #E1910) on a Tecan Spark 10M with injection control. Transfection efficiency was controlled for by dividing Firefly luminescence by Renilla luminescence, and final activity was normalized to a negative control.

RNA was extracted with the RNeasy Plus kit (Qiagen, #74134) according to the manufacturer's instructions. Reverse transcription of mRNA was performed using 1 μg RNA with the iScript Reverse Transcription Supermix (Bio-Rad, #1708841). For real-time PCR analysis, 5 μL of cDNA (50 ng of starting RNA) was amplified per reaction using the iTaq Universal SYBR Green Supermix (Bio-Rad, #1725124) and the Bio-Rad CFX96 qPCR system.

The experiment was performed as described (27) with the following modifications. Nuclei were treated with 1000 U EcoRI (NEB, #R3101S) at 37°C overnight. 100 U T4 enzyme (NEB, #M0202S) was added to digested nuclei and incubated at 16°C for 4 hours. Another 100 U T4 enzyme was added, and nuclei were incubated with rotation at 4°C overnight. Ligated DNA (150 ng) was quantified in triplicate by TaqMan real-time PCR using the PrimeTime Gene-Expression Supermix (IDT, #1055772). Control 3C template was generated by using two bacterial artificial chromosomes (BAC) encompassing the entire EGFR gene, RP11-159M24 and RP11-148P17, were purchased from the Children's Hospital Oakland Research Institute (CHORI). Equimolar of the two BACs were digested with EcoRI and ligated. The ligation product from BAC control was used for normalization. The relative interaction frequency was calculated as: 2^{Ct (BAC)-Ct (3C)}.

Guide RNAs were designed using the MIT CRISPR Design website (http://crispr.mit.edu). To minimize potential off-target effects of guides, only high-score guide RNAs (score >80) were used. Guide RNAs were annealed and diluted 1:200 in ddH2O and used for downstream applications.

Guide RNAs were cloned into pX330-BFP (from Dr. Tim Fenton) for upstream guides or pX458-GFP (Addgene, plasmid #48138) for downstream guides. Products were sequence verified by Sanger sequencing with the hU6-F primer. Constructs were cotransfected with Lipofectamine 3000 (Thermo Fisher, #L3000015) into 4.5 × 10⁵ cells in a 6-well plate. After 24 hours, cells were collected and the top 1% of BFP⁺/GFP⁺ cells were sorted using the SH800S Cell Sorter (Sony Biotechnology). Single cells were plated in 96-well plates and grown for 2 to 3 weeks. Single clones were screened using PCR, and a minimum of 2 homozygous clones were mixed at equal ratios and used for downstream analysis.

SF767 and HN12 cells were transduced with lenti-dCas9-KRAB-blast (from Dr. Paul Mischel) and selected with 10 μg/mL blasticidin (Gibco, #A1113903) for 72 hours after transduction. Guide RNAs were cloned into the lentiGuide-Puro vector (Addgene, plasmid #52963) and transformed into Stbl3 bacteria. Constructs were confirmed by Sanger sequencing using the hU6-F primer. Constructs were transduced individually into cells expressing dCas9-KRAB and selected with 1 μg/mL puromycin. After assessing EGFR transcript levels by RT-qPCR, one highly effective CE1 guide-expressing cell line was selected for double gRNA expression and transduced with the complementary CE2 guide.

HN12, U87, or Cal27 cells (5 × 10²) or SF767 cells (1 × 10⁴) were seeded in 5 replicate black, clear bottom 96-well plate in 6 replicate wells in complete media. After 24 hours, complete media were removed and 100 μL of 10 μg/mL blasticidin and 1 μg/mL puromycin in DMEM + 0.5% FBS was added to each well. Baseline luminescence was measured at day 1 with the ATPlite 1 step Luminescence Assay System (PerkinElmer, #6016731) on a Tecan Spark 10M. Luminescence measurements were obtained at every other day for 9 days and plotted using GraphPad Prism.

HN12 parental (1 × 10⁶; n = 5) or CE1^−/−CE2^−/− (n = 4) cells were injected into the right flank of nude mice. Tumors were measured by caliper every 4 days until visible tumors formed and then measured every day until the appearance of necrosis or the volume reached 500 mm³. The animal studies are approved by UCSD IACUC according to the NIH guidelines.

SF767 (1 × 10⁵) or HN12 (5 × 10⁵) cells were seeded in 12-well plates and grown overnight. siRNAs were transfected into each well with Lipofectamine 2000 in serum-free and antibiotic-free DMEM. Media were changed to complete media 6 hours later. Samples were collected in SDS sample buffer 48 to 72 hours later. siRNAs used for this study include BRD2 (Ambion, #s12070), BRD3 (Ambion, #s15544), BRD4 (Ambion, #s23902), and scramble control (Invitrogen, #12935-300).

SF767 or HN12 cells were treated with 0.5 μmol/L JQ1 dissolved in DMSO (MedChemExpress, #HY-13030) for 24 hours. Vehicle control samples were treated with equal volume DMSO for 24 hours.

Protein samples were collected in SDS sample buffer, separated using gel electrophoresis, and transferred via wet transfer onto a PVDF membrane. The membrane was blocked with 5% milk in TBST and probed with primary antibodies at 1:1,000 dilution overnight at 4°C and secondary HRP antibodies at 1:2,000 for 1 hour at RT. Signal was assessed via chemiluminescence with the SuperSignal West Pico PLUS substrate (Thermo Fisher, #34580) and visualized on a ChemiDoc MP system (Bio-Rad). Anti-Cas9 antibody (Cell Signaling Technology, #14697), anti-β-actin (Sigma, #A3854), anti-c-Jun (Cell Signaling Technology, #9165S), anti-c-Fos (Santa Cruz Biotechnology, #sc-52), anti-HA-HRP (Santa Cruz Biotechnology, #sc-805), anti-BRD4 (Active Motif, #39909), anti-BRD3 (Santa Cruz Biotechnology, #sc-515729), anti-BRD2 (Cell Signaling Technology, #5848S), c-Myc (Cell Signaling Technology, #9402), and anti-EGFR (BD Biosciences, #610017) were used for analysis.

Chromatin immunoprecipitation was performed as described previously (28) with the following modifications. Chromatin was sheared in diluted lysis buffer to 200 to 500 bp using a Covaris M220 Focused-Ultrasonicator with the following parameters: 10 minutes, peak incident power 75, duty factor 10%, 200 cycles/burst. Antibodies for ChIP were obtained from commercially available sources: anti-H3K27Ac (Active Motif, #39133), anti-BRD4 (Active Motif, #39909), anti-c-Jun (Cell Signaling Technology, #9165T), anti-c-Fos (Santa Cruz Biotechnology, #sc-166940), anti-BRD3 (Santa Cruz Biotechnology, #sc-515729), and anti-BRD2 (Cell Signaling Technology, #5848S). Five percent of the chromatin was not exposed to antibody and was used as control (input). For ChIP-qPCR analysis, DNA quantity for each ChIP sample was normalized against input DNA. For ChIP-seq samples, after DNA purification ChIP-seq DNA libraries were prepared with either the TruSeq ChIP Library Prep Kit (Illumina, #IP-202-1012) or the Accel-NGS 2S Plus DNA Library kit (Swift Bioscience, #21024) and sequenced using 75 bp single-end sequencing on an Illumina Hi-seq 4000.

Approximately 50,000 permeabilized nuclei were transposed using Tn5 transposase (Illumina, #FC-121-1030) as described previously (29). Libraries were amplified using NEBNext High-Fidelity 2X PCR Master Mix (NEB, #M0541) with primer extension at 72°C for 5 minutes, denaturation at 98°C for 30 seconds, followed by 8 cycles of denaturation at 98°C for 10 seconds, annealing at 63°C for 30 seconds, and extension at 72°C for 60 seconds. Each library was size selected and sequenced on an Illumina NextSeq500 or HiSeq4000 to a depth of ≥ 20 million usable reads pairs. Sequencing runs that did not meet the read pair threshold were sequenced again, and all replicates were pooled for analysis. Tracks shown in Fig. 1A and Supplementary Fig. S2 are pooled from two biological replicates.

Correlative statistical analysis was performed on the relationship between ChIP-seq read count data and relative luciferase expression, EGFR expression, or relative interaction frequency using Spearman rank order correlation. Data from chromosomal interaction, CRISPR/Cas9 deletion, cell proliferation, and ChIP-qPCR enrichment were compared, and statistical analysis was performed using a Student t test.

All DNA oligo (primer, gRNA) sequences are listed in Supplementary Table S2.

All raw and processed sequencing data generated in this study have been submitted to the NCBI Gene-Expression Omnibus (GEO; http://www.ncbi.nlm.nih.gov/geo/) under accession number GSE128275. Accession numbers for publicly available data accessed are listed in Supplementary Table S1.

To gain further insight into the mechanisms responsible for EGFR transcriptional control in HNSCC and GBM, we performed ChIP-seq and ATAC-seq in 7 GBM and 3 HNSCC cell lines with nonamplified EGFR copy numbers and varying levels of EGFR activation. Overlay of IGV tracks of all 10 cell lines showed conservation of H3K27Ac intensity and open chromatin throughout intron 1, indicating the presence of enhancers in these regions (Fig. 1A). Because superenhancers (SE) are defined as large clusters of transcriptional enhancers that drive expression of genes that define cell identity, and are often found at oncogenes (30) we identified SEs using the ROSE (rank ordering of super-enhancers) algorithm (refs. 6, 31; Supplementary Fig. S1A). Using this algorithm, we discovered cell line–specific SEs in the first intron of EGFR, many of which rank highly in several of our cell lines (Supplementary Fig. S1B). Interestingly, the location and size of these SEs varied and were dependent upon the local enrichment of H3K27Ac ChIP-seq signal.

To determine if the presence of putative enhancers in EGFR intron 1 was predictive of EGFR expression, we first measured EGFR transcript in each GBM and HNSCC cell line (Fig. 1B). Due to cell line–specific presence and location of predicted EGFR super enhancers, we used the total number of intron-1–mapped H3K27Ac ChIP-seq reads as a measure of enhancer presence and plotted these values against the fold change in EGFR expression. Analyzing the relationship by Spearman correlation showed a highly significant correlation (Fig. 1C). Together, these results identify regions containing characteristics of SEs in the first intron of EGFR and suggest that presence of these putative enhancers is important for high levels of EGFR transcript.

It has been reported previously that SEs are congregations of active constituent enhancers (CE; ref. 32). To determine which CEs of the identified SEs are active, we segmented regions that exhibited highly conserved H3K27Ac ChIP-seq and ATAC-seq signals into 2 kb segments (Supplementary Fig. S2). Each segment was then measured for enhancer activity by in vitro bioluminescence (Fig. 2A). Regions that exhibited conserved luciferase expression included 1, 3, and 16 to 19. Interestingly, activity of segments 8 and 9 was HNSCC specific, whereas segments 1 and 3 were more glioma specific. Segments 16 to 19 consistently enhanced luciferase activity in both tumor models. Combining H3K27Ac presence by ChIP-seq, open chromatin accessible regions by ATAC-seq, and functional activity as defined by our luciferase system in both tumor types, we define two distinct CEs. Specifically, we combined segments 1 to 3 into an approximately 6-kb region which we have termed constituent enhancer 1 (CE1), and combined segments 16 to 19 into an approximately 8-kb region which we have termed constituent enhancer 2 (CE2; Fig. 2A). CE1 and CE2 reside in regions of conserved high levels of H3K27Ac; however, other segments do not have associated H3K27Ac enrichment. We hypothesized that the specific enrichment of H3K27Ac at an enhancer segment would correlate to luciferase expression in the matched cell line. Indeed, plotting the normalized luciferase intensity against the average H3K27Ac read intensity (Fig. 2B and C) reveals a highly significant relationship (P < 0.001, Spearman correlation). These results define two conserved putative CE's within EGFR intron 1 and establish a relationship between histone acetylation at these regions and enhancer activity as measured by in vitro bioluminescence.

Enhancers typically make contact with one or more gene promoters through long-range interactions (30). To test if the predicted enhancer regions CE1 or CE2 interact with the EGFR promoter we performed a chromosome conformation capture (3C) assay. Primers were designed around the EcoRI sites in CE1 (F1 to F5) and CE2 (F6 to F10) and nearby the EGFR promoter (Fig. 3A). Because we hypothesized that an increased interaction frequency would correlate with increased transcript levels, we chose the cell lines with the highest EGFR expression levels in both tumor types (HN12/SF767) and compared them against a cell line with virtually no EGFR expression (TS576). Compared with TS576, in SF767 cells low levels of interaction with CE1 and CE2 were found; however, significantly stronger interactions were identified at F2 and F4 of CE1 and 4 out of 5 regions of CE2 in the HN12 cell line (Fig. 3B). These primers reside in highly acetylated regions in HN12 cells; thus, we hypothesized there is a correlation between H3K27Ac intensity and interaction frequency. Correspondingly, we identified a significant correlation between H3K27Ac peak intensity and interaction frequency when comparing all three cell lines (Fig. 3C). Together, the data thus far indicate that CE1 and CE2 have characteristics ascribing them to active enhancers: surrounded by nucleosomes with high H3K27Ac, open chromatin, transcriptional enhancement, and interaction with a promoter.

To directly assess if CE1 and CE2 were essential for EGFR expression, we used the CRISPR/Cas9 system to delete CE1 (hg 38, chr7: 55,060,994-55,066,815) and CE2 (hg38, chr7: 55,127,646-55,135,347). To test how enhancer loss effects cells with different EGFR expression levels we made these deletions in two GBM cell lines and two HNSCC cell lines which had either high (SF767, HN12) or low (U87, Cal27) relative EGFR expression. Single enhancer deletions were generated with a dual-guide deletion strategy (Fig. 4A) and the deletion was confirmed by genotyping PCR (Fig. 4B). Compared with parental cell lines, EGFR transcript and protein was significantly decreased in each deletion (Fig. 4C). In cell lines that express EGFR highly, CE1 deletion provided the strongest knockdown of EGFR transcript and protein; however, U87 had the strongest effect in the CE2 knockout cell line. Cell type–specific differences in transcript levels between CE1 and CE2 indicate there may be differential utilization of either CE1 or CE2 in different cell lines, possibly due to cell line–specific expression of TFs critical for EGFR enhancer activity.

Previous studies have indicated partial redundant control of a gene by multiple CEs within a single SE (32). To evaluate if there was a compensatory effect on EGFR transcription by either CE1 or CE2, we performed a second round of CRISPR/Cas9 editing on the homozygous edited populations (Fig. 4D). Compared with parental cell lines, EGFR transcript levels were significantly decreased with loss of both enhancers (Fig. 4E). Notably, the amplitude of EGFR transcript loss was greater in CE1^−/− + CE2^−/− when compared against CE1^−/− or CE2^−/− alone. EGFR transcript and protein loss corresponded to significant proliferation deficiencies in each cell line measured, notably most significantly in double-deleted cell lines (Fig. 4F). In HN12 cells, this proliferation difference between parental and double knockout cells translates to a significant repression of tumor growth in a subcutaneous tumor model (Fig. 4G). Restoring EGFR protein by lentiviral transduction of wild-type EGFR in CE1^−/− + CE2^−/− cells results in significant rescue of proliferation in HN12 cells (Fig. 4H). Incomplete rescue of growth in these cells indicates CE1 or CE2 may be enhancing other genes important for cell growth. Together, these results demonstrate both cell type–specific CE utilization as well as a cooperative relationship between the CE1 and CE2, whereby double enhancer deletion results in more significant deleterious effects than single deletions alone.

To eliminate the possibility of structural variation being the root cause of EGFR expression loss in CRISPR/Cas9 deleted clones, we hypothesized that histone deacetylation would be sufficient for EGFR transcriptional repression. dCas9-KRAB is known to recruit endogenous chromatin modifying complexes to deacetylate histones (33); therefore, to test our hypothesis, we targeted a nuclease deactivated Cas9 (dCas9) protein fused to the Krüppel-associated box (KRAB) domain of Kox1 (34) with four CE1 and five CE2 gRNAs in HN12 (Fig. 5A) and SF767 (Supplementary Fig. S3A) cell lines.

To ensure successful targeting of dCas9-KRAB, we confirmed deacetylation at the targeted regions. Compared with an enhancer off-target (O-T) gRNA, the targeting of dCas9-KRAB resulted in significant decreases in H3K27Ac in HN12 cells (Fig. 5B). Additionally, we measured EGFR transcript and protein levels and observed significant decreases in transcript in 80% (8/10) of enhancer-specific targeted regions (Fig. 5C). Interestingly, targeting CE2 had an overall stronger repressive effect on EGFR transcript levels, with repression by gRNA targeting CE2.4 achieving a greater than 3-fold decrease in expression. We performed double dCas9-KRAB repression by adding a second gRNA (CE1.4 + CE2.4) and observed the most significant decrease in EGFR transcript and protein levels, achieving an 8-fold decrease in RNA expression compared with the off-target gRNA (Fig. 5C). In SF767 cells, dCas9-KRAB was expressed to lower levels than that of HN12 (Supplementary Fig. S3A) and correspondingly we saw less significant reductions in H3K27Ac (Supplementary Fig. S3B). Importantly, we still saw significant repression of EGFR transcript levels in 50% (5/10) of enhancer-specific targeted regions in SF767 (Supplementary Fig. S3C). In these cells there was no observable preference for CE2; however, the strongest repression was again observed by combining CE1.4 and CE2.4.

Finally, the effect of dCas9-KRAB-mediated EGFR repression on cell proliferation was assessed by ATPlite assay. At low (0.5%) serum, the relative proliferation of the indicated EGFR-repressed cells was significantly inhibited (Fig. 5D; Supplementary Fig. S3D). Importantly, cell lines with stronger EGFR repression exhibited reduced proliferation over time. Reintroducing high levels of wild-type EGFR by lentiviral transduction resulted in a rescue of the proliferation ability of double knockdown cells. These results indicate that EGFR transcriptional changes in enhancer-deleted regions (Fig. 4) are not due solely to structural alteration within the first intron. Moreover, these results demonstrate that loss of H3K27Ac at the identified EGFR enhancers is sufficient for significant decreases in EGFR protein and transcript levels and this reduction in EGFR is sufficient to reduce cell growth.

To identify TFs important for enhancer activity, we further analyzed our H3K27Ac ChIP-seq and ATAC-seq data. To eliminate nonenhancer regulatory regions, we intersected ATAC-seq peaks with enhancer peaks from H3K27Ac and kept the TSS-distal (±2.5 kb) ATAC-seq peaks that mapped within an enhancer. Performing de novo motif analysis on these peaks in EGFR-expressing cells (SF767 and HN12) identified an AP-1 transcription factor motif as the most significantly enriched motif (Supplementary Fig. S4A).

To validate the TF motifs identified in silico, we examined AP-1 family transcription factor ChIP-seq data deposited by the Encyclopedia of DNA Elements (ENCODE) consortium (Supplementary Table S1). Using this approach, multiple c-Jun and c-Fos peaks were identified within the CE1 and CE2 of EGFR in HN12 and SF767 cells (Fig. 6A). Importantly, H3K27Ac ChIP-seq in HeLa-S3 and HUVEC also shows high levels of enhancer marks in the AP-1 marked regions (Supplementary Fig. S4B). For further analysis, we chose c-Jun and c-Fos as the prototype AP-1 factors that form heterodimers that bind to the consensus 5′-TGA(C/G)TCA motif-3′ (35). We validated and quantified c-Fos and c-Jun enrichment at the CE1 and CE2 regions in EGFR-expressing (SF767 and HN12) and nonexpressing (TS576) cells and identified significant fold enrichment of both factors in EGFR-expressing cells (Fig. 6B). Additionally, we identified significantly increased binding of c-Fos at the CE1-AP1-3 and c-Jun in CE1-AP1.3 and CE2-AP1.4 subregions in HN12 cells, indicating these regions may be important for the increased EGFR expression levels in these cells.

To validate the role of AP-1 TFs in EGFR transcription, we utilized a dominant-negative version of c-Jun (JunDN; ref. 36). JunDN can dimerize with other AP-1 family members and bind DNA; however, the transcriptional activation capability is eliminated. HN12 and U87 cells transduced with JunDN showed decreased EGFR protein levels, thus supporting a role for c-Jun heterodimers in the regulation of EGFR transcription (Fig. 6C; Supplementary Fig. S4C). Additionally, we detected a decrease in c-Jun levels when JunDN was present, likely due to autoregulation of the Jun promoter by c-Jun heterodimers (37). To confirm the effect of JunDN was due to reduced c-Jun heterodimer activity, we utilized a luciferase reporter containing a trimerized AP-1 binding motif. Luciferase reporter assays showed significant decreases in activity when JunDN is present compared with the empty vector control, confirming the reduced transcriptional activation potential of JunDN (Fig. 6D; Supplementary Fig. S4D). Taken together, these data suggest that AP-1 family members are critical for fine-tuned regulation of EGFR expression and specifically bind to EGFR enhancer regions. Perturbation of this AP-1 transactivation effect by expression of a dominant-negative results in a significant repression of c-Jun heterodimer targets including EGFR and c-Jun, confirming the role of this family of TFs in intron 1–mediated EGFR expression.

Previous research from our lab has shown treatment of mice harboring GBM neurosphere PDX models with the BET protein inhibitor JQ1 significantly prolongs survival, and combination of JQ1 with anti-EGFR therapy further increases this effect (38). To determine if this effect of JQ1 is partially attributable to downregulation of EGFR transcription through BET proteins, we treated HN12 and SF767 cells with JQ1 and measured the effects on EGFR expression. Interestingly, after 24 hours of 0.5 μmol/L JQ1 EGFR protein and transcript levels were decreased in both HN12 (Fig. 7A and B) and SF767 (Supplementary Fig. S5A and S5B) cell lines. Protein expression of BET proteins BRD2, BRD3, and BRD4 increased in response to JQ1; however, protein level of the known BRD4 target c-Myc was downregulated (Fig. 7A).

BRD4 is known to co-occupy enhancers with AP-1 family members (9, 10) and is enriched at enhancers (6), and BRD2/3 have been shown to bind to hyperacetylated regions and allow for the activity of RNA polymerase II (39). JQ1 is a pan-BET inhibitor (40); thus, to determine if JQ1 treatment was affecting EGFR levels through reduced activity of select BET family members, we first looked at expression levels of BET proteins in our cell line panel and found their expression to be highly variable but present in most cell lines (Supplementary Fig. S5C). We then looked for evidence of binding of these factors to CE1 and CE2. Recent data in the liposarcoma cell line LPS141 (41) show the presence of H3K27Ac in CE1 and CE2, and has binding of BET family members BRD2, BRD3, and BRD4 in those regions (Supplementary Fig. S5D). To confirm binding of these factors in GBM and HNSCC and to interrogate their relationship with AP-1, we performed ChIP-qPCR for c-Fos, c-Jun, BRD2, BRD3, BRD4, and H3K27Ac at regions of open chromatin in CE1 and CE2. In HN12 cells, treatment with JQ1 significantly reduces occupancy of H3K27Ac at all measured regions, and significantly reduces binding of BET and AP-1 family TFs to CE1 and CE2 (Fig. 7C). Interestingly, in contrast to steady state (Fig. 6C), which suggests CE2-AP1.4 is a critical c-Fos and c-Jun binding site, treatment with JQ1 only affects binding of BRD4 at that region (Fig. 7C). In SF767 cells, significant reductions in TF occupancy were observed primarily in CE2, with only BRD4 showing a significant reduction in binding to CE1 (Supplementary Fig. S5E). Importantly, JQ1 treatment also impedes EGFR transcription by inhibiting the interaction between CE1/CE2 and the EGFR promoter as measured by 3C (Fig. 7D; Supplementary Fig. S5F).

To determine the specific BET protein critical for the maintenance of EGFR expression, we performed siRNA knockdowns of BRD2/3/4 individually as well as in combination. In HN12 cells, single knockdown of BRD2 or BRD4, but not BRD3, resulted in downregulation of the EGFR protein. Complete loss of the BRD2/3/4 protein by combination siRNA treatment resulted in the strongest downregulation of the EGFR protein (Fig. 7D). In SF767 cells, single knockdown of BRD2 produced the strongest downregulation of EGFR protein (Supplementary Fig. S5F), indicating BRD2 and BRD4 activity on EGFR expression may be cell type specific.

To further interrogate the relationship between EGFR, AP-1, and BET clinically, we utilized H3K27Ac ChIP-seq and matched RNA-seq in 44 patient-derived glioblastoma stem cells (GSC) and 50 primary tumors (24). We organized these samples by their EGFR expression by RNA-seq and found that tumors that have high EGFR also have high levels of H3K27Ac in intron 1 (Supplementary Fig. S6A). We then plotted EGFR RNA expression versus expression of JUN, FOS, BRD2, and BRD4 and found significant correlations between expression of these TFs and EGFR levels (Supplementary Fig. S6B). Finally, to determine if the activity of the EGFR enhancers are tumor specific, we compared the expression of these genes between tumor and normal samples in the TCGA database (ref. 42; Supplementary Fig. S6C). We found in GBM tumor samples transcript levels for EGFR, FOS, JUN, BRD2, and BRD4 are significantly upregulated versus normal samples. In HNSCC tumor samples, transcript levels for EGFR, BRD2, and BRD4 are significantly upregulated versus normal samples. Taken together, these results implicate a role for BRD2 and BRD4 cooperating with AP-1, in the maintenance of EGFR expression in GBM and HNSCC.

In this study, we identify regions of epigenetic regulation within the first intron of the EGFR gene in GBM and HNSCC, characterizing DNA regions that are cell type specific in their H3K27Ac deposition, but contain conserved regions of open chromatin and histone acetylation within EGFR-expressing cells. These regions pass the threshold to be considered superenhancers and contain individual constituents that demonstrate functional attributes of active enhancers, including transcriptional enhancement in reporter assays, 3D interactions with the EGFR promoter, and negative regulation of their target gene when removed or repressed. We identify the presence and activity of AP-1 TFs in the CE1 and CE2, and when the activity of these TFs is eliminated significant effects are seen in expression of target genes including EGFR and JUN, indicating direct AP-1 dependency of these genes. Pharmacologic disruption of the transcription factor complexes at these enhancers has significant effects on EGFR expression, providing a mechanism by which this transcriptional control mechanism may be targeted. Although identification of EGFR superenhancers in other tumor types has been published previously (8, 10, 24–26), this study provides significant advances in our understanding of the most critical EGFR enhancer regions by undertaking functional genomic and pharmacologic approaches to directly perturb their activity and demonstrating functional consequences on tumor cell proliferation both in vitro and in vivo.

Elevated EGFR is a well-established therapeutic target; however, responses to EGFR tyrosine kinase inhibitors (TKI) are sporadic. The mechanisms behind resistance to EGFR TKI are unique to each tumor type, and secondary resistance mutations are common. In HNSCC, high EGFR copy numbers are statistically associated with cetuximab and gefitinib resistance (43), and although rare, kinase domain mutations may be associated with altered responses to EGFR inhibitors (44). In GBM, resistance mechanisms are less well understood with prevailing theories including intratumoral heterogeneity (45), EGFR-amplified extrachromosomal DNA (ecDNA; ref. 46), and loss of PTEN (47). Another tumor type that has a high prevalence of EGFR overexpression is non–small cell lung cancer (NSCLC), the resistance mechanisms for which are largely kinase domain mutations (48, 49). With these challenges in mind, this study presents a kinase domain–independent mechanism by which EGFR expression and activity can be prevented. Recent studies have shown targeted transcription factor blockade can overcome EGFR TKI resistance (50); thus, this study presents an additional pathway that can be targeted alone or in combination with EGFR TKIs to treat EGFR-positive tumors. Further studies should interrogate the CE1 and CE2 in EGFR-mutated models, specifically in NSCLC, where TKI-resistant EGFR mutations are present in as much as 63% of tumors (48).

This study focuses exclusively on nonamplified EGFR; however, significant fractions of both HNSCC and GBM tumors have high copy numbers of the gene (16, 51). In GBM, EGFR can be amplified both as a homogeneously staining region (HSR) or on ecDNA (46), both of which include the entire EGFR gene and surrounding regions. Amplified enhancer regions maintain their enhancer signatures (52) and focal amplifications of superenhancers have been shown to drive transcription (53). Because amplified enhancers have standard enhancer features, we hypothesize that targeting amplified EGFR directly (Fig. 5) or pharmacologically (Fig. 7) will have similar inhibitory effects. Future studies should explore this hypothesis to further broaden the clinical implications of targeting EGFR intron 1 enhancers.

In tumor types that frequently alter EGFR, structural variations (SVs) are common. One such SV, EGFR variant III (EGFRvIII), is an extracellular domain mutation that shows constitutive tyrosine kinase activity (54). The mutant receptor is highly tumor specific (55) and relatively common in GBM (51). EGFRvIII protein has been identified in other tumor types, including HNSCC (56), although genomic deletion of exons 2–7 has not been detected. In GBM, every incidence of EGFRvIII has unique genomic breakpoints within intron 1 and intron 8, removing exons 2–7 and large portions of intron 1 depending on breakpoint location, often including regions that we have identified as enhancers (ref. 57; Fig. 1A). It is important to note the majority of identified EGFRvIII breakpoints occur nearer the 3′ end of EGFR intron 1, often removing the region that we have identified as CE2 (57). Interestingly, in some of our cell line models, the data indicate a more significant role for CE2 in interaction with the EGFR promoter by 3C (Fig. 3) and influence on the expression of EGFR by dCas9-KRAB (Fig. 5). In other models, the utilization of each CE is either equal (Supplementary Fig. S3) or demonstrates a CE1 bias (Fig. 4C, HN12/SF767/Cal27). But, targeting both CEs always produces a combinatorial effect (Figs. 4 and 5). These data demonstrate an intricate regulatory system that is cell type–specific and makes predicting the functional consequences of intron 1-loss on EGFRvIII expression difficult. Further research should interrogate directly how the extent of intron 1-loss affects EGFRvIII through analysis of patient-derived and/or CRISPR/Cas9-generated EGFRvIII⁺ models.

Many potent oncogenes are associated with superenhancers, including MYC and HER2 (30, 58). These genes and others have been successfully targeted in cancers that express them by inhibiting the BET protein BRD4, a hallmark factor involved in super enhancer identity (6). The BET bromodomain inhibitor JQ1 has been shown to inhibit BRD4 at superenhancers (6) and additionally sensitizes EGFR-amplified GBM (38) and HNSCC (59) cells to EGFR TKI. Our data suggest that this sensitization may also be due in part to JQ1-mediated inhibition of AP-1 and other BET proteins at EGFR intron 1 enhancers (Fig. 7). Although site-specific inhibition of TF binding is cell type dependent, whether targeting TF binding pharmacologically (Fig. 7A–D; Supplementary Fig. S5A–S5G) or through RNAi (Fig. 7E; Supplementary Fig. S5F), the effects on EGFR are consistently deleterious. These results support the global targeting of AP-1 and BET rather than site-specific repression (e.g., dCas9-KRAB). Additionally, our study shows significant reduction of EGFR protein and transcript after JQ1 treatment independent of gene expression (Fig. 7A and B; Supplementary Fig. S5A and S5B), indicating a broader effect of JQ1 on tumor models that express EGFR at varying levels. Indeed, EGFR transcript levels are reduced to similar levels in SF767 cells, which express less EGFR and exhibit less significant TF occupancy differences in response to JQ1 treatment (Supplementary Fig. S5E). Further, clinical benefit for patients with EGFR-expressing GBMs is supported by the presence of superenhancers in EGFR-positive GSC and primary tumor samples (Supplementary Fig. S6A). These tumors exhibit positive correlations between critical TF expression and EGFR mRNA, and suggest limited off-target effects in the brain due to significant differences in TF expression between normal and tumor tissues (Supplementary Fig. S6B and S6C). These data combined further support the combination of EGFR TKI and JQ1 as treatment for EGFR-positive malignancies.

In conclusion, we found that EGFR expression is maintained in part through presence and activity of critical enhancers present in intron 1 of the gene. Characterization of CE1 and CE2 in multiple cell line systems identified a novel role for BET transcriptional coactivators and AP-1 TFs in these enhancers, and provided the rationale for therapeutic targeting of EGFR through perturbation of BET and AP-1 in EGFR-positive malignancies.

No potential conflicts of interest were disclosed.

Conception and design: N.M. Jameson, F. Furnari

Development of methodology: N.M. Jameson, J. Benitez

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): N.M. Jameson, J. Ma, A. Izurieta, J.Y. Han, R. Mendez

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): N.M. Jameson, A. Parisian, F. Furnari

Writing, review, and/or revision of the manuscript: N.M. Jameson, J. Ma, J. Benitez, A. Izurieta, A. Parisian, F. Furnari

Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): N.M. Jameson

Study supervision: N.M. Jameson, J. Benitez, F. Furnari

We thank Dr. Bing Ren for sharing the ChIP protocol; Dr. Paul Mischel and Dr. Sudhir Chowdry for sharing the dCas9-KRAB-Blast and lentiGuide-Puro vectors; Dr. David Gorkin for advising on the ATAC-seq data generation and analysis; Dr. Silvio Gutkind for sharing the HNSCC cell lines; Dr. Roman Sasik for providing the ChIP-seq analysis pipeline.

This work was supported by grants from the NIH (NS080939; to F. Furnari) and the UCSD Cancer Biology, Informatics & Omics Training Program (T32-CA067754-22; to N.M. Jameson).

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.