A major hurdle to the application of precision oncology in pancreatic cancer is the lack of molecular stratification approaches and targeted therapy for defined molecular subtypes. In this work, we sought to gain further insight and identify molecular and epigenetic signatures of the Basal-like A pancreatic ductal adenocarcinoma (PDAC) subgroup that can be applied to clinical samples for patient stratification and/or therapy monitoring. We generated and integrated global gene expression and epigenome mapping data from patient-derived xenograft models to identify subtype-specific enhancer regions that were validated in patient-derived samples. In addition, complementary nascent transcription and chromatin topology (HiChIP) analyses revealed a Basal-like A subtype-specific transcribed enhancer program in PDAC characterized by enhancer RNA (eRNA) production that is associated with more frequent chromatin interactions and subtype-specific gene activation. Importantly, we successfully confirmed the validity of eRNA detection as a possible histologic approach for PDAC patient stratification by performing RNA-ISH analyses for subtype-specific eRNAs on pathologic tissue samples. Thus, this study provides proof-of-concept that subtype-specific epigenetic changes relevant for PDAC progression can be detected at a single-cell level in complex, heterogeneous, primary tumor material.

Implications:

Subtype-specific enhancer activity analysis via detection of eRNAs on a single-cell level in patient material can be used as a potential tool for treatment stratification.

This article is featured in Selected Articles from This Issue, p. 865

Pancreatic ductal adenocarcinoma (PDAC) is the most common form of pancreatic cancer with a complex molecular landscape and highly heterogeneous nature. Despite the decline in overall cancer death rates, the mortality rate of pancreatic cancer remains high due to the lack of proper early detection methods and effective therapeutic options beyond surgery. Consequently, pancreatic cancer is projected to become the second leading cause of cancer-related deaths by 2030 (1). Accordingly, there is an urgent need to elucidate the molecular mechanisms underlying tumor development, therapeutic resistance, and metastasis of this lethal disease.

A number of large-scale studies involving transcriptomic, genetic and epigenetic analyses of PDAC tumor samples (2–6) consistently identified two principally distinct molecular subtypes referred to as “classical-” and “Basal-like”, which display differences in phenotypic behavior, outcome and therapeutic responsiveness. Subsequently, a more detailed transcriptomic analysis including samples from advanced-stage PDAC tumors (7) revealed that the previously identified major subtypes can be further segregated into two subclusters each based on distinct gene expression signatures. Among the identified subtypes, Basal-like A PDAC was shown to associate with particularly aggressive and chemoresistant behavior. While the origin of Basal-like cell identity in PDAC is still poorly understood, it has been closely linked to the expression of the the ΔN isoform of the transcription factor TP63 (ΔNp63) (3, 8, 9).

ΔNp63 is a master regulator that plays an essential role in the development of stratified epithelial tissues such as epidermis, breast and prostate (10). A number of studies have revealed the oncogenic potential of ΔNp63, which was found to be specifically overexpressed in the “Basal” or “squamous” subtypes of diverse epithelial cancers that display a more aggressive clinical course and poor prognosis (11, 12). While ΔNp63 expression was previously linked to Basal-like PDAC (3, 13, 14), its functional relevance in PDAC was unclear until recently. ΔNp63 directs subtype-specific gene expression in PDAC through the dynamic epigenetic regulation of enhancer landscapes (8, 9, 15) and promotes inflammatory changes in the tumor microenvironment (16).

As a lineage-defining TF, ΔNp63 occupies a large number of binding sites at distal enhancers compared with promoters (8, 9). Enhancers are distal regulatory elements which have a prominent role in driving the expression of lineage-specific genes upon binding of TFs and cofactors, which in turn mediate long-range interactions with target gene promoters (17). Growing evidence suggests that tumor-associated molecular aberrations result in a rewiring of the enhancer networks, thereby directing cancer cell identity (18, 19). Notably, in addition to epigenetic marking, some enhancers are transcribed and produce short-lived noncoding RNA molecules referred as enhancer RNAs (eRNA; refs. 20, 21). Although the functional relevance of eRNAs remains to be definitively established, multiple studies have uncovered an activating role of eRNAs in transcriptional regulation (22–26). Spatio-temporal control of gene expression is achieved by the dynamic interaction of distal enhancers with the target promoters. Interestingly, eRNAs have been reported to contribute to the maintenance and stabilization of enhancer-promoter loop structures (18, 27, 28). In fact, eRNAs can directly associate with TFs and stabilize their binding to the cognate enhancers (27–30). Hence, based on the current body of evidence, we hypothesized that eRNAs can serve as a precise readout of cell-type specific enhancer activity. Consistently, enhancer activity profiling across breast cancer cell lines was used to identify the TFs driving subtype-specific biology of breast cancer (31).

In this study, we performed comprehensive multi-omics analyses to characterize the molecular determinants of ΔNp63-dependent transcription in PDAC. For the first time, we report a Basal-like A subtype-specific transcribed enhancer program (B-STEP) which is essential for subtype-specific gene expression in Basal-like A PDAC. Moreover, utilization of a CRISPR-dCas9-KRAB epigenome editing approach confirmed the essential role of several B-STEP enhancers. Importantly, we show that the B-STEP program and in vivo growth of Basal-like patient-derived xenograft (PDX) can be inhibited by pharmacologic targeting of Bromodomain- and Extra-Terminal (BET) domain proteins. Finally, we validated the detectability of B-STEP eRNAs in formalin-fixed paraffin-embedded (FFPE) tissue samples from xenograft and patient samples. Overall, this study uncovers a unique repertoire of transcribed, subtype-specific enhancers in PDAC and displays the potential utility of their detection for clinical stratification and/or accurate monitoring of therapeutic response.

PDX

PDX models generated in Göttingen, Germany (GöPDX) were prepared as previously described (32). Informed and written consent was obtained from all patients. Most PDX samples used for mRNA sequencing (mRNA-seq) were F1 generation with exception of GöPDX5 (F4), GöPDX12 (F2), and GöPDX13 (F3). Animal experiments were performed following the protocols approved by the Institutional Animal Care and Use Committee (IACUC; 33.9–42502–04–17/2407). The generation and utilization of PDX models was approved by the ethical review board of the University Medical Center Göttingen (UMG; 70112108).

Mouse xenografts

L3.6pl xenografts were generated during a previous study (33). HT29 xenografts were generated by injecting one million cells under sterile conditions subcutaneously into the right abdominal flank of SCID mice. Tumors were allowed to grow for 21 days. Mouse weight and tumor size were measured daily. Animal experiments were performed following the protocols approved by the IACUC (33.9–42502–04–15/2039).

Cell culture

L3.6pl (RRID:CVCL_0384; ref. 34) were kindly provided by Elisabeth Hessmann (UMG, Germany) and AsPC1 cells (RRID:CVCL_0152) – by Frauke Alves (35) [Max-Planck-Institute for Multidisciplinary Scienes Göttingen, Germany]. L3.6pl cells were maintained and cultured in Eagle Minimum Essential Medium (Gibco, Life Technologies) and AsPC1 cells - in RPMI medium (Gibco, Life Technologies). Cells were regularly tested for Mycoplasma using the MycoAlert Mycoplasma Detection Kit (Lonza, #LT07–318). Cells were maintained in culture for a maximum of 10 to 15 passages on average before thawing new stock. Although cell line authentication was not performed, cell lines were regularly monitored for phenotypic appearance and behavior. Transfection of siRNAs was achieved using Lipofectamine RNAiMAX (Invitrogen) according to the manufacturer's instructions. Control non-targeting siRNA #5 (NT5) (siGENOME, #D-001210–05–20, target sequence: UGGUUUACAUGUCGACUAA) and TP63 siRNA pool (siGENOME; #M-003330–01–0005 SMARTpool consisting of D-003330–05 CAUCAUGUCUGGACUAUUU; D-003330–06 CAAACAAGAUUGAGAUUAG; D-003330–07 GCACACAGACAAAUGAAUU; D-003330–08 CGACAGUCUUGUACAAUUU) were purchased from Dharmacon, ThermoScientific. For JQ1 experiments, 24 hours after seeding, cells were treated with 500 nmol/L JQ1 (dissolved in DMSO) or an equivalent volume of DMSO for 24 hours.

For proliferation assay, 6,000 cells were seeded in each well of a 48-well plate and left to grow overnight. Next day, cells were treated with different concentrations of pan BET inhibitor CC-90010 (Trotabresib, BMS-986378) (MedChemExpress) and transferred to Incucyte S3 Live Cell Analysis Instrument (Sartorius; RRID:SCR_023147) for live cell imaging. Cell growth inhibition was assessed for total duration of 96 hours while treatment was replenished every 48 hours. Confluence was measured by Incucyte 2022B Rev2 software. IC50 was calculated from dose-response curves (three parameters) using GraphPad Prism 9.5.1 (RRID:SCR_002798).

Establishment of stable dCas9-KRAB expressing cell lines

HEK293T (RRID:CVCL_0063) cells were transfected with pMD2.G (Addgene, #12259, RRID:Addgene_12259), psPAX2 (Addgene, #12260, RRID:Addgene_12260), and Lenti-dCas9-KRAB-blast (Addgene, #89567, RRID:Addgene_89567) using Lipofectamine 3000 (Invitrogen) according to the manufacturer's instructions. Viral supernatant was pooled by collecting the media at 24, 48, and 72 hours post-transfection and centrifuged at 500 × g for 10 minutes to remove any cellular material. 500,000 L3.6pl cells were transduced by spin-oculation at 1,000 × g for 2 hours at 30°C using 8 μg/mL polybrene in a dilution of viral supernatant and normal culture media. After 48 hours, the cells were resuspended in appropriate culture medium supplemented with 10 μg/mL Blasticidin S (EMD Millipore Corp). Selection was maintained throughout all experiments.

For epigenome editing, 400,000 L3.6pl-dCas9-KRAB cells per well (in a 6-well plate) were reverse transfected using Lipofectamine RNAiMax (Invitrogen) with custom designed synthetic guide-RNAs (Integrated DNA Technologies) targeting individual super-enhancers (SE) or the SNAI2 transcriptional start site (Supplementary Table S1) following the manufacturer's instructions. RNA was harvested 48 hours post-transfection and gene expression was assessed via qRT-PCR.

IHC staining

IHC staining was performed in the Clinical Pathology Core of the Mayo Clinic (Rochester, MN) with the Leica Bond platform (Leica) using antibodies against p40 (ΔNp63, Biocare, #ACI 3066 A, C; clone: BC28; RRID:AB_2858274), HNF1B (Sigma, #HPA002083, Polyclonal; RRID:AB_1080232), KRT5 (Leica, #NCL-L-CK5; clone: XM26) and KRT7 (Dako, # M7018; clone: OV-TL 12/30) after heat-induced antigen retrieval. Studies on human samples were approved by the Institutional Review Board of the Mayo Clinic (#19–005710).

RNAScope

RNA-ISH staining was performed at the Mayo Clinic Pathology Research Core (Rochester, MN) using the Leica Bond RX stainer (Leica). RNA-ISH staining was performed on-line. Slides were dewaxed and baked prior to being exposed to the retrieval solution, Epitope Retrieval 2 (EDTA; Leica) for 15 minutes at 95°C and the Protease III pretreatment (ACD) for 15 minutes at 40°C. The slides were incubated for 2 hours at 40°C with either a positive control probe (Hs-PPIB probe, ACD), negative control (DapB, ACD) or custom designed probes against the eRNA-producing regions contained within the SEs associated with SNAI2 (Hs-eSNAI2–2, ACD) or TP63 (Hs-eTPRG1, ACD). After hybridization, the slides were incubated in a series of amplification reagents (RNAscope 2.5 LSx Reagent Kit-Brown, ACD) according to the manufacturer's protocol. Slides were incubated for 20 minutes in DAB for visualization and counterstained in hematoxylin for 5 minutes, followed by a bluing agent for 2 minutes. Slides were removed from the stainer and rinsed in tap water, dehydrated in increasing concentrations of ethanol, and cleared 3 times in xylene prior to mounting with coverslips in xylene-based medium.

Chromatin immunoprecipitation

Chromatin immunoprecipitation (ChIP) was performed as previously described (32). PDX tissue samples were homogenized on dry ice using the Sciences BioMascher II Micro Tissue homogenizer. Powdered tissue and cell pellets from L3.6pl cells were fixed in 1% formaldehyde in PBS for 20 minutes and quenched with 125 mmol/L glycine for 5 minutes at room temperature. Sonicated chromatin extract was immunoprecipitated overnight with 1 μg H3K27ac antibody (Diagenode, #C15410196, RRID:AB_2637079) or 1 μg control rabbit IgG (Diagenode, #C15410206, RRID:AB_2722554). Protein A sepharose was used to capture immune complexes. The ChIP efficiency was verified via qRT-PCR and normalized to input DNA.

Length extension chromatin run-on followed by sequencing

A detailed protocol for length extension chromatin run-on followed by sequencing (leChRO-seq) was kindly provided by Charles G. Danko and was performed as previously described (36).

H3K27ac HiChIP

HiChIP was performed as previously described (37). Cells were crosslinked with 1% formaldehyde. For IP, 6 μg of H3K27ac antibody Diagenode, #C15410206, RRID:AB_2722554) were used. Sequencing libraries were prepared using the KAPA HyperPrep (Roche) kit. The libraries were sequenced by the Genome Analysis Core at Mayo Clinic, Minnesota.

Statistical analysis

For gene expression analyses, statistical comparisons were performed using Student t test or one-way ANOVA followed by Dunnett post test in GraphPad Prism (GraphPad Software, RRID:SCR_002798). Individual P values and type of statistical test are indicated in the figure legends. Data is represented by the mean and standard deviation.

Data availability

The sequencing data generated in this study have been deposited into the Gene Expression Omnibus (RRID:SCR_005012) under GSE213397 superseries. The additional data used in this study were obtained from ArrayExpress (http://www.ebi.ac.uk/arrayexpress) under accession codes E-MTAB-5632 (PDX H3K27ac ChIP-seq), E-MTAB-5639 (PDX mRNA-seq; ref. 38); E-MTAB-7034 (L3.6pl TP63 & BRD4 ChIP-seq; ref. 9).

Further methodologic details including tumor growth analysis following JQ1 treatments of PDX-bearing mice, library preparation and next-generation sequencing and data analysis can be found in Supplementary Data.

PDX models faithfully recapitulate PDAC subtype-specific gene expression signatures

A recent study showed that the two major PDAC molecular subtypes “Basal-like” and “classical” PDAC can be further classified into further subgroups based on distinct gene expression patterns (7). To determine whether these molecular subtypes are faithfully reflected in experimental model systems, we established an independent cohort of PDX samples from tumor specimens resected at the University Medical Center Göttingen, Germany (GöPDX). Application of signature gene analysis with subsequent clustering of these specimens identified two PDX samples (GöPDX4 and GöPDX15) clearly classifiable as Basal-like A, expressing higher levels of “Signature 2” genes (7), which define the most aggressive Basal-like A PDAC subtype in patient samples (Fig. 1A, red box). Of the 13 samples analyzed, eight (GöPDX5, 8, 9, 11, 12, 19, 21, 23) clearly expressed “Signature 1 and 6” genes and were hence classified as “classical A/B” (Fig. 1A, blue box). We verified subtype identity by examining the expression of several signature genes in PDX samples classified as Basal-like A (GöPDX4 and GöPDX15) and “classical A/B” (GöPDX5 and GöPDX8; Supplementary Fig. S1A). Importantly, IHC analyses of the Basal, Signature 2 marker KRT5 and the classical, Signature 6 marker HNF1B conducted in primary PDAC tissue from the donor patients matching these PDX tumors supported the respective molecular subtype (Fig. 1B). Hence, these PDX models indeed recapitulate molecular features of the original tumors and thus represent a powerful preclinical tool to understand the underlying molecular biology and targeting thereof.

Figure 1.

ΔNp63 defines the Basal-like A-specific enhancer landscape in PDAC. A, Hierarchical clustering of PDX samples from patient cohort generated in Göttingen, Germany using the gene signatures (1,2,6 & 10) from Chan-Seng-Yue M and colleagues's study. The heat maps show the raw z-score values generated on log2 transformed normalized counts from mRNA-seq. B, Primary tumor samples from patients GöPat4, GöPat15, GöPat5, and GöPat8 were subjected to IHC detection of KRT5, KRT7, and HNF1B protein expression. Scale bar = 50 μm. C, Heat maps depict ΔNp63 binding (left) and H3K27ac (right) occupancy (normalized to 1x) in Basal-like A GöPDX4 centered ±5 KB around the ΔNp63 summits. D, Ranking of enhancers based on H3K27ac signal intensity in GöPDX4 generated by ROSE algorithm. In total 551 SEs were identified. E, 551 GöPDX4 –specific SE identified in (D) were overlapped with ΔNp63 summits to evaluate the distribution of ΔNp63 binding in SE regions. Then, 338 ΔNp63 bound SE were used in downstream analysis to identify the closest nearby genes. Of these 288 genes, those associated with “Signature 2” are shown. F, IGV tracks show normalized ΔNp63 and H3K27ac occupancy in reads per genomic content (RGC) over the 2 SE regions indicated in (D) in GöPDX4 as well as the H3K27ac occupancy in Basal-like A PDX samples 1.033, AO.IP and foie8 and “classical A/B” PDX samples 1.064, 2.045 and 2.116 from Lomberk G. et al. study.

Figure 1.

ΔNp63 defines the Basal-like A-specific enhancer landscape in PDAC. A, Hierarchical clustering of PDX samples from patient cohort generated in Göttingen, Germany using the gene signatures (1,2,6 & 10) from Chan-Seng-Yue M and colleagues's study. The heat maps show the raw z-score values generated on log2 transformed normalized counts from mRNA-seq. B, Primary tumor samples from patients GöPat4, GöPat15, GöPat5, and GöPat8 were subjected to IHC detection of KRT5, KRT7, and HNF1B protein expression. Scale bar = 50 μm. C, Heat maps depict ΔNp63 binding (left) and H3K27ac (right) occupancy (normalized to 1x) in Basal-like A GöPDX4 centered ±5 KB around the ΔNp63 summits. D, Ranking of enhancers based on H3K27ac signal intensity in GöPDX4 generated by ROSE algorithm. In total 551 SEs were identified. E, 551 GöPDX4 –specific SE identified in (D) were overlapped with ΔNp63 summits to evaluate the distribution of ΔNp63 binding in SE regions. Then, 338 ΔNp63 bound SE were used in downstream analysis to identify the closest nearby genes. Of these 288 genes, those associated with “Signature 2” are shown. F, IGV tracks show normalized ΔNp63 and H3K27ac occupancy in reads per genomic content (RGC) over the 2 SE regions indicated in (D) in GöPDX4 as well as the H3K27ac occupancy in Basal-like A PDX samples 1.033, AO.IP and foie8 and “classical A/B” PDX samples 1.064, 2.045 and 2.116 from Lomberk G. et al. study.

Close modal

ΔNp63 defines the Basal-like A specific enhancer program

Subtype-specific gene expression patterns in cancer are controlled by the epigenetic enhancer landscape and the tissue-specific TFs that define them (39). The oncogenic function of ΔNp63 in PDAC has been linked with the activation of a “squamous” gene expression program that defines the aggressive clinical course of Basal-like PDAC (7–9). To further explore the direct role of ΔNp63 in preclinical PDX models, we performed ChIP-seq analysis for ΔNp63 in a Basal-like A PDX (GöPDX4). Consistent with a role in gene activation, ΔNp63-occupied genomic regions in GöPDX4 were marked by the active epigenetic mark H3K27ac (Fig. 1C). As a key driver of cell fate specification, ΔNp63 exerts its function via establishment of squamous-specific enhancer clusters, i.e., SEs, which have been strongly implicated in oncogenic transcriptional pattern (40). In fact, ΔNp63 displayed around 18% promoter-specific binding pattern, whereas around 40% of ΔNp63 occupancy occurred at distal intergenic regions (Supplementary Fig. S1B). Consistently, of 551 SEs identified in GöPDX4 (Fig. 1D), ca. 60% were occupied by ΔNp63 (Fig. 1E; Supplementary Table S2). Moreover, among the genes associated with ΔNp63-bound SEs (Supplementary Table S3) 23 belonged to Basal-like A “Signature 2” as indicated in (Fig. 1E). To verify the “Basal-like” nature of identified SEs we used data from an independent published cohort of PDX PDAC samples (38). Consistent with previous classification, three PDX samples (1.033; AO.IP; foie8) that were characterized as squamous based on ΔNp63-driven epigenetic and gene expression patterns (9), displayed a particularly high expression of genes from “Signature 2” genes (Supplementary Fig. S1C, red box). The remaining three PDX tumors classified as Basal in the initial publication (37) expressed genes of both “Basal” and “classical” signatures to differing degrees, probably more closely representing the “hybrid” subtype of PDAC. Nine PDX specimens consistently expressed “classical A” and “classical B” genes, with no clear separation into either of the two classical subclusters. Thus, we classified these as “classical A/B” (Supplementary Fig. S1C, blue box). Consistent with a key role in cellular identity specification, the ΔNp63-enriched SEs from GöPDX4 displayed a subtype-specific activation pattern that was not marked by H3K27ac in the classical PDX samples (Fig.1F; Supplementary Fig. S1D). In summary, the identification of Basal-like A-specific enhancers supports a role for ΔNp63 in shaping the enhancer architecture of the most aggressive PDAC subtype.

ΔNp63 activates a transcribed enhancer program to maintain oncogenic transcription

To further delineate the underlying mechanisms of ΔNp63-dependent enhancer activation in shaping the Basal-like A specific disease phenotype, we sought to use a cell culture model where ΔNp63 expression could be modulated. Consistent with our previous identification of ΔNp63 as a central component defining L3.6pl cellular identity (9), an unbiased principal component analysis of all H3K27ac-enriched regions revealed that L3.6pl cells and GöPDX4 were both similar to the Basal-like A PDX samples from Lomberk, and colleagues (Fig. 2A). Interestingly, upon loss of ΔNp63 in L3.6pl cells (Fig. 2B; Supplementary Fig. S2A) whereas some of “Signature 2” genes were downregulated, the expression of the three other signatures was mainly enriched (Fig. 2C; Supplementary Fig. S2B), thereby confirming the central role of ΔNp63 in driving the squamous gene expression program and suppressing other subtype identities.

Figure 2.

ΔNp63 defines the Basal-like A-specific oncogenic expression via transcribed enhancer program. A, principal component analysis (PCA) of H3K27ac occupancy in L3.6pl, Basal-like A PDX samples 1.033, AO.IP, foie8 and GöPDX4 as well as “classical A/B” PDX samples 1.064, 2.045, and 2.116 on all H3K27ac regions. B, The volcano plot shows ΔNp63-regulated genes following TP63 depletion in L3.6pl cells. Red color indicates significantly up- or downregulated genes with abs log2fold change |L2FC| ≥ 1 and Padj ≤ 0.05; the blue dots – significantly (Padj ≤ 0.05) regulated genes with less than abs|L2FC| < 1; the green dots – insignificantly regulated genes with abs|L2FC| ≥ 1 and Padj > 0.05 and the black dots – unregulated genes with Padj > 0.05 and abs|L2FC| < 1. C, Gene set enrichment analysis was performed to compare the enrichment of the “Signature 2, 10,1 and 6” in L3.6pl cells upon TP63 knockdown. Normalized enrichment score (NES) and FDRs are shown within each plot. D, The schematic explaination how putative ΔNp63-bound enhancers in L3.6pl were defined. These regions were subjected to DiffBind analysis using H3K727ac signal change upon TP63 depletion in L3.6pl cells. The MA plots shows the change of H3K27ac occupancy on these enhancer regions with red points representing differentially bound sites. In total 718 ΔNp63-bound enhancer regions show significant change (FDR < 0.05) in H3K27ac signal upon TP63 knockdown, namely 601 sites–losing the signal (DN) and 117 gaining it (UP). E, Genomic distribution of ΔNp63-dependent enhancers identified in (D). F, Heat maps depict ΔNp63, BRD4, H3K4me3, and H3K4me1 occupancy (normalized to 1x) in control as well as H3K27ac occupancy in control and TP63-depleted L3.6pl cells on ΔNp63-dependent enhancer regions centered ±5 KB around ΔNp63 summits. G, Average leChRO-seq signal normalized to counts per million (CPM) depicted at ΔNp63-dependent enhancer regions centered ±1 KB around ΔNp63 summits in control (NT5) and TP63-depleted (siTP63) L3.6pl cells.

Figure 2.

ΔNp63 defines the Basal-like A-specific oncogenic expression via transcribed enhancer program. A, principal component analysis (PCA) of H3K27ac occupancy in L3.6pl, Basal-like A PDX samples 1.033, AO.IP, foie8 and GöPDX4 as well as “classical A/B” PDX samples 1.064, 2.045, and 2.116 on all H3K27ac regions. B, The volcano plot shows ΔNp63-regulated genes following TP63 depletion in L3.6pl cells. Red color indicates significantly up- or downregulated genes with abs log2fold change |L2FC| ≥ 1 and Padj ≤ 0.05; the blue dots – significantly (Padj ≤ 0.05) regulated genes with less than abs|L2FC| < 1; the green dots – insignificantly regulated genes with abs|L2FC| ≥ 1 and Padj > 0.05 and the black dots – unregulated genes with Padj > 0.05 and abs|L2FC| < 1. C, Gene set enrichment analysis was performed to compare the enrichment of the “Signature 2, 10,1 and 6” in L3.6pl cells upon TP63 knockdown. Normalized enrichment score (NES) and FDRs are shown within each plot. D, The schematic explaination how putative ΔNp63-bound enhancers in L3.6pl were defined. These regions were subjected to DiffBind analysis using H3K727ac signal change upon TP63 depletion in L3.6pl cells. The MA plots shows the change of H3K27ac occupancy on these enhancer regions with red points representing differentially bound sites. In total 718 ΔNp63-bound enhancer regions show significant change (FDR < 0.05) in H3K27ac signal upon TP63 knockdown, namely 601 sites–losing the signal (DN) and 117 gaining it (UP). E, Genomic distribution of ΔNp63-dependent enhancers identified in (D). F, Heat maps depict ΔNp63, BRD4, H3K4me3, and H3K4me1 occupancy (normalized to 1x) in control as well as H3K27ac occupancy in control and TP63-depleted L3.6pl cells on ΔNp63-dependent enhancer regions centered ±5 KB around ΔNp63 summits. G, Average leChRO-seq signal normalized to counts per million (CPM) depicted at ΔNp63-dependent enhancer regions centered ±1 KB around ΔNp63 summits in control (NT5) and TP63-depleted (siTP63) L3.6pl cells.

Close modal

After validating the Basal-like A identity of L3.6pl cells and their high degree of similarity to patient-derived tumor samples, we sought to identify the ΔNp63-depedent enhancers that control the identity of these cells. To specifically elucidate the function of ΔNp63-bound distal regulatory elements, we focused on distal ΔNp63 peaks not contained within the transcribed region of any active genes (Fig. 2D) and performed differential H3K27ac occupancy analysis on these upon ΔNp63 depletion. In total, 718 ΔNp63-bound H3K27ac regions displayed significant change, with the 83% of these sites losing H3K27ac occupancy upon ΔNp63 depletion (Fig. 2D and E), supporting the key function of ΔNp63 in establishing the epigenetic status of the distal regulatory sites. Moreover, these ΔNp63-dependent distal sites (Supplementary Table S4) displayed typical active enhancer features characterized by BRD4 and H3K4me1 occupancy (Fig. 2F).

A number of recent studies have revealed that enhancer activity strongly correlates with noncoding eRNA production (41). To identify to what degree the ΔNp63-dependent enhancer repertoire was associated with active transcription, we performed leChRO-seq (36) that allows mapping of RNAPII-driven nascent transcription. Whereas nascent transcription remained unaffected on ΔNp63-independent H3K27ac bound regions (Supplementary Fig. S2C and S2D), bidirectional RNA transcription was substantially reduced on ΔNp63-dependent enhancer regions following loss of ΔNp63 (Fig. 2G), thereby further supporting a central regulatory role of ΔNp63 in controlling a Basal-like A transcribed enhancer program.

Transcribed enhancers display frequent higher order chromatin interactions with target genes

Enhancers regulate gene expression via long-range chromosomal interactions with target gene promoters, and eRNAs have been shown to contribute to stabilization of these contacts (27, 28). To identify genes directly associated with ΔNp63-dependent enhancers we performed HiChIP in L3.6pl cells, using H3K27ac to identify genomic interactions with active enhancer regions. In total, 626 genes were found to engage with ΔNp63-dependent enhancers (Fig. 3A; Supplementary Table S5) including 76 ΔNp63-dependent genes. In a recent work we identified a subgroup of highly interacting enhancers associated with chemoresistance in PDAC that displayed particularly high rates of eRNA production (42). In addition, high transcriptional activity was found to be associated with subset of enhancers forming dense interactome with genes defining the cell identity (43). Thus, we hypothesized that ΔNp63-controlled eRNA transcription may be associated with enhancer-gene interaction frequency as well. We classified ΔNp63-dependent enhancers into 5 groups based on the number of interactions originating from them (Fig. 3B) and observed that eRNA production positively correlated with interaction frequency (Fig. 3C). This finding is consistent with a potential role of eRNAs and/or the active RNA Polymerase II in mediating and/or stabilizing the enhancer-promoter looping. Consistently, the 338 GöPDX4 ΔNp63-bound SEs identified in (Fig. 1D) displayed similar ΔNp63-dependent active enhancer features in L3.6pl cells based on H3K27ac and nascent transcription profile changes (Fig. 3D-F; Supplementary Fig. S3A and S3B).

Figure 3.

B-STEP enhancers engage in long-range chromatin interactions to control target gene expression. A, Genes interacting with ΔNp63-dependent enhancers were identified from H3K27ac HiChIP and overlapped with genes that are downregulated in ΔNp63-depleted L3.6pl cells (log2fold change <−1; Padj < 0.05). B, The 601 ΔNp63-dependend enhancer regions identified in (Fig. 2D) were categorized on the basis of number of interactions identified through H3K27ac HiChIP into different groups as stated in the figure. Box plots shows average leChro-seq signal intensity [log2 values of counts per million (CPM)] within each interaction group. C, Shmechatic representation of transcriptional activity of enhancers correlating with interaction frequency originating from these enhancers. D–F, IGV tracks show H3K27ac HiChIP interactions as well as ΔNp63 and H3K27ac occupancy in reads per genomic content (RGC) and leChRO-seq signal in CPM around the ΔNp63-dependent enhancers of MIR205HG, FAT2, and SNAI2 genes in L3.6pl cells and in GöPDX4. G, Relative expression of ΔNp63, SNAI2, and FAT2 genes was assessed after CRISPR dCas9-KRAB mediated silencing of TPRG1, SNAI2, and FAT2 SEs using target specific guide RNAs (gRNA). Silencing of SNAI2 TSS was used as positive control. tracrRNA, transactivating crRNA without the target specificity. Mean ± SD, n = 3. One-way ANOVA followed by Dunnett test were performed for statistical analysis where ***, P < 0.001; **, P < 0.01; *, P < 0.05.

Figure 3.

B-STEP enhancers engage in long-range chromatin interactions to control target gene expression. A, Genes interacting with ΔNp63-dependent enhancers were identified from H3K27ac HiChIP and overlapped with genes that are downregulated in ΔNp63-depleted L3.6pl cells (log2fold change <−1; Padj < 0.05). B, The 601 ΔNp63-dependend enhancer regions identified in (Fig. 2D) were categorized on the basis of number of interactions identified through H3K27ac HiChIP into different groups as stated in the figure. Box plots shows average leChro-seq signal intensity [log2 values of counts per million (CPM)] within each interaction group. C, Shmechatic representation of transcriptional activity of enhancers correlating with interaction frequency originating from these enhancers. D–F, IGV tracks show H3K27ac HiChIP interactions as well as ΔNp63 and H3K27ac occupancy in reads per genomic content (RGC) and leChRO-seq signal in CPM around the ΔNp63-dependent enhancers of MIR205HG, FAT2, and SNAI2 genes in L3.6pl cells and in GöPDX4. G, Relative expression of ΔNp63, SNAI2, and FAT2 genes was assessed after CRISPR dCas9-KRAB mediated silencing of TPRG1, SNAI2, and FAT2 SEs using target specific guide RNAs (gRNA). Silencing of SNAI2 TSS was used as positive control. tracrRNA, transactivating crRNA without the target specificity. Mean ± SD, n = 3. One-way ANOVA followed by Dunnett test were performed for statistical analysis where ***, P < 0.001; **, P < 0.01; *, P < 0.05.

Close modal

One example of a “Signature 2” gene associated with a ΔNp63-dependent eRNA-producing SE is the SNAI2 gene (Fig. 3E). SNAI2 was not previously recognized as a direct target of ΔNp63-dependent SE due to the limitations of the 2D genome-association rules typically used (nearest genes), but was identified via our HiChIP analyses and thereby revealed a mechanistic link between ΔNp63 expression and SNAI2-driven epithelial-to-mesenchymal (EMT) transition (44–46). To assess the specific role that ΔNp63-dependent, eRNA-producing enhancers play in controlling target gene expression, we performed CRISPR/dCas9-KRAB–mediated epigenome editing to inactivate the SE associated with the FAT2, SNAI2 and TP63 genes (Supplementary Fig. S3C; Fig. 3G). In case of the latter, TPRG1 was shown to harbor a SE that directly regulates ΔNp63 expression (15). Notably, dCas9-KRAB-mediated inactivation of SE activity resulted in the significant and specific downregulation of all 3 target genes. Furthermore, consistent with the master regulatory function of ΔNp63, silencing of the TPRG1 SE upstream of the TP63 gene not only attenuated expression of ΔNp63, but also resulted in downregulation of the ΔNp63-dependent, SE-driven genes SNAI2 and FAT2. In contrast, inactivation of the SE associated with SNAI2 or FAT2 specifically downregulated the expression of the associated target gene, but not the other SE-associated genes investigated (Fig. 3G). Altogether, these results indicate that the ΔNp63-established B-STEP is a defining feature of the most aggressive subtype of pancreatic cancer.

Subtype-specific eRNA detection

The highly specific expression pattern of eRNAs suggests that they might serve as specific and sensitive markers of enhancer activity and therefore aid in the selection of therapy and monitoring of therapeutic efficacy. To test the potential clinical utility of detecting B-STEP eRNAs, we evaluated the detectability of TPRG1 and SNAI2 eRNAs in FFPE xenograft samples from the Basal-like PDAC cell line L3.6pl and a negative control colorectal cancer cell line (HT29) by using the RNAScope RNA-ISH approach (Fig. 4A and B). Remarkably, one to two individual foci could be observed for TPRG1 and SNAI2 in many nuclei of the L3.6pl xenograft sample, but not in the colorectal cancer xenograft. To confirm the feasibility of this approach on clinical material, we used histologic tumor material obtained from a Mayo Clinic patient following resection, which we confirmed to be of the Basal subtype based on the expression of both ΔNp63 (“p40”) and its downstream target KRT5 (Fig. 4C). Notably, positively stained nuclear foci could be observed for both TPRG1 and SNAI2 eRNAs in individual cells (Fig. 4D).

Figure 4.

PDAC subtypes can be defined on the basis of eRNA detection. A and B, RNAscope detection of RNA in (A) HT29 and (B) L3.6pl xenografts using probes against SNAI2 and TPRG1 eRNAs. PPIB (Peptidylprolyl isomerase B) probe was used as positive control, whereas probe against the bacterial gene dapB – as negative control. Nuclei were counterstained with hematoxylin. Scale bar = 20 μm for HT29 and 10 μm for L3.6pl xenograft images. C, IHC detection of ΔNp63 protein (p40) and KRT5 in primary tumor material from a patient with PDAC. Scale bar = 100 μm. D, RNAscope detection of RNA in primary FFPE tissue from PDAC patient from (C) using the probes against the TPRG1 and SNAI2 eRNAs. Detected RNA molecules in form of single dots are indicated with arrows. Scale bar = 20 μm. E, Gene set enrichment analysis was performed to compare the enrichment of “ΔNp63-dependent genes” in L3.6pl cells following 24 hours of 500 nmol/L JQ1 treatment. Normalized enrichment score (NES) and FDRs are shown within the plot. F and G, The relative expression analysis of Basal-like A marker (F) genes and (G) eRNAs was evaluated via qRT-PCR in L3.6pl cells following 24-hour 500 nmol/L JQ1 treatment. Mean ± SD, n = 3. Student t test were performed for statistical analysis where ***, P < 0.001; **, P < 0.01; *, P < 0.05. H, Cell proliferation analysis of Basal-like A L3.6pl cells (left) and classical AsPC1 cells (right) following treatment with different concentrations of pan-BET inhibitor CC-90010. Mean ± SD, n = 3. One-way ANOVA followed by Dunnett test were performed for statistical analysis where ***, P < 0.001; **, P < 0.01; *, P < 0.05. I, Relative tumor growth was monitored in mice transplanted with patient-derived Basal-like PDAC biopsies (Bo85PDX) over the time of 24 days following 50 mg/kg i.p. JQ1 treatment. P < 0.0001 (estimated using two-way ANOVA).

Figure 4.

PDAC subtypes can be defined on the basis of eRNA detection. A and B, RNAscope detection of RNA in (A) HT29 and (B) L3.6pl xenografts using probes against SNAI2 and TPRG1 eRNAs. PPIB (Peptidylprolyl isomerase B) probe was used as positive control, whereas probe against the bacterial gene dapB – as negative control. Nuclei were counterstained with hematoxylin. Scale bar = 20 μm for HT29 and 10 μm for L3.6pl xenograft images. C, IHC detection of ΔNp63 protein (p40) and KRT5 in primary tumor material from a patient with PDAC. Scale bar = 100 μm. D, RNAscope detection of RNA in primary FFPE tissue from PDAC patient from (C) using the probes against the TPRG1 and SNAI2 eRNAs. Detected RNA molecules in form of single dots are indicated with arrows. Scale bar = 20 μm. E, Gene set enrichment analysis was performed to compare the enrichment of “ΔNp63-dependent genes” in L3.6pl cells following 24 hours of 500 nmol/L JQ1 treatment. Normalized enrichment score (NES) and FDRs are shown within the plot. F and G, The relative expression analysis of Basal-like A marker (F) genes and (G) eRNAs was evaluated via qRT-PCR in L3.6pl cells following 24-hour 500 nmol/L JQ1 treatment. Mean ± SD, n = 3. Student t test were performed for statistical analysis where ***, P < 0.001; **, P < 0.01; *, P < 0.05. H, Cell proliferation analysis of Basal-like A L3.6pl cells (left) and classical AsPC1 cells (right) following treatment with different concentrations of pan-BET inhibitor CC-90010. Mean ± SD, n = 3. One-way ANOVA followed by Dunnett test were performed for statistical analysis where ***, P < 0.001; **, P < 0.01; *, P < 0.05. I, Relative tumor growth was monitored in mice transplanted with patient-derived Basal-like PDAC biopsies (Bo85PDX) over the time of 24 days following 50 mg/kg i.p. JQ1 treatment. P < 0.0001 (estimated using two-way ANOVA).

Close modal

Our previous work revealed an essential role of the Bromodomain and Extraterminal (BET) histone reader protein BRD4 in controlling eRNA expression in breast cancer (47). Moreover, BET family member bromodomain testis-specific protein (BRDT) was found to co-localize and interact with ΔNp63 to drive a unique transcriptional program and modulate cell phenotype in esophageal squamous cell carcinoma (37). Given the central role of BET proteins as a defining feature of SE activity (48, 49) in cancer, we hypothesized that small molecule inhibition of BET proteins may inhibit B-STEP-dependent transcription of the Basal-like A-specific gene expression pattern. Consistently, the BET inhibitor JQ1 downregulated ΔNp63-dependent (Fig. 4E) and Basal-like A-defining “Signature 2” genes (Fig. 4F). Moreover, similar to the effects of ΔNp63 loss (Supplementary Fig. S3D and S3E), JQ1 treatment abrogated expression of B-STEP eRNAs (Fig. 4G). Importantly, based on cell growth inhibition in L3.6pl cells ΔNp63-program conferred higher sensitivity to trotabresib, a reversible and orally active pan-BET inhibitor CC-90010 (50) currently under phase-I clinical trials (NCT03220347), when compared with cells of classical origin AsPC1 (ref. 8; Fig. 4H). Finally, in line with the anti-cancer effect of BETi, JQ1 treatment of mice implanted with a Basal-like A PDX (Supplementary Fig. S3F) substantially impaired tumor growth (Fig. 4I).

In summary, we conclude that B-STEP eRNAs can be detected on a single-cell level in histologic and primary patient material, are highly sensitive to BETi treatment and can therefore likely serve as predictive markers for therapy stratification.

A number of recent studies have examined the transcriptomic and genomic diversity within PDAC and uncovered distinct subtypes of PDAC with different clinical courses and therapeutic responsiveness (2–5, 51, 52). In this work, we have uncovered important epigenetic principles controlling the most aggressive Basal-like A subtype of PDAC. Growing evidence indicates that epigenetic reprogramming of the enhancer landscape is a key regulatory event controlling tumor cell identity. However, while technological advances have enabled the detailed interrogation of the mechanistic and epigenetic underpinnings of tumor subtype-specific gene expression in cell culture, the application of these findings to the clinical realm has been halted by the impracticality and difficulty of applying complex high throughput procedures to clinical material. In this study, we used a multi-omics approach to characterize the epigenetic and transcriptional landscape of clinically relevant model systems and identified a novel B-STEP we have termed B-STEP. Our findings demonstrate that eRNAs can function as a robust proxy of subtype-specific enhancer activity in PDAC.

To date, the application of genomic methodologies for the identification of subtype-specific enhancers in patients has been fraught with difficulties including technical limitations due to low quantities of primary biopsy material, extensive desmoplasia and high tumor heterogeneity. For the first time we not only identified tumor subtype-specific eRNAs, but also demonstrated the ability to detect these in histologic tumor material. We propose that detection of eRNAs in histologic material may serve as a novel diagnostic approach for patient stratification on a single-cell level. eRNA expression reflects mRNA expression kinetics more closely than H3K27ac occupancy at enhancers (53) and alterations in enhancer transcription precede the successive wave of changes in transcription (54), suggesting that analysis of eRNAs underlying complex cell identity programs in patient tumors and/or circulating tumor cells can be alternatively used as a tool for patient stratification and the assessment of therapeutic responses.

A major hurdle to identifying mechanism-based, clinically meaningful biomarkers is the availability of experimental model systems that faithfully recapitulate patient tumor characteristics. Established cancer cell lines often significantly diverge from primary tumor samples due to the high selective pressure of 2D cell culture conditions and lack microenviromental signals. Our results show that PDX models largely reflect the diversity in gene expression programs observed in clinical PDAC samples. Moreover, consistent with our previous work (9), our results demonstrate that the highly metastatic L3.6pl cell line closely mirrors the transcriptional and epigenetic enhancer program present in Basal-like A PDX models.

Functionally, eRNAs likely promote target gene activation and/or promoter/enhancer looping by interacting with transcriptional activators and co-activators (23, 27–29). Given their short length and unstable nature, discovery of eRNAs has proven challenging, especially in patient-derived clinical samples with substantial RNA degradation. However, this limitation was recently overcome by the development of leChRO-seq (36). Application of this technique in L3.6pl cells and a Basal-like A PDX enabled us to identify the B-STEP. Moreover, consistent with the recent data by Hamdan and colleagues (42) and Briend and colleagues (43), we similarly observed a positive correlation between contact frequency and the transcriptional activity of B-STEP enhancers, which also supports a role of eRNAs in promoting enhancer-promoter interactions.

Notably, while epigenetic activation of enhancers has frequently been correlated with cell state-specific gene expression, the direct interaction of enhancers with target genes and their necessity for the expression of specific genes has rarely been investigated. Our HiChIP analyses enabled us to identify direct targets of B-STEP enhancers. Moreover, we confirmed the essential role of these enhancers via epigenome editing. Intriguingly, one B-STEP ΔNp63-dependent target gene is the EMT TF SNAI2, which is one of the “Signature 2” genes (7). Our results reveal that SNAI2 expression in Basal-like A PDAC is under the direct control of an upstream ΔNp63-dependent SE. The aggressive and highly invasive nature of the Basal-like A PDAC subtype has been partially explained by its higher enrichment of EMT gene signatures (7, 55) and recently SNAI2 was shown to contribute to a more aggressive behavior of PDAC by inducing EMT in response to nutrient deprivation (46). Although the role of the ΔNp63-SNAI2 axis in other cancers has been reported before (44, 45), we describe the first mechanistic link explaining ΔNp63 and SNAI2 co-expression in Basal-like A PDAC. Moreover, we confirm that ΔNp63 plays a key master regulatory function in establishing the transcribed enhancer landscape of Basal-like A PDAC subtype. Furthermore, perturbation of the TPRG1 SE that controls ΔNp63 expression (15) is sufficient to block downstream ΔNp63-driven enhancer activity and target gene expression.

Another important finding from our work is the discovery that the B-STEP program is particularly amenable to modulation by small molecule BET inhibitor treatment both in vitro and in vivo. This finding is particularly relevant for the clinical application of BETi for PDAC treatment. A major limitation to the clinical utility of BETi is the occurrence of treatment-emergent adverse events (TEAE) such as thrombocytopenia. To date, BETi have only been marginally effective as a monotherapy in clinical trials for solid tumors, probably due to the limited tumor cell cytotoxicity that they invoke. On the basis of our findings we propose a new subtype-specific mechanism-based strategy to BETi application and dosage in PDAC. Specifically, we propose that B-STEP eRNA transcription may be a better indicator of on-target effects and can serve as a highly sensitive and accurate readout for clinical dosing and monitoring of BETi treatment in Basal-like A PDAC. Importantly, we believe this approach may enable lower, more accurate BETi dosing that will help to minimize TEAE. On the basis of the findings that Basal-like A PDAC is generally refractory to standard of care chemotherapy, we anticipate that transient inhibition of the B-STEP program and hence of the Basal-like A phenotype by BETi treatment can serve as an ideal inductive therapy to sensitize Basal-like A PDAC tumors to chemotherapy treatment. According to our model, verification of B-STEP eRNA expression may serve as an ideal biomarker for early phase clinical trial dosing studies and later for patient stratification (i.e., inclusion criteria for BETi clinical trials) as well as for monitoring of on-target therapeutic responsiveness. We envision that therapeutic strategies targeting cell identity programs will significantly benefit from the utilization of improved mechanism-based biomarkers and enable a better, more accurate identification and evaluation of therapy-induced effects.

Limitations of the study

To validate the potential utility of subtype-specific eRNA expression as a marker for patient/treatment stratification, more patient and PDX samples need to be included in the study. However, given the low level of eRNA expression and their high instability it is technically challenging to obtain suitable samples for analysis. Very often eRNA integrity is highly compromised in archived FFPE patient material with varying fixation conditions. Moreover, given the fact that most Basal A tumors are diagnosed at an advanced stage (i.e., using biopsies) fresh material is not readily available and not in sufficient quantity. While we agree that large numbers would strengthen our findings, we do believe that the data presented in the study demonstrate that eRNAs can be detected under suitable conditions, both in xenograft and primary patient material.

P. Hermann reports grants from Deutsche Krebshilfe (Max Eder Program), Hector Foundation; and grants from Deutsche Forschungsgemeinschaft (SFB1279) during the conduct of the study. J.T. Siveke reports grants and personal fees from Bristol-Myers Squibb/Celgene, Roche/Genentech, AstraZeneca,; grants from Eisbach Bio, Abalos Therapeutics; personal fees from Immunocore, Bayer, MSD Sharp Dohme, Novartis, Servier; and other support from Pharma15 outside the submitted work. No disclosures were reported by the other authors.

X. Wang: Data curation, software, formal analysis, investigation, visualization, writing–original draft, writing–review and editing. A.P. Kutschat: Formal analysis, investigation, visualization, methodology, writing–original draft, writing–review and editing. J. Aggrey-Fynn: Investigation, visualization, writing–review and editing. F.H. Hamdan: Investigation, writing–review and editing. R.P. Graham: Methodology, writing–review and editing. A.Q. Wixom: Investigation, writing–review and editing. Y. Souto: Validation, investigation, writing–review and editing. S. Ladigan-Badura: Resources, methodology, writing–review and editing. J.A. Yonkus: Methodology, writing–review and editing. A.M. Abdelrahman: Resources, methodology, writing–review and editing. R. Alva-Ruiz: Resources, methodology, writing–review and editing. J. Gaedcke: Resources, methodology, writing–review and editing. P. Ströbel: Resources, methodology, writing–review and editing. R.L. Kosinsky: Resources, methodology, writing–review and editing. F. Wegwitz: Resources, methodology, writing–review and editing. P. Hermann: Funding acquisition, writing–review and editing. M.J. Truty: Resources, methodology, writing–review and editing. J.T. Siveke: Funding acquisition, writing–review and editing. S.A. Hahn: Resources, funding acquisition, methodology, writing–review and editing. E. Hessmann: Resources, funding acquisition, methodology, writing–review and editing. S.A. Johnsen: Conceptualization, formal analysis, supervision, funding acquisition, investigation, writing–original draft, project administration, writing–review and editing. Z. Najafova: Conceptualization, supervision, investigation, visualization, methodology, writing–original draft, project administration, writing–review and editing.

This project was funded by the Robert Bosch Stiftung (to S.A. Johnsen, Z. Najafova, R.L. Kosinsky), the Deutsche Krebshilfe (70112505; PIPAC Consortium; to S.A. Johnsen, E. Hessmann, S.A. Hahn, J.T. Siveke, P. Hermann, F. Wegwitz, J. Gaedcke, S. Ladigan-Badura, and Z. Najafova), a grant from NCI (CA 102701; to S.A. Johnsen), a grant from National Institute of Diabetes and Digestive and Kidney Diseases (T32 DK07198; to A.Q. Wixom), grants from the Deutsche Forschungsgemeinschaft (DFG, JO 815/3–2 and KFO 5002, to S.A. Johnsen and E. Hessmann, P. Strobel, J. Gaedcke, respectively). Z. Najafova was partially supported by the Dorothea Schlözer Program (University of Göttingen). This work was supported by the Optical Microscopy Core of the Mayo Clinic Center for Cell Signaling in Gastroenterology (P30DK084567). Work in the lab of J.T. Siveke is supported by the German Cancer Consortium (DKTK), Deutsche Forschungsgemeinschaft (DFG; 405344257/SI1549/3–2 and SI1549/4–1); the German Federal Ministry of Education and Research (BMBF; 01KD2206A/SATURN3). The authors thank the Genome Analysis Core at the Mayo Clinic for performing the sequencing and the Pathology Research Core at Mayo Clinic for performing IHC analysis and RNAScope.

The publication costs of this article were defrayed in part by the payment of publication fees. Therefore, and solely to indicate this fact, this article is hereby marked “advertisement” in accordance with 18 USC section 1734.

Note: Supplementary data for this article are available at Molecular Cancer Research Online (http://mcr.aacrjournals.org/).

1.
Rahib
L
,
Smith
BD
,
Aizenberg
R
,
Rosenzweig
AB
,
Fleshman
JM
,
Matrisian
LM
.
Projecting cancer incidence and deaths to 2030: the unexpected burden of thyroid, liver, and pancreas cancers in the United States
.
Cancer Res
2014
;
74
:
2913
21
.
2.
Collisson
EA
,
Sadanandam
A
,
Olson
P
,
Gibb
WJ
,
Truitt
M
,
Gu
S
, et al
.
Subtypes of pancreatic ductal adenocarcinoma and their differing responses to therapy
.
Nat Med
2011
;
17
:
500
3
.
3.
Bailey
P
,
Chang
DK
,
Nones
K
,
Johns
AL
,
Patch
AM
,
Gingras
MC
, et al
.
Genomic analyses identify molecular subtypes of pancreatic cancer
.
Nature
2016
;
531
:
47
52
.
4.
Moffitt
RA
,
Marayati
R
,
Flate
EL
,
Volmar
KE
,
Loeza
SGH
,
Hoadley
KA
, et al
.
Virtual microdissection identifies distinct tumor- and stroma-specific subtypes of pancreatic ductal adenocarcinoma
.
Nat Genet
2015
;
47
:
1168
78
.
5.
Puleo
F
,
Nicolle
R
,
Blum
Y
,
Cros
J
,
Marisa
L
,
Demetter
P
, et al
.
Stratification of pancreatic ductal adenocarcinomas based on tumor and microenvironment features
.
Gastroenterology
2018
;
155
:
1999
2013
.
6.
Cancer Genome Atlas Research Network
.
Integrated genomic characterization of pancreatic ductal adenocarcinoma
.
Cancer Cell
2017
;
32
:
185
203
.
7.
Chan-Seng-Yue
M
,
Kim
JC
,
Wilson
GW
,
Ng
K
,
Figueroa
EF
,
O'Kane
GM
, et al
.
Transcription phenotypes of pancreatic cancer are driven by genomic events during tumor evolution
.
Nat Genet
2020
;
52
:
231
40
.
8.
Somerville
TDD
,
Xu
Y
,
Miyabayashi
K
,
Tiriac
H
,
Cleary
CR
,
Maia-Silva
D
, et al
.
TP63-mediated enhancer reprogramming drives the squamous subtype of pancreatic ductal adenocarcinoma
.
Cell Rep
2018
;
25
:
1741
55
.
9.
Hamdan
FH
,
Johnsen
SA
.
DeltaNp63-dependent super-enhancers define molecular identity in pancreatic cancer by an interconnected transcription factor network
.
Proc Natl Acad Sci USA
2018
;
115
:
E12343
52
.
10.
Barbieri
CE
,
Pietenpol
JA
.
p63 and epithelial biology
.
Exp Cell Res
2006
;
312
:
695
706
.
11.
Westfall
MD
,
Pietenpol
JA
.
p63: molecular complexity in development and cancer
.
Carcinogenesis
2004
;
25
:
857
64
.
12.
DeYoung
MP
,
Ellisen
LW
.
p63 and p73 in human cancer: defining the network
.
Oncogene
2007
;
26
:
5169
83
.
13.
Basturk
O
,
Khanani
F
,
Sarkar
F
,
Levi
E
,
Cheng
JD
,
Adsay
NV
.
DeltaNp63 expression in pancreas and pancreatic neoplasia
.
Mod Pathol
2005
;
18
:
1193
8
.
14.
Ito
Y
,
Takeda
T
,
Wakasa
K
,
Tsujimoto
M
,
Sakon
M
,
Matsuura
N
.
Expression of p73 and p63 proteins in pancreatic adenocarcinoma: p73 overexpression is inversely correlated with biological aggressiveness
.
Int J Mol Med
2001
;
8
:
67
71
.
15.
Andricovich
J
,
Perkail
S
,
Kai
Y
,
Casasanta
N
,
Peng
W
,
Tzatsos
A
.
Loss of KDM6A activates super-enhancers to induce gender-specific squamous-like pancreatic cancer and confers sensitivity to BET inhibitors
.
Cancer Cell
2018
;
33
:
512
26
.
16.
Somerville
TD
,
Biffi
G
,
Daßler-Plenker
J
,
Hur
SK
,
He
XY
,
Vance
KE
, et al
.
Squamous trans-differentiation of pancreatic cancer cells promotes stromal inflammation
.
Sawyers
CL
,
Ojala
PM
,
Oliver
T
,
editors
.
eLife
.
eLife Sciences Publications, Ltd;
2020
;
9
:
e53381
.
17.
Lewis
MW
,
Li
S
,
Franco
HL
.
Transcriptional control by enhancers and enhancer RNAs
.
Transcription
2019
;
10
:
171
86
.
18.
Chen
H
,
Li
C
,
Peng
X
,
Zhou
Z
,
Weinstein
JN
,
Cancer Genome Atlas Research Network
;
Liang
H
.
A pan-cancer analysis of enhancer expression in nearly 9000 patient samples
.
Cell
2018
;
173
:
386
99
.
19.
Sur
I
,
Taipale
J
.
The role of enhancers in cancer
.
Nat Rev Cancer
2016
;
16
:
483
93
.
20.
De Santa
F
,
Barozzi
I
,
Mietton
F
,
Ghisletti
S
,
Polletti
S
,
Tusi
BK
, et al
.
A large fraction of extragenic RNA pol II transcription sites overlap enhancers
.
PLoS Biol
2010
;
8
:
e1000384
.
21.
Kim
TK
,
Hemberg
M
,
Gray
JM
,
Costa
AM
,
Bear
DM
,
Wu
J
, et al
.
Widespread transcription at neuronal activity-regulated enhancers
.
Nature
2010
;
465
:
182
7
.
22.
Zhao
Y
,
Wang
L
,
Ren
S
,
Wang
L
,
Blackburn
PR
,
McNulty
MS
, et al
.
Activation of P-TEFb by androgen receptor-regulated enhancer RNAs in castration-resistant prostate cancer
.
Cell Rep
2016
;
15
:
599
610
.
23.
Schaukowitch
K
,
Joo
JY
,
Liu
X
,
Watts
JK
,
Martinez
C
,
Kim
TK
.
Enhancer RNA facilitates NELF release from immediate early genes
.
Mol Cell
2014
;
56
:
29
42
.
24.
Lai
F
,
Orom
UA
,
Cesaroni
M
,
Beringer
M
,
Taatjes
DJ
,
Blobel
GA
, et al
.
Activating RNAs associate with mediator to enhance chromatin architecture and transcription
.
Nature
2013
;
494
:
497
501
.
25.
Rahnamoun
H
,
Lee
J
,
Sun
Z
,
Lu
H
,
Ramsey
KM
,
Komives
EA
, et al
.
RNAs interact with BRD4 to promote enhanced chromatin engagement and transcription activation
.
Nat Struct Mol Biol
2018
;
25
:
687
97
.
26.
Bose
DA
,
Donahue
G
,
Reinberg
D
,
Shiekhattar
R
,
Bonasio
R
,
Berger
SL
.
RNA binding to CBP stimulates histone acetylation and transcription
.
Cell
2017
;
168
:
135
49
.
27.
Hsieh
CL
,
Fei
T
,
Chen
Y
,
Li
T
,
Gao
Y
,
Wang
X
, et al
.
Enhancer RNAs participate in androgen receptor-driven looping that selectively enhances gene activation
.
Proc Natl Acad Sci USA
2014
;
111
:
7319
24
.
28.
Li
W
,
Notani
D
,
Ma
Q
,
Tanasa
B
,
Nunez
E
,
Chen
AY
, et al
.
Functional roles of enhancer RNAs for estrogen-dependent transcriptional activation
.
Nature
2013
;
498
:
516
20
.
29.
Sigova
AA
,
Abraham
BJ
,
Ji
X
,
Molinie
B
,
Hannett
NM
,
Guo
YE
, et al
.
Transcription factor trapping by RNA in gene regulatory elements
.
Science
2015
;
350
:
978
81
.
30.
Melo
CA
,
Drost
J
,
Wijchers
PJ
,
van de Werken
H
,
de Wit
E
,
Vrielink
JAFO
, et al
.
eRNAs are required for p53-dependent enhancer activity and gene transcription
.
Mol Cell
2013
;
49
:
524
35
.
31.
Franco
HL
,
Nagari
A
,
Malladi
VS
,
Li
W
,
Xi
Y
,
Richardson
D
, et al
.
Enhancer transcription reveals subtype-specific gene expression programs controlling breast cancer pathogenesis
.
Genome Res
2018
;
28
:
159
70
.
32.
Kutschat
AP
,
Hamdan
FH
,
Wang
X
,
Wixom
AQ
,
Najafova
Z
,
Gibhardt
CS
, et al
.
STIM1 mediates calcium-dependent epigenetic reprogramming in pancreatic cancer
.
Cancer Res
2021
;
81
:
2943
55
.
33.
Mishra
VK
,
Wegwitz
F
,
Kosinsky
RL
,
Sen
M
,
Baumgartner
R
,
Wulff
T
, et al
.
Histone deacetylase class-I inhibition promotes epithelial gene expression in pancreatic cancer cells in a BRD4- and MYC-dependent manner
.
Nucleic Acids Res
2017
;
45
:
6334
49
.
34.
Bruns
CJ
,
Harbison
MT
,
Kuniyasu
H
,
Eue
I
,
Fidler
IJ
.
In vivo selection and characterization of metastatic variants from human pancreatic adenocarcinoma by using orthotopic implantation in nude mice
.
Neoplasia
1999
;
1
:
50
62
.
35.
Napp
J
,
Pardo
LA
,
Hartung
F
,
Tietze
LF
,
Stühmer
W
,
Alves
F
.
In vivo imaging of tumor xenografts with an antibody targeting the potassium channel Kv10.1
.
Eur Biophys J
2016
;
45
:
721
33
.
36.
Chu
T
,
Rice
EJ
,
Booth
GT
,
Salamanca
HH
,
Wang
Z
,
Core
LJ
, et al
.
Chromatin run-on and sequencing maps the transcriptional regulatory landscape of glioblastoma multiforme
.
Nat Genet
2018
;
50
:
1553
64
.
37.
Wang
X
,
Kutschat
AP
,
Yamada
M
,
Prokakis
E
,
Böttcher
P
,
Tanaka
K
, et al
.
Bromodomain protein BRDT directs ΔNp63 function and super-enhancer activity in a subset of esophageal squamous cell carcinomas
.
Cell Death Differ
2021
;
28
:
2207
20
.
38.
Lomberk
G
,
Blum
Y
,
Nicolle
R
,
Nair
A
,
Gaonkar
KS
,
Marisa
L
, et al
.
Distinct epigenetic landscapes underlie the pathobiology of pancreatic cancer subtypes
.
Nat Commun
2018
;
9
:
1978
.
39.
Hnisz
D
,
Abraham
BJ
,
Lee
TI
,
Lau
A
,
Saint-André
V
,
Sigova
AA
, et al
.
Super-enhancers in the control of cell identity and disease
.
Cell
2013
;
155
:
934
47
.
40.
He
Y
,
Long
W
,
Liu
Q
.
Targeting super-enhancers as a therapeutic strategy for cancer treatment
.
Front Pharmacol
2019
;
10
:
361
.
41.
Li
W
,
Notani
D
,
Rosenfeld
MG
.
Enhancers as noncoding RNA transcription units: recent insights and future perspectives
.
Nat Rev Genet
2016
;
17
:
207
23
.
42.
Hamdan
FH
,
Abdelrahman
AM
,
Kutschat
AP
,
Wang
X
,
Ekstrom
TL
,
Jalan-Sakrikar
N
, et al
.
Interactive enhancer hubs (iHUBs) mediate transcriptional reprogramming and adaptive resistance in pancreatic cancer
.
Gut
2023
;
72
:
1174
85
.
43.
Briend
M
,
Rufiange
A
,
Moncla
L-HM
,
Mathieu
S
,
Bossé
Y
,
Mathieu
P
.
Connectome and regulatory hubs of CAGE highly active enhancers
.
Sci Rep
2023
;
13
:
5594
.
44.
Srivastava
K
,
Pickard
A
,
Craig
SG
,
Quinn
GP
,
Lambe
SM
,
James
JA
, et al
.
ΔNp63γ/SRC/Slug signaling axis promotes epithelial-to-mesenchymal transition in squamous cancers
.
Clin Cancer Res
2018
;
24
:
3917
27
.
45.
Dang
TT
,
Westcott
JM
,
Maine
EA
,
Kanchwala
M
,
Xing
C
,
Pearson
GW
.
ΔNp63α induces the expression of FAT2 and Slug to promote tumor invasion
.
Oncotarget
2016
;
7
:
28592
611
.
46.
Recouvreux
MV
,
Moldenhauer
MR
,
Galenkamp
KMO
,
Jung
M
,
James
B
,
Zhang
Y
, et al
.
Glutamine depletion regulates Slug to promote EMT and metastasis in pancreatic cancer
.
J Exp Med
2020
;
217
:
e20200388
.
47.
Nagarajan
S
,
Hossan
T
,
Alawi
M
,
Najafova
Z
,
Indenbirken
D
,
Bedi
U
, et al
.
Bromodomain Protein BRD4 is required for estrogen receptor-dependent enhancer activation and gene transcription
.
Cell Rep
2014
;
8
:
460
9
.
48.
Lovén
J
,
Hoke
HA
,
Lin
CY
,
Lau
A
,
Orlando
DA
,
Vakoc
CR
, et al
.
Selective inhibition of tumor oncogenes by disruption of super-enhancers
.
Cell
2013
;
153
:
320
34
.
49.
Shi
J
,
Vakoc
CR
.
The mechanisms behind the therapeutic activity of BET bromodomain inhibition
.
Mol Cell
2014
;
54
:
728
36
.
50.
Moreno
V
,
Sepulveda
JM
,
Vieito
M
,
Hernández-Guerrero
T
,
Doger
B
,
Saavedra
O
, et al
.
Phase I study of CC-90010, a reversible, oral BET inhibitor in patients with advanced solid tumors and relapsed/refractory non-Hodgkin's lymphoma
.
Ann Oncol
2020
;
31
:
780
8
.
51.
Raphael
BJ
,
Hruban
RH
,
Aguirre
AJ
,
Moffitt
RA
,
Yeh
JJ
,
Stewart
C
, et al
.
Integrated genomic characterization of pancreatic ductal adenocarcinoma
.
Cancer Cell
2017
;
32
:
185
203
.
52.
Aung
KL
,
Fischer
SE
,
Denroche
RE
,
Jang
GH
,
Dodd
A
,
Creighton
S
, et al
.
Genomics-driven precision medicine for advanced pancreatic cancer: early results from the COMPASS trial
.
Clin Cancer Res
2018
;
24
:
1344
54
.
53.
Tyssowski
KM
,
DeStefino
NR
,
Cho
JH
,
Dunn
CJ
,
Poston
RG
,
Carty
CE
, et al
.
Different neuronal activity patterns induce different gene expression programs
.
Neuron
2018
;
98
:
530
46
.
54.
Arner
E
,
Daub
CO
,
Vitting-Seerup
K
,
Andersson
R
,
Lilje
B
,
Drabløs
F
, et al
.
Transcribed enhancers lead waves of coordinated transcription in transitioning mammalian cells
.
Science
2015
;
347
:
1010
4
.
55.
Mueller
S
,
Engleitner
T
,
Maresch
R
,
Zukowska
M
,
Lange
S
,
Kaltenbacher
T
, et al
.
Evolutionary routes and KRAS dosage define pancreatic cancer phenotypes
.
Nature
2018
;
554
:
62
8
.