Cell culture
Cells were cultured in a humidified 37°C incubator with 5% CO2. K562 and Jurkat cells were grown in RPMI-1640 medium supplemented with 10% fetal bovine serum and 1X penicillin streptomycin antibiotic. For the PARP1 inhibition experiment, cells were treated with either 10 µM Olaparib, initially diluted in DMSO, or DMSO only, for 90 min before immediately crosslinked for Micro-C or lysed in NUN buffer (0.3 M NaCl, 1 M Urea, 1% NP-40, 20 mM HEPES, pH 7.5, 7.5 mM MgCl2, 0.2 mM EDTA, 1 mM DTT, 20 units per ml SUPERase In RNase Inhibitor (Life Technologies, AM2694), 1X Protease Inhibitor Cocktail (Roche, 11 873 580 001)) for ChRO-seq.
mECSs harboring a homozygous endogenous NELFB-FKBP12F36V fusion protein were cultured on 0.1% gelatin (Millipore) in PBS+/+ coated tissue-culture grade plates. For routine culture, cells were grown in Serum/LIF conditions: DMEM (Gibco), supplemented with 2 mM L-glutamine (Gibco), 1x MEM non-essential amino acids (Gibco), 1 mM sodium pyruvate (Gibco), 100 U/ml penicillin and 100 U/ml streptomycin (Gibco), 0.1 mM 2-mercaptoethanol (Gibco), 15% Fetal Bovine Serum (Gibco), and 1000 U/ml of recombinant leukemia inhibitory factor (LIF).
To induce NELFB degradation, dTAG-13 (Bio-Techne) was reconstituted in DMSO (Sigma) at 5 mM. dTAG-13 was diluted in maintenance medium to 500 nM and added to cells with medium changes for the specified amounts of time. For dTAG washes, the cells were washed 4 times, twice with PBS +/+ and twice with maintenance medium following the treatment time to ensure complete removal of the dTAG ligand. At the end of each dTAG-13 treatment time point, cells were detached using Trypsin-EDTA (0.05%) (Gibco) and counted before crosslinking for Micro-C.
Micro-C
Micro-C for K562, Jurkat and mESCs was performed by following the published protocol for mammalian Micro-C21,22,83. Cells were crosslinked with 1 ml per million cells of 1% formaldehyde for 10 minutes at room temperature and quenched by 0.25 M Glycine for 5 min. After spin-down for 5 minutes at 300Xg at 4°C, cells were washed at a density of 1 ml per million cells in ice cold PBS. Cells were crosslinked a second time, with 1 ml per 4 million cells of 3 mM disuccinimidyl glutarate (DSG) (ThermoFisher Scientific, 20593) for 40 min at room temperature and quenched by 0.4 M Glycine for 5 min. Following two washes with ice cold PBS, cells were flash-frozen and kept at -80°C until further use. For MNase digestion, cells were thawed on ice for 5 min, incubated with 1ml MB#1 buffer (10 mM Tris-HCl, pH 7.5, 50 mM NaCl, 5 mM MgCl2, 1 mM CaCl2, 0.2% NP-40, 1x Roche cOmplete EDTA-free (Roche diagnostics, 04693132001)) and washed twice with MB#1 buffer. MNase concentration for each cell type was predetermined using MNase titration experiments exploring 2.5-20U of MNase per million cells. We selected the MNase concentration that gives ~ 90% mononucleosomes. Chromatin was digested with MNase for 10 min at 37°C and digestion was stopped by adding 8 ul of 500 mM EGTA and incubating at 65°C for 10 min.
Following dephosphorylation with rSAP (NEB #M0371) and end polishing using T4 PNK (NEB #M0201), DNA polymerase Klenow fragment (NEB #M0210) and biotinylated dATP and dCTP (Jena Bioscience #NU-835-BIO14-S and #NU-809-BIOX-S, respectively), ligation was performed in a final volume of 2.5 ml for 3h at room temperature using T4 DNA ligase (NEB #M0202). Dangling ends were removed by a 5 min incubation with Exonuclease III (NEB #0206) at 37°C and biotin enrichment was done using 20 ul Dynabeads™ MyOne™ Streptavidin C1 beads (Invitrogen #65001). Libraries were prepared with the NEBNext Ultra II Library Preparation Kit (NEB #E7103). Samples were sequenced on a combination of Illumina’s NovaSeq 6000 and HiSeq 2500 at Novogene.
Micro-C data mapping and visualization
All Micro-C mapping was done using the mirnylab/distiller-nf: v0.3.3 pipeline84. Raw data were mapped to the hg38 human genome assembly (K562 and Jurkat) or mm10 mouse genome assembly (mESCs). For analysis of contacts in the MYC locus, data was mapped to hg19 human genome assembly due to a large gap present in this locus when mapping K562 sequencing data to hg38. For data visualization by contact maps, multi cool (mcool) files, balanced by iterative correction and eigenvector decomposition (ICE) for resolutions of 200 bp to 10 Mb were generated from contacts with both ends having a mapq score ≥ 30. Micro-C data visualization as contact maps in genome-browser shots with available PRO-seq, dREG, CRISPRi and histone marks tracks was done using the HiCExplorer tool85 and pyGenomeTracks86. Virtual 4C tracks were prepared as described previously49. 1D signal near enhancer and promoter TSSs (Fig. S2) was calculated based on the distiller-nf output pairs files, filtered for intra-chromosomal with mapq ≥ 30. Contacts assigned to the 5’ of single reads were shifted 75bp downstream, based on their orientation, to the probable center of the nucleosome.
ChRO-seq
For chromatin isolation from cultured K562 cells, we added 1 ml of 1X NUN buffer and vigorously vortexed the samples for 1 min. An additional 500 µl of NUN buffer was added and samples were vigorously vortexed for an additional 30 s. The samples were then incubated on ice for 30 min with a brief vortex every 10 min and centrifuged at 12,500g for 30 min at 4°C, after which the NUN buffer was removed from the chromatin pellet. The chromatin pellet was washed 3 times with 1 ml of 50 mM Tris-HCl, pH 7.5, supplemented with 40 U/ml of RNase inhibitor, centrifuged at 10,000Xg for 5 min at 4°C, and the supernatant discarded. After washing, 100 µl of chromatin storage buffer (50 mM Tris-HCl, pH 8.0, 25% glycerol, 5 mM Mg(CH3COO)2, 0.1 mM EDTA, 5 mM DTT, 40 U/ml Rnase inhibitor) was added to each sample. The samples were loaded into a Bioruptor and sonicated with the power setting on high, with a cycle time of 10 min with cycle durations of 30 s on and 30 s off. The sonication was repeated up to three times, as needed, to get the chromatin pellet into suspension. Samples were stored at − 80°C.
ChRO-seq library preparation was performed following a published protocol24, with minor modifications. Chromatin from 1 million K562 cells in 25 µl chromatin storage buffer was mixed with 25 µl of 2X chromatin run-on buffer (10 mM Tris-HCl, pH 8.0, 5 mM MgCl2,1 mM DTT, 300 mM KCl, 10 µM Biotin-11-ATP (Perkin Elmer, NEL544001EA), 40 µM Biotin-11-CTP (Perkin Elmer, NEL542001EA), 10 µM Biotin-11-GTP (Perkin Elmer, NEL545001EA), 40 µM Biotin-11-UTP (Perkin Elmer, NEL543001EA), 2ng/µl Yeast tRNA (VWR, 80054-306), 0.8 U/µl RNase inhibitor, 1% Sarkosyl (Fisher Scientific, AC612075000)). The run-on reaction was incubated in a thermomixer at 37°C for 5 min at 750 RPM. The reaction was stopped by adding 300 µl of Trizol LS (Life Technologies, 10296-010). To clean up the reaction prior to base hydrolysis, 40 µl of 1-bromo-3-chloropropane (BCP) (Sigma, B9673) was added to the samples and samples were vortexed for 20 sec, incubated for 3 min at RT and centrifuged at 17,000Xg at 4°C for 5 min. ~250 µl of aqueous phase was transferred to a new 1.5 ml tube and samples were pelleted at 17,000Xg at 4°C for 15 min in 650 ml ice-cold 100% EtOH with GlycoBlue (2.5 µl) (Ambion, AM9515) to visualize the RNA pellet. The pellet was further washed with 0.5 ml ice-cold 75% by vortex-mixing and centrifugation at 17,000Xg at 4°C for 5 min. Any residual EtOH was removed by air-drying the pellet for 5 min in RT. The RNA pellet was resuspended in 20 µl of diethylpyrocarbonate (DEPC)-treated water and heat denatured at 65°C for 40 sec. For base hydrolysis, 5 µl of 1N NaOH was added to the RNA sample to get 0.2N NaOH, and samples were incubated on ice for 8 min. Base hydrolysis was stopped by adding 25 µl of Tris-HCl pH 6.8 and gently mixing. 3′ Adapter ligation was done using T4 RNA Ligase 1 (NEB, M0204L). A first binding to streptavidin beads (NEB, S1421S) and washed as described24. RNA was removed from beads by Trizol (Life Technologies, 15596-026) and followed by a 5′ decapping with RNA 5′ pyrophosphohydrolase (RppH, NEB, M0356S). The 5′ end was phosphorylated using T4 polynucleotide kinase (NEB, M0201L). A second bead binding was performed, followed by an on-beads 5′ adapter ligation. After bead washes, RNA was removed from beads by adding 300 µl of TRI reagent (MRC, TR 118). 40 µl of BCP were added to the samples and samples were vortexed for 20 sec, incubated for 3 min at room temperature and centrifuged at 17,000Xg at 4°C for 5 min. ~180µl of aqueous phase was transferred to a new 1.5 ml tube and samples were pelleted at 17,000Xg at 4°C for 15 min in 450 ml ice-cold 100% EtOH with GlycoBlue (2.5 µl) to visualize the RNA pellet. The pellet was further washed with 0.4 ml ice-cold 75% EtOH by vortex-mixing and centrifugation at 17,000Xg at 4°C for 5 min. Any residual EtOH was removed by air-drying the pellet for 5 min in RT. A reverse transcription reaction using Superscript III Reverse Transcriptase (Life Technologies, 18080-044) was used to generate cDNA. cDNA was then amplified using Q5 High-Fidelity DNA Polymerase (NEB, M0491L) to generate the ChRO-seq libraries following the protocol provided by NEB. Libraries were sequenced using Illumina’s HiSeq 4000 at Novogene.
ChRO-seq, PRO-seq and GRO-seq data processing and analysis
Processing for ChRO-seq data (as well as PRO-seq and GRO-seq available raw data) in this study was done using the Proseq2.0 pipeline available from GitHub (https://github.com/Danko-Lab/proseq2.0)87. Differential expression analyses between K562 and Jurkat cells for pausing signal and gene body transcription levels was performed by DEseq288 either on signal between the TSS and 250bp downstream (pause signal) or signal downstream to the first 250 bp through the annotated (GENCODE V29) polyadenylation cleavage site (gene body signal). For visualization of the changes, fold change in expression following NELFB-dTAG in mESCs or PARP1 inhibition in K562 was performed by deepTools bigwigCompare command at 1bp resolution, using 0.25 as pseudocount. For the NELFB-dTAG PRO-seq visualization, fold-change and normalized PRO-seq signal matrices were calculated in a stranded manner, followed by a contacennation of the two strands’ matrices to generate a single, stranded matrix (Fig. 4B).
Definition of TIRs, enhancers, promoters and transcription start sites
For mESCs we first defined TIRs genome-wide as detected by dREG30,31 using available GRO-seq data from mESCs47,89. To finely and unbiasedly define the position of transcription initiation at each of these TIRs, we used the position with the most 5’ mapped GRO-seq reads within the dREG peak (maxTSN). For the analyses of K562 cells, we first called TIRs using dREG from available PRO-seq data31. The center of these TIRs was defined as the center of enhancers and promoters for the analysis comparison of contacts between functional and nonfunctional enhancer-promoter pairs, based on CRISPRi data. For any further analyses the center of enhancers and promoters was defined as the maxTSNs, called using the data from coPRO with enrichment for 5’ capping (coPRO-capped)36. For the comparison between K562 and Jurkat cell lines, we called TIRs in both K562 and Jurkat using PRO-seq data31 and dREG and determined maxTSN based on coPRO-capped from K562 cells. We defined promoters based on the existence of any known human (K562 and Jurkat) of mouse (mESCs) stable 5’ mapped transcripts from CAGE90 within 5kb away in the direction of maximum initiation. In analyses including Jurkat and K562 cell lines, we considered only shared promoters based on proximity to the best transcription start site defined by the deconvolution of expression for nascent RNA-sequencing data (DENR), based on GENCODE V29 annotations, in both cell types. We used a combined set of enhancers from TIRs detected in both cell lines to define enhancers. Since promoters make a relatively small fraction of all TIRs found in the data and can act as enhancers for other distal genes91 we included promoters under the definition for enhancers whenever we calculated enhancer-promoter contacts genome-wide.
Micro-C contact normalization and comparison between CRISPRi-defined functional and nonfunctional pairs
Comparison between functional and nonfunctional enhancer-promoter pairs was based either on CRISPRi genetic screens for enhancer function either in the MYC locus, based on cell viability28 or based on expression from single-cell RNA sequencing analysis27. In addition to the seven active enhancers in the K562 MYC locus, we also detected 54 other TIRs that were marked by DNase-I hypersensitivity sites (DHSs) and histone modifications, were located in the same TAD, and were tested by CRISPRi, but which did not affect the growth rate of K562 cells. These TIRs were considered nonfunctional and compared to the seven functional MYC enhancers. For the genome wide-analysis based on data from27, we defined functional pairs as pairs where CRISPR inhibition of the enhancer resulted in a loss of a minimum of 10% of gene expression, with an empirical p-value < 0.05. nonfunctional pairs were defined as pairs with less of 5% effect on gene expression and empirical p-value > 0.9. To remove possible confounding effect, we filtered the functional and nonfunctional pairs to have similar distributions of enhancer-promoter contacts (limited between 5kb and 1Mb), accessibility (by ATAC-seq) and gene body transcription levels (by PRO-seq) in the target gene (Fig. S1C). All CRISPRi-targeted enhancers and target promoters were reassigned to their nearest dREG peak center, within 5kb, on the same strand. Enhancer-promoter contacts were defined as contacts that map to a 5kb window near the promoter on one end and the enhancer on the other hand. Contacts were normalized to the expected based on a non-parametric LOWESS smoothing of the contacts-by-distance function in a region corresponding to a 1Mb in the orientation of the promoter, relative to the enhancer (Fig. S1A). Observed over expected ratios were then compared between functional and nonfunctional pairs (Fig. 1B,C, Fig. S1B). Differences in contacts between CRISPRi-defined functional and nonfunctional pairs were calculated based on cell-by-cell differences between APA matrices for all functional and all nonfunctional pairs, normalized for the number of pairs. The differences were calculated as the medians of the differences based on 1000 bootstrapping iterations of the functional and nonfunctional pairs, to remove outlier background. These differences were presented as the number of contact differences per 1000 pairs (Fig. 1E and F). The APA matrices were centered on the coPRO-based maxTSN as the TSS assigned for each TIR.
Aggregated Peak Analysis (APA) and APA matrices normalization for comparison between samples
We expected significant changes in chromatin after manipulating Pol II transcription70,71. As such, not only are enhancer-promoter contacts expected to change, but the background contacts with at least one end originating at enhancer- and promoter- regions may be affected between conditions. As APAs are often used to characterize contacts7,22,46, we devised an APA that normalizes enhancer-promoter contacts to all contacts associated with enhancers and promoters. We calculated changes in aggregated 1D signal mapped to the same windows around enhancers and promoters, as used in the APA. We calculated the expected change in each pixel based on the ratios between the sample-specific sum of the 1D signal in each treatment condition (Fig. 2A). For all other, intra-sample APA analysis we used the aggregated raw counts. For APAs calculated at windows of 20kb around the anchors, we considered all possible anchor pairs within a genomic distance of 25-150kb. For higher resolution APAs with 2kb window around enhancer and promoter TSSs (Fig. 1D and Fig. S3), we considered all possible enhancer-promoter pairs within a genomic distance of 5-100kb.
Individual contact comparison between samples and treatments
We also devised an alternative normalization scheme which compares the number of contacts between enhancer-promoter pairs to the local background near each enhancer and promoter anchor. We calculated the number of contacts between each pair of anchors (enhancer-promoter or CTCF binding sites) using a 5kb window around each anchor. As a background, we counted the number of contacts between each anchor (in a 5kb window) and regions 10–150 kb from the second anchor (Fig. 2C). The ratio between the anchor-to-anchor contacts and background contacts was presented in scatterplots or used to calculate the statistical significance of changes between treatments and samples. To avoid the impact of noise, we analyzed only contacts that met a minimum baseline of anchor-to-anchor contacts (at least 8 contacts per billion contacts (CPB)) in one of the treatment conditions. Since TSS calling data (PRO-seq and coPRO-capped) was more abundant for K562 than Jurkat, when comparing K562 and Jurkat libraries we considered enhancer-promoter pairs with at least 8 CPB in both cell lines, to avoid ascertainment bias. The distribution of ratios between enhancer-promoter and background contacts in treated samples (Olaparib\TRP\FLV or dTAG treated cells) was compared to the median ratio in the respective control samples (Figs. 2D-E, 4C-D, 5B and S6A). For comparison between cell lines, enhancer-promoter contacts at promoters with increased gene body transcription and\or Pol II pausing signal in one cell line, were compared to their median at the other cell line (Figs. 3D-E and S4A,C).
Definition of CTCF binding sites
Contacts between CTCF binding sites were used as a control to determine whether the effects of a treatment were specific to enhancer-promoter contacts. We defined pairs of CTCF binding sites as CTCF motifs that were shown to bind CTCF based on ENCODE ChIP-seq data, within the same minimum and maximum allowed genomic distances as for enhancers and promoters. We focused only on CTCF sites that show no overlap with any dREG-defined TIR within 5kb.
PARP1 CUT&Tag
24h prior to the experiment, cells were split at a final concentration of 0.5 million cells/mL in fresh RPMI supplemented with 10% FBS and 1X penicillin streptomycin antibiotic. PARP1 CUT&Tag was done in two replicates with 250,000 K562 cells each for baseline control and anti-PARP1 CUT&Tag. Experiments were done according to the bench top CUT&Tag V.3 protocol92, using CUTANA™ Concanavalin A Conjugated Paramagnetic Beads (EpiCypher, SKU: 21-1401) and pAG-Tn5 for CUT&Tag (EpiCypher, SKU: 15-1017), with the following exceptions: (A) cells were not cross-linked and (B) NE1 buffer was used to extract and permeabilize nuclei. For anti-PARP1 CUT&Tag, we used an antibody against the N-terminal of PARP1 (Active Motif, AB_2793257). An anti-IgG antibody (Cell signaling, Normal Rabbit IgG #2729) was used as baseline control.
CUT&Tag data processing
fastq files were trimmed from Nextera adaptor-associated sequences using cutadapt93 and mapped to hg38 human genome assembly using bowtie294. We then used bedtools genomecov to generate bedgraph files and calculated background-normalized PARP1 binding using macs2 bdgcmp command95.
Statistical analysis
Throughout the manuscript, Mann-Whitney U test is used for independent samples, such as comparison of changes between different sets of genomic loci or pairs. Wilcoxon signed-rank test is used for paired samples, usually being the same loci\pairs compared between samples\conditions. For assessment of trends in our data, such as the changes in the difference in contacts between functional and nonfunctional pairs, or assessing the effect of changes in paused Pol II occupancy on changes in enhancer-promoter contacts, we used Pearson’s correlation coefficient (R). The confidence intervals for the medians throughout the manuscript were calculated using 1000 iterations of bootstrap.