ARTICLE https://doi.org/10.1038/s41467-022-31192-9 OPEN Cohesin couples transcriptional bursting probabilities of inducible enhancers and promoters Irene Robles-Rebollo 1,2,3, Sergi Cuartero 1,2,7, Adria Canellas-Socias1,2,8, Sarah Wells1,2, Mohammad M. Karimi 1,2,9, Elisabetta Mereu4, Alexandra G. Chivu 1,2,10, Holger Heyn 4, Chad Whilding1,2, Dirk Dormann1,2, Samuel Marguerat 1,2, Inmaculada Rioja5, Rab K. Prinjha5, Michael P. H. Stumpf 3,6, Amanda G. Fisher 1,2 & Matthias Merkenschlager 1,2 ✉ 1234567890():,; Innate immune responses rely on inducible gene expression programmes which, in contrast to steady-state transcription, are highly dependent on cohesin. Here we address transcrip- tional parameters underlying this cohesin-dependence by single-molecule RNA-FISH and single-cell RNA-sequencing. We show that inducible innate immune genes are regulated predominantly by an increase in the probability of active transcription, and that probabilities of enhancer and promoter transcription are coordinated. Cohesin has no major impact on the fraction of transcribed inducible enhancers, or the number of mature mRNAs produced per transcribing cell. Cohesin is, however, required for coupling the probabilities of enhancer and promoter transcription. Enhancer-promoter coupling may not be explained by spatial proxi- mity alone, and at the model locus Il12b can be disrupted by selective inhibition of the cohesinopathy-associated BET bromodomain BD2. Our data identify discrete steps in enhancer-mediated inducible gene expression that differ in cohesin-dependence, and suggest that cohesin and BD2 may act on shared pathways. 1 MRC London Institute of Medical Sciences, Institute of Clinical Sciences, Faculty of Medicine, Imperial College London, Du Cane Road, London W12 0NN, UK. 2 Institute of Clinical Sciences, Faculty of Medicine, Imperial College London, Du Cane Road, London W12 0NN, UK. 3 Theoretical Systems Biology Group, Imperial College London, London, UK. 4 CNAG-CRG, Centre for Genomic Regulation (CRG), The Barcelona Institute of Science and Technology (BIST), Barcelona, Spain. 5 Epigenetics RU, GlaxoSmithKline Medicines Research Centre, Stevenage, UK. 6 School of BioSciences and School of Mathematics and Statistics, University of Melbourne, Parkville, VIC 3010, Australia. 7Present address: Josep Carreras Leukaemia Research Institute, Badalona, Spain. 8Present address: IRB, Institute for Research in Biomedicine, Barcelona, Spain. 9Present address: Comprehensive Cancer Centre, School of Cancer & Pharmaceutical Sciences, Faculty of Life Sciences & Medicine, King’s College London, London, UK. 10Present address: Department of Molecular Biology and Genetics, Cornell University, New York, USA. ✉email: matthias.merkenschlager@lms.mrc.ac.uk NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications 1 ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 C ohesin cooperates with CTCF in organising the genome of the response to LPS (Fig. 1a and Supplementary Fig. 1b). into self-interacting domains that each contain a small smRNA-FISH counts were highly consistent between replicates number of genes and enhancers, thus potentially facil- (Supplementary Fig. 1c). The copy numbers of Il12b transcripts itating enhancer–promoter cooperation1–6. However, the acute increased with time after LPS activation (Supplementary Fig. 1d). depletion of CTCF or cohesin from mouse or human cell lines Transcript copy numbers of Il12b and other inducible genes were deregulates only a limited number of genes at steady-state, correlated with transcript abundance as determined by RNA-seq notably direct targets of CTCF promoter binding7, and boundary- (Fig. 1b). proximal genes that show altered transcriptional bursting8, respectively. By contrast to steady-state gene expression, loss of Regulation of inducible immune gene transcription. To identify cohesin severely disrupts inducible transcriptional programmes in the transcriptional parameters that underlie the transcriptional macrophages9 and neurons10, as well as oestrogen-responsive induction of immune genes in response to LPS we used nascent genes11. smRNA-FISH with probes complementary to transcribed introns. Inducible genes in macrophages are paradigmatic for the core Nascent smRNA-FISH captures the fraction of actively tran- transcriptional response to microbial signals12. Mechanistically, scribing alleles and the intensity of probe signal at transcription the engagement of anti-microbial pattern recognition receptors start sites (Burst intensity, Fig. 1c). We found that the inducible drives the recruitment of transcription factors, chromatin remo- expression of Il12b (Fig. 1d) and other inducible genes tested by dellers, and cofactors, including histone acetyltransferases. BET smRNA-FISH (Egr2, Prdm1, Ifnb1, Peli1 and Sertad2, Supple- (Bromo- and Extra-Terminal domain) epigenetic reader proteins mentary Fig. 2) was associated primarily with an increase in the interact with acetylated lysine residues via N-terminal tandem fraction of actively transcribing alleles. We compared the fraction bromodomains termed BD1 and BD2 (ref. 13) and recruit Med- of actively transcribing alleles, as directly quantified by smRNA- iator and RNA polymerase to initiate eRNA transcription12,14. FISH, with the burst frequency as inferred mathematically from Enhancer activation culminates in target gene transcription12 by the moments (mean and variance) of mature mRNA facilitating recruitment of the transcriptional machinery to distribution18,23,24,29 (Fig. 1e, left). We found good agreement, promoters14 and by increasing the frequency of promoter which indicates that burst fraction and burst frequency are cor- bursting15. Transcription is typically discontinuous and can be related in LPS-activated macrophages. In contrast, burst intensity described in terms of the frequency and the size of transcriptional as quantified by smRNA-FISH poorly correlated with burst size bursts16–19. Both the size20–24 and the frequency22,23,25–28 of as inferred by the moment of mature mRNA copy number dis- transcriptional bursts have been linked to inducible gene tribution (Fig. 1e, right). expression in different experimental systems. How cohesin con- To ask whether the expression of inducible genes was regulated tributes to the regulation of transcriptional bursting during by the frequency of transcript-expressing cells and/or the level of immune gene activation remains unknown. transcripts per cell, we performed scRNA-seq at 0, 2 and 8 h of Here we address the transcriptional parameters that underlie the LPS stimulation (Fig. 2). LPS increased the frequency of cells that activation of inducible genes and enhancers in macrophages and expressed inducible transcripts for each of 6 previously identified determine steps within the process that are cohesin-dependent or classes of inducible primary30 (Fig. 2, Bhatt classes A1/2, B) as well -independent. We use single-molecule RNA fluorescence in situ as secondary response genes30 (Fig. 2, Bhatt classes C, D, E, F). In hybridisation (smRNA-FISH) to quantitatively assess the fraction contrast to the frequency of transcript-expressing cells, the and intensity of enhancer and target gene transcription as well as number of transcripts per cell increased only for class C inducible the coordination between the transcriptional activity of enhancers genes at the early time point. Class C to F transcript numbers and promoters at the single-molecule level. smRNA-FISH indicates accumulated at the late time point (Fig. 2, right). Taken together, that inducible immune genes and the activity of associated these data indicate that the regulation of inducible immune genes enhancers are regulated primarily by an increase in the fraction of is mediated primarily by increasing the fraction of actively actively transcribing alleles, and that inducible enhancer activation transcribing alleles and the frequency of cells that express is tightly coupled to an increased probability of target gene tran- inducible transcripts. scription in wild-type cells. We find that cohesin is dispensable for activation-induced increases in the probability of enhancer bursting, and that the number of mature mRNAs produced per transcribed Regulation of inducible immune gene-enhancer eRNA tran- target gene allele is largely independent of cohesin. However, scription. Inducible immune gene expression is regulated by cohesin is critical for coupling the probabilities of enhancer and inducible enhancers12,31,32. To explore the cis-regulatory control target gene bursting. Finally, we present evidence that the role of of inducible immune genes we designed smRNA-FISH probes for cohesin in enhancer–promoter coupling may not be explained by inducible immune gene enhancers that were positioned between spatial proximity alone and demonstrate potential convergence of 10 and 200 kb from their target genes within the same topolo- cohesin and BD2, the second BET protein bromodomain, in the gically associating domains and in the same—or immediately activation of inducible genes. adjacent—contact domains in macrophage Hi-C maps: the functionally validated Il12b HSS1 enhancer33 and Ifnb1 L2 enhancer34, as well as enhancers associated with Egr2, Peli1, and Results Prdm1 (ref. 31). These enhancers were actively transcribed in LPS- To validate the role of transcription in LPS-inducible macrophage activated macrophages (Fig. 3a and Supplementary Fig. 3, GRO- gene expression we compared GRO-seq (as a measure of tran- seq), showed binding of constitutive (PU.1) as well as signal- scription) with RNA-seq (as a measure of transcript levels). RNA- activated transcription factors (STAT2 and IRF3), and inducible seq and GRO-seq were highly correlated (r = 0.83, P < 10e-16; acetylation of histone H3 lysine 27 (Fig. 3a and Supplementary Supplementary Fig. 1a). This indicates that LPS-induced macro- Fig. 3, H3K27ac, STAT2, IRF3). We used smRNA-FISH to phage genes are regulated at the level of transcription. To char- visualise transcriptional eRNA bursting of inducible macrophage acterise the transcription of inducible genes such as Il12b we enhancers (Fig. 3b). Quantification of eRNA burst parameters adopted a high throughput smRNA-FISH approach that allows indicated that activation-induced transcription of inducible the quantification of transcript copy numbers and transcriptional immune gene enhancers occurred primarily through modulation bursts in thousands of individual macrophages during the first 2 h of burst fraction (Fig. 3c), not burst intensity (Fig. 3d). 2 NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 ARTICLE a Untreated + LPS c Il12b > > Exon probes Intron probe RNA-seq GRO-seq chr11: (kb) 44,215 44,220 44,225 44,230 44,215 44,220 44,225 44,230 Il12b exon probe DAPI Il12b intron probe DAPI Il12b exon Il12b intron DAPI Il12b exon probe DAPI Il12b exon probe DAPI Il12b Il12b d 0.25 P = 0.012 20 P = 0.24 Burst fraction 0.20 TSS intensity P = 0.87 P = 1.6e-123 15 0.25 P = 0.27 10 smRNA-FISH 0.10 0.05 P = 6.7e-17 5 0.00 0 0 60 90 120 0 60 90 120 Time after LPS (min) Time after LPS (min) b e Burst frequency (moments) 8 r = 0.857 Burst size (moments) 80 log2 mean counts P = 2.26e-11 2.5 r = -0.25 Stat1 smRNA-FISH 6 P = 0.33 Il12b 2.0 60 Il12b Actb 4 Ifit1 Cd40 1.5 40 Peli1 Atf3 2 1.0 Cxcl10 Tnf r = 0.89 20 Ifnb1 0 0.5 P = 1.7e-06 0 2.5 5.0 7.5 10 0.0 0.2 0.4 0.6 3 4 5 6 7 log2 mean FPKM + 1 Burst fraction (smRNA-FISH) Burst intensity (smRNA-FISH) RNA-seq Fig. 1 Inducible gene transcription is regulated primarily by the fraction of actively transcribing alleles. a Induction of Il12b transcripts by macrophage activation as assessed by RNA-seq (top), GRO-seq (middle) and smRNA-FISH (bottom). LPS treatment times: 1 h for GRO-seq and 2 h for RNA-seq and smRNA-FISH. Representative of three independent biological replicates for RNA-seq, two independent biological replicates for GRO-seq and >3 independent biological replicates for smRNA-FISH. b Correlation between transcript abundance by RNA-seq (log2 FPKM) and log2 mean absolute transcript copy numbers by smRNA-FISH. r: Pearson correlation; P: p value in this and all subsequent figures. Statistical test: two-sided Pearson's product- moment correlation (Null hypothesis: correlation = 0). Il12b, Stat1 and Actb: three independent biological replicates, n = 8899 cells. Cd40, Atf3 and Tnf: 3 independent biological replicates, n = 12,372 cells. c Representative images illustrating the use of intron probes to localise transcriptional bursts. Greater than three independent biological replicates. d Burst fraction and burst intensity in resting and activated wild-type macrophages. In this and all subsequent box plot figures, boxes show the median and the upper and lower quartiles, whiskers show 1.5 of the interquartile range. P: p value. In this and all subsequent figures, we used a two-sided Tukey HSD test (ANOVA design formula ~ Sample + Replicate), adjusted for multiple testing for burst intensities and a two-sided Cochran–Mantel–Haenszel test (replicate as strata) with Bonferroni correction for multiple testing for burst fractions. N = 24,160 cells, four independent biological replicates per transcript. e Left: correlation between the fraction of actively transcribing alleles measured by smRNA-FISH and transcriptional burst frequencies inferred from the distribution of mature transcripts (see methods) in macrophages stimulated with LPS for 90 min (Il12b, Ifnb1, Peli1), 120 min (Cxcl10) or 180 min (Ifit1). Right: correlation between burst intensity measured by smRNA-FISH and burst size inferred from the distribution of mature transcripts (see methods). Intron probes were used for Il12b and Peli1, exon probes were used for Ifit1, Cxcl10, and for the intron-less Ifnb1 gene. r: Pearson correlation. P value derived from two-sided Pearson's product-moment correlation test as implemented in the cor.test function in R against the null hypothesis correlation is equal to 0. N = 18,494 cells, three to four independent biological replicates per transcript. Error bands: 95% confidence interval. Functional coupling of enhancer and promoter bursting. We bursting in individual cells, as illustrated for Il12b in Fig. 4d. To used smRNA-FISH to visualise transcriptional bursting of both ask whether eRNA and mRNA bursting is coordinated, we enhancers and target genes in response to LPS. Consistent with determined the frequencies of simultaneous eRNA and mRNA previous reports that enhancer activation can precede the tran- transcription. Simultaneous eRNA and mRNA bursts in response scription of target genes in macrophages32,35 and other cell to LPS activation were significantly more frequent than expected types11,36,37, the LPS-induced transcriptional burst fraction at the if enhancers and promoters were to burst independently (Fig. 4e). Il12b HSS1 enhancer peaked earlier than Il12b mRNA (Fig. 4a, Since mature macrophages are quiescent diploid cells, the P = 0.0036, two-way ANOVA). We recorded eRNA and mRNA occurrence of two enhancers or two promoter bursts in a single transcriptional bursts for five inducible enhancers and associated cell indicates that bursts occur on different chromosomes in genes in a total of 38,093 cells (Supplementary Fig. 4). Burst trans. The great majority (>95%) of such bi-allelic bursts were fractions of inducible enhancers and promoters were strongly >1 μm apart (Fig. 4f). We, therefore, used 1 μm as a threshold correlated (Fig. 4b). The number of cells that contained two below which we could assign enhancer and promoter bursts to enhancer or two promoter bursts was greater than expected based the same allele with >95% confidence (red line in Fig. 4f). In on random co-occurrence38 (Fig. 4c), suggesting that certain cells contrast to pairs of enhancer-enhancer and promoter-promoter were more likely than others to show enhancer or gene activity. bursts, pairs of enhancer–promoter bursts often occurred at smRNA-FISH probes for enhancers and their targets allowed distances <1 μm (Fig. 4f). Enhancer and promoter bursting on the exploration of the relationship between enhancer and promoter same allele occurred significantly more frequently than expected NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications 3 ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 3.3e-11 1.1e-7 9.0e-9 8.2e-9 8.8e-7 5.0e-5 4.5e-13 2.2e-10 3.2e-13 1.4e-11 1.0 1.0 LPS 6.4e-7 0.25 1.0 1.1e-6 7.8e-13 1.4e-11 10.0 9.5e-11 5.5e-9 - 2h 8h 1.0 1.0 log10 transcript expr. (CPM + 1) 8.2e-12 6.7e-8 4.7e-13 2.3e-10 1.3e-11 0.002 1.4e-11 Fraction of cells with transcripts 1.0 9.5 0.02 1.0 2.3e-11 0.06 1.0 3.2e-13 9.0 0.8 3.4e-8 8.5 LPS Genes - 0.6 8.0 2h 7.5 8h 0.4 7.0 0.2 Cells 6.5 0.0 6.0 0 2 4 6 8 10 12 ln(CPM+1) A1+2 B C D E F A1+2 B C D E F Primary Secondary Primary Secondary Fig. 2 Inducible gene transcription is driven primarily by the fraction of expressing cells. Left: heatmap of transcript levels (CPM) of inducible macrophage primary and secondary response genes30 detected by scRNA-seq in at least 20% of cells. Centre: Fraction of cells with detectable transcripts for each of six previously defined classes of inducible genes30 at 0, 2 and 8 h after LPS activation. Right: Transcript levels (CPM+1) per transcript-positive cell is shown for the same classes of inducible genes at 0, 2 and 8 h after LPS activation. Numbers represent adjusted P values for 0 versus 2 h and 0 versus 8 h. Two-sided Wilcoxon signed-rank test with Bonferroni correction for multiple testing. N = 682 cells. by chance (Fig. 4g). Taken together, these data show that transcripts per cell was similar for wild-type and Rad21−/− enhancer burst fractions are efficiently coupled to the burst macrophages after 2 h of LPS stimulation. fractions of target promoters, both at the cellular and the allelic To examine the transcriptional parameters associated with the level. reduced fraction of Rad21−/− macrophages expressing inducible genes, we used intron smRNA-FISH for Egr2, Prdm1, Il12b, Ifnb1 Impaired expression of inducible genes in cohesin-deficient and Peli1. In Rad21−/− macrophages, inducible genes displayed macrophages is the result of reduced transcription. To examine reduced burst fractions, while burst intensities were only mildly the role of the chromatin-associated cohesin complex in the affected (Fig. 5d). These results suggest that the reduced transcriptional bursting of enhancers and promoters, we deleted expression of inducible immune genes in Rad21−/− macrophages floxed Rad21 (Rad21lox/lox) alleles encoding the cohesin subunit at the population level was due primarily to a reduced probability RAD21 by 4-OHT-mediated activation of ERt2Cre in bone of transcription, which led to a reduced fraction of cells marrow-derived mature macrophages9. This approach removed expressing inducible immune genes. 90–99% of Rad21lox/lox alleles (Supplementary Fig. 5a) and ~90% of RAD21 protein (Supplementary Fig. 5b, c). Loss of RAD21 Transcriptional bursting is productive in the absence of protein was uniform across individual cells, as demonstrated by cohesin. To examine whether cohesin controls the number of immunofluorescence staining and confocal microscopy (Supple- mature transcripts, we compiled mRNA counts for cytoplasmic mentary Fig. 5d, e). Mature macrophages are suitable for transcripts and burst parameters for the inducible genes Cxcl10, mechanistic studies on the role of cohesin in gene regulation Ifit1, Ifnb1 and Il12b, and compared the log2 burst measures because their quiescent state avoids interference with essential cell against the log2 mRNA mean counts for each replicate in wild- cycle functions of cohesin39. The conditional deletion of the type and in Rad21−/− macrophages. This analysis comprised a cohesin subunit Rad21 in quiescent primary macrophages results total of 77,482 cells from 76 samples, 38 for wild-type and 38 for in a broad reduction in the expression of immune response genes Rad21−/− macrophages (Fig. 6a). Mature transcript copy num- as assessed by RNA-seq (ref. 9). Comparison of RNA-seq with bers were highly correlated with burst fractions both in wild-type GRO-seq showed that this reduction in inducible immune gene and Rad21−/− macrophages (Fig. 6a). expression occurs predominantly at the transcriptional level in We next addressed whether wild-type and Rad21−/− macro- both resting and LPS-activated macrophages (Fig. 5a). phages produced the same number of Cxcl10, Ifit1, Ifnb1 and Il12b transcripts for a given burst fraction over the range of burst Transcript copy numbers scale with burst fraction. We used fractions observed (Fig. 6b). For each gene, burst fraction and smRNA-FISH with exon probes to quantify the fraction of cells transcript copy number were highly correlated for both wild-type expressing inducible transcripts and the number of transcripts per (black) and Rad21−/− macrophages (red) (Fig. 6b). There were cell, as exemplified for Il12b in Fig. 5b. Rad21−/− macrophages no significant differences in the relationship between burst showed a reduced faction of transcript-expressing cells compared fraction and transcript copy number for any of the genes to wild-type macrophages (Fig. 5b). In contrast, the number of examined (Fig. 6b). This analysis showed that even though transcripts per cell was similar for Rad21−/− and wild-type Rad21−/− macrophages showed reduced transcriptional burst macrophages (Fig. 5b). This indicates that the reduced expression fractions, they produced mature mRNAs in numbers propor- of inducible transcripts in Rad21−/− macrophages is primarily tional to the observed fraction of transcribing alleles. We further due to a reduced frequency of cells expressing Il12b transcripts. compared burst parameters observed in resting and LPS-activated Consistent with these smRNA-FISH results, analysis by wild-type and cohesin-deficient cells with transcriptional activity scRNA-seq showed that the frequency of individual cells (GRO-seq, Supplementary Fig. 6a), and with the abundance of expressing previously defined inducible transcripts31 was reduced mature transcripts (RNA-seq, Supplementary Fig. 6b). Both in Rad21−/− compared to wild-type macrophages, both at transcriptional activity (GRO-seq) and the abundance of mature baseline and after LPS stimulation for 2 or 8 h (Fig. 5c). In transcripts (RNA-seq) strongly correlated with the fraction of contrast to the fraction of transcript-containing cells, the level of transcriptional bursting in both wild-type and Rad21−/− 4 NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 ARTICLE a Genes 1.5 chr10: 67,020 67,040 kb chr11: 44,210 44,220 kb chr4: 88,172 88,182 kb chr11: 2.0 20,980 21,020 kb chr10:44,200 0.75 44,300 kb 0.75 0.5 Peli1 Prdm1 Egr2 Il12b Ifnb1 GRO-seq resting -1.5 -0.75 -0.5 -2.0 -0.75 1.5 0.75 0.5 2.0 0.75 GRO-seq 1h LPS -1.5 -0.75 -0.5 -2.0 -0.75 RNA-FISH probes Egr2 Enh HSS1 Enh Il12b Ifnb1 L2 Enh Enh Peli1 Prdm1 Enh 400 150 400 400 400 PU.1 resting 400 400 PU.1 1h LPS 400 150 400 30 30 20 30 30 Stat2 resting 30 30 20 30 30 Stat2 1h LPS 100 50 150 100 100 Irf3 resting 100 50 150 100 100 Irf3 1h LPS 100 150 50 100 100 H3K27ac resting 100 150 50 100 100 H3K27ac 1h LPS Mφ enhancers Egr2 Enh. DAPI Il12b HSS1 Enh. DAPI Ifnb1 L2 Enh. DAPI Peli1 Enh. DAPI Prdm1 Enh. DAPI b Untreated + LPS 10 µm 10 µm c Egr2 Enh P = 5.2e-109 Il12b HSS1 Enh P = 4.8e-87 Ifnb1 L2 Enh P = 3.9e-11 Peli1 Enh Prdm1 Enh 0.10 0.50 P = 0.0 0.4 P = 2.5e-41 Burst Fraction 0.5 0.2 0.05 0.25 0.2 0.0 0.0 0.00 0.00 0.0 0 60 0 90 0 90 0 120 0 60 Time after LPS (min) Time after LPS (min) Time after LPS (min) Time after LPS (min) Time after LPS (min) d 60 P = 7.1e-10 10 P = 1.0 40 P = 0.43 10 P = 0.99 15 P = 2.5e-6 Burst 40 10 intensity 5 20 5 20 5 0 0 0 0 0 0 60 0 90 0 90 0 120 0 60 Time after LPS (min) Time after LPS (min) Time after LPS (min) Time after LPS (min) Time after LPS (min) Fig. 3 Inducible enhancer activation is linked to increased eRNA burst fractions. a Inducible enhancers associated with inducible genes in macrophages with GRO-seq reads from resting and 1 h LPS-activated macrophages, ChIP-seq for PU.1 and STAT2, IRF3 and H3K27ac. Macrophage enhancers were defined previously31. The positions of smRNA-FISH probes are magenta. b smRNA-FISH Images of eRNA bursts in untreated and LPS-activated macrophages with probes for the Egr2 enhancer (0 and 60 min LPS), the Il12b HSS1 and Ifnb1 L2 (0 and 90 min LPS), and the Peli1 enhancer (0 and 120 min LPS). smRNA-FISH probes are in magenta. At least three independent biological replicates per gene and condition, n = 64,682 cells. c eRNA burst fractions in resting and LPS-activated macrophages. At least three independent biological replicates per gene and condition, n = 64682 cells. d eRNA burst intensities in resting and LPS-activated macrophages. P values were determined by a two-sided Tukey HSD test (ANOVA design formula ∼ Sample + Replicate) and adjusted for multiple testing. At least three independent biological replicates per gene and condition, n = 64,682 cells. macrophages (Supplementary Fig. 6). Taken together, these data transcription factors are dysfunctional in cohesin-deficient indicate that cohesin is required for the initiation of transcrip- macrophages9. This is illustrated here for the Ifnb1 L2 enhancer tional bursting, but, once initiated, the number of transcripts (Fig. 7a). However, of 1112 intergenic enhancers31 that were produced is a reflection of burst fraction, both in the presence and detectably transcribed by GRO-seq and showed significant the absence of cohesin. induction after 1 h of LPS treatment in wild-type macrophages, the great majority (1048 of 1112 or 94.2%) were transcribed to at Cohesin couples transcriptional burst fractions of inducible least wild-type levels in Rad21−/− macrophages (adj. P > 0.05). enhancers and their target genes. A subset of inducible enhan- Enhancers associated with Il12b, Prdm1, Peli1 and Egr2 were cers characterised by the binding of Irf- and Stat- family efficiently induced in Rad21−/− as well as wild-type macrophages NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications 5 ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 a Il12b HSS1 enhancer bursts Il12b intron bursts b 0 0.35 P = 5.1e-10 0.10 P = 0.23 −2 log2 gene burst fraction 0.30 P = 1.4e-63 P = 1.7e-41 0.08 −4 0.25 Burst fraction Burst fraction P = 8.2e-21 0.20 0.06 −6 0.15 0.04 P = 1.9e-27 −8 0.10 0.02 −10 r = 0.609 0.05 P = 0.000167 0.00 0.00 −12 0 60 90 120 0 60 90 120 −6 −4 −2 0 Time after LPS (min) Time after LPS (min) log2 enhancer burst fraction c Il12b HSS1enh. bursts Il12b intron bursts Egr2 Il12b Untreated + LPS with 2 prom. bursts P = 4.31e-7 with 2 enh. bursts P = 6.39e-25 Prdm1 Fraction of cells Fraction of cells 0.10 0.02 Ifnb1 Peli1 P = 3.7e-5 0.05 0.01 NaN P = 3.9e-6 0.00 0.00 P = 3.4e-11 0 90 Observed 0 90 f 15 P = 0.95 Expected n = 73 n = 261 n = 139 Distance in μm 10 d Enhancer probe mRNA probe Merge 5 0 Il12b intr. HSS1 enh. HSS1 enh. Il12b intr. HSS1 enh. Il12b intr. 10 µm <1 μm: 1.2% 4.1% 46.8% Il12b HSS1 DAPI Il12b intr. DAPI Il12b HSS1 Il12b intr. DAPI e Il12b HSS1- Il12b Ifnb1 L2 - Ifnb1 Prdm1 Enh - Prdm1 P = 2.2e-10 P = 0.0001 P = 0.0018 0.06 0.025 0.015 P = 0.15 0.020 P = 0.0004 Co-occurrence of Observed 0.04 0.010 0.015 Expected enh-prom bursts 0.02 0.010 at the cell level NaN 0.005 0.005 0.00 0.000 0.000 0 90 0 90 0 60 Time after LPS (min) Time after LPS (min) Time after LPS (min) g Il12b HSS1- Il12b Ifnb1 L2 - Ifnb1 Prdm1 Enh - Prdm1 0.03 P = 1.48e-9 P = 0.003 P = 0.004 0.006 0.010 P = 0.26 0.008 P = 2.2e-5 Co-occurrence of 0.02 0.004 0.006 Observed enh-prom bursts 0.004 Expected at the allelic level 0.01 0.002 (<1 μm) NaN 0.002 0.00 0.000 0.000 0 90 0 90 0 60 Time after LPS (min) Time after LPS (min) Time after LPS (min) (Fig. 7a). The Il12b HSS1 enhancer, Prdm1 enhancer and Peli1 Prdm1 and Peli1 were markedly reduced in Rad21−/−compared enhancer displayed higher burst fractions in Rad21−/− compared to wild-type macrophages (Fig. 7b). The promoter burst fraction to wild-type macrophages (Fig. 7a). The Egr2 enhancer was sig- of Egr2 was more strongly reduced than the Egr2 enhancer burst nificantly induced in Rad21−/− macrophages, albeit to lower fraction in Rad21−/− macrophages. levels than in wild-type macrophages (Fig. 7a). In contrast to To quantify the observed defects in promoter activation we enhancer burst fractions, the promoter burst fractions of Il12b, analysed the relationship between transcriptional burst fractions 6 NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 ARTICLE Fig. 4 Functional coupling of enhancer and promoter bursting. a Transcriptional burst fractions of the Il12b HSS1 enhancer (left, 3–4 independent biological replicates, n = 14,759 cells) and Il12b (right, four independent biological replicates, n = 17,038 cells). b log2-transformed burst fractions of the indicated enhancers and promoters in resting (UT) and LPS-activated macrophages (LPS). Two-way ANOVA formula (BF_Gene ~ BF_Enhancer*Gene). P value (P) for whether there is a relationship between the bursting of enhancers and genes (frequency_enhancer) = 1.94e-09. At least three independent biological replicates per gene and condition, n = 38093 cells. Error band: 95% confidence interval. c Bi-allelic enhancer (left) and promoter bursts (right) occur in individual cells more often than expected based on the burst fractions of enhancers and promoters (3–4 independent biological replicates, n = 14759 cells for Il12b HSS1 and n = 17038 cells for Il12b). Cochran–Mantel–Haenszel test (replicate as strata) with Bonferroni correction for multiple testing, two-sided considering a random population of the same size where the sampled frequency equals the expected frequency in case of bursting events happening independently. d Simultaneous application of enhancer (magenta) and target gene probes (green) illustrated for the Il12b locus. The position of probes is shown in Fig. 3a. More than eight biological and four technical replicates. e Experimentally observed (observed) co-occurrence of eRNA and mRNA transcriptional bursts over co-occurrence expected by chance (expected) at the level of individual cells. Cochran–Mantel–Haenszel test (replicate as strata) with Bonferroni correction for multiple testing, two-sided. At least three independent biological replicates per gene and condition, n = 38093 cells. f The distance distribution of Il12b intron - Il12b intron pairs, Il12b HSS1 enhancer—Il12b HSS1 enhancer pairs and Il12b intron—Il12b HSS1 enhancer pairs is plotted. A horizontal red line indicates a distance of 1 μm. The percentages of Il12b intron— Il12b intron pairs (1.2%), Il12b HSS1 enhancer— Il12b HSS1 enhancer pairs (4.1%) and Il12b intron—Il12b HSS1 enhancer pairs below 1 μm (46.8%) is indicated. Tukey HSD test (ANOVA design formula distance ~ pair + replicate), two-sided. g Box plots of experimentally observed (observed) co-occurrence of eRNA and mRNA transcriptional bursts over co-occurrence expected by chance (expected) at the level of individual alleles (i.e. at distances below 1 μm). Cochran–Mantel–Haenszel test (replicate as strata) with Bonferroni correction for multiple testing, two-sided. At least three independent biological replicates per gene and condition, n = 38,093 cells. of genes and enhancers in wild-type and cohesin-deficient cells. challenges: direct enhancer–promoter looping is neither strictly In wild-type macrophages, enhancer activation was linked to required, nor is it sufficient for the activation of target genes40–44, efficient promoter activation (black in Fig. 7c), indicative of and many loci that undergo enhancer–promoter looping during functional coupling between enhancers and promoters. Rad21−/− active transcription can do so in the absence of cohesin across macrophages showed a blunted increase of target gene transcrip- significant genomic distances of ~100 kb (refs. 10,45–49). Based on tional burst fractions, despite elevated burst fractions of the Il12b these and other considerations discussed more fully below, we set HSS1, Prdm1, Peli1 and Egr2 enhancers in response to activation out to test the validity of the spatial proximity model of cohesin signals (red in Fig. 7c), indicating that enhancer–promoter action in our experimental system. coupling was substantially weakened (Il12b, Peli1, Egr2) or lost Analysis of smRNA-FISH data showed that enhancer and completely (Prdm1) in the absence of cohesin (Fig. 7d). promoter bursts occurred over a wide range of distances (Fig. 8, We next asked whether cohesin has a role in coupling the top). Even in the sub-micron range that contains enhancers and activity of inducible enhancers to target gene transcription promoters on the same allele (see above), cumulative genome-wide. To this end, we analysed the LPS-induced enhancer–promoter distances showed no—or at most minor— transcription of canonical macrophage enhancers31 by GRO-seq detectable differences between wild-type and cohesin-deficient (enhancer log2 fold-change, Fig. 7e) and the expression of cells (Fig. 8, bottom). The Egr2 enhancer was slightly closer to the protein-coding genes (RNA-seq, gene log2 fold-change, Fig. 7e) Egr2 target gene in cohesin-deficient macrophages (genomic that were associated by 3D genome organisation at the level of distance ~61 kb), and vice versa for Peli1 (genomic distance TADs or contact domains. Wild-type macrophages showed a ~43 kb). No significant differences were found for the Prdm1 and correlation between the transcription of enhancers and genes Il12b enhancer-gene pairs (~198 and ~10 kb, respectively). These located in the same TADs (Fig. 7e left, r = 0.22 P = 4.09e-05) as data show that transcriptional bursts can occur across a range of well as contact domains (Fig. 7e right, r = 0.32, P = 5.67e-05). In different enhancer–promoter distances. Furthermore, within the contrast, as observed for the genes investigated by smRNA-FISH, resolution limits of our experimental setup, enhancer–promoter cohesin-deficient macrophages showed a reduced correlation distances are similar in wild-type and cohesin-deficient macro- between the transcription of enhancers and genes located in the phages, suggesting that enhancer–promoter distances alone may same TADs (Fig. 7e left, r = 0.12, P = 0.028) or contact domains not explain the role of cohesin in the functional coupling of (Fig. 7e right, r = 0.14, P = 0.01). The difference in the correlation enhancers and promoters. between enhancer and gene transcription between wild-type and cohesin-deficient macrophages was significant at the level of both Non-additive effects of cohesin and BD2 perturbations. We TADs (Fig. 7e left, P = 0.0006) and contact domains (Fig. 7e considered a potential involvement of BET proteins, since (i) human right, P = 1.99e-05). This analysis indicates that cohesin couples mutations in cohesin genes and in BRD4 cause related the activity of inducible macrophage enhancers to inducible pathology50–52, and (ii) either BET inhibition or targeted degrada- macrophage genes genome-wide. Taken together with the tion of BRD4 abrogate transcription but leave enhancer–promoter smRNA-FISH data, these results show that the loss of cohesin contacts intact53. We used the BD1-selective small molecule inhibitor uncouples transcriptional burst fractions of inducible genes from GSK778, which largely phenocopies pan-BET inhibitors, as well as transcriptional burst fractions of inducible enhancers in innate the BD2-selective inhibitor GSK046, which has more limited effects immune cells and that enhancer-mediated gene activation on steady-state transcription54. BD1 and BD2 inhibitors primarily involves discrete steps characterised by selective reliance on affected burst fractions rather than burst intensities (Fig. 9a). Selec- cohesin. tive inhibition of either BD1 or BD2 significantly blunted the LPS- induced increase of Egr2 enhancer and Egr2 promoter burst fractions Transcriptional bursts occur across a range of enhancer– in wild-type macrophages (Fig. 9a). Interestingly, selective inhibition promoter distances. A straightforward interpretation of these of BD2 allowed for full LPS induction of burst fractions at the Il12b results would be that cohesin mediates enhancer–promoter cou- HSS1 enhancer, but significantly inhibited Il12b promoter burst pling simply by decreasing the proximity between enhancers and fractions (Fig. 9a). Selective inhibition of BD2 therefore in effect promoters in 3D nuclear space. Although appealing, the view of uncoupled the Il12b HSS1 enhancer from the Il12b promoter, cohesin as a mere provider of spatial proximity faces significant reminiscent of the uncoupling of enhancer and promoter bursting NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications 7 ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 a b Wild-type DAPI Il12b Rad21-/-DAPI Il12b Fraction of Il12b+ cells 3.2e-58 0.6 log2 FC GRO-seq Rad21-/- vs WT Inducible genes 0.5 1.2e-24 Untreated 2 0.4 1 0.3 1.0 0 0.2 0.07 0.1 −1 0.0 −2 0 60 90 120 −3 Time after LPS (min) r = 0.72 # of Il12b transcripts 200 −4 P = 2.3e-50 0.89 175 0.99 −5 150 LPS −8 −6 −4 −2 0 2 4 125 log2 FC RNA-seq Rad21-/- vs WT 100 0.89 75 50 1.0 25 Wild-type 0 0 60 90 120 Rad21-/- Time after LPS (min) c Untreated 2h LPS 8h LPS 1.0 6.7e-33 2.7e-28 9.0 4.0e-52 Transcript levels 2.1e-9 0.8 8.5 (CPM + 1) Fraction of 5.8e-39 7.0e-30 Genes cells with 0.6 8.0 transcripts 0.4 7.5 0.2 7.0 Wild-type Rad21-/- 0.0 6.5 Cells 0 2 8 0 2 8 0 2 4 6 8 10 12 Time after LPS (h) Time after LPS (h) ln(CPM+1) d Egr2 P = 0.0 0.10 Prdm1 Il12b Ifnb1 1.0 Peli1 P = 0.0 P = 1.2e-13 P = 2.4e-21 P = 9.2e-8 P = 3.3e-89 0.6 0.08 P = 4.5e-6 0.8 P = 0.0007 0.10 0.10 P = 0.0 0.4 P < 2.2e-16 0.06 0.6 Burst fraction 0.04 0.05 0.05 0.4 0.2 P = 0.025 0.02 P = 1.0 0.2 0.0 0.00 0.00 0.00 0.0 0 60 0 60 0 90 0 90 0 120 100 20 40 20 P = 0.036 80 P = 5.9e-9 P = 0.93 P = 0.015 P = 0.78 10 15 30 P = 0.95 15 Burst P = 0.94 P = 1.0 60 P = 0.65 P = 1.0 P = 1.0 P = 2.7e-9 intensity 10 20 10 40 5 20 5 10 5 WT 0 0 0 0 0 Rad21-/- 0 60 0 60 0 90 0 90 0 120 Time after LPS (min) Time after LPS (min) Time after LPS (min) Time after LPS (min) Time after LPS (min) Fig. 5 Deregulation of inducible genes in the absence of cohesin is due primarily to a reduced probability of transcription. a Loss of cohesin alters the expression of inducible genes at the transcriptional level. One hour LPS for GRO-seq, 2 h LPS for RNA-seq. r Pearson correlation. Two-sided Pearson's product-moment correlation (null hypothesis: correlation = 0). Three biological replicates for RNA-seq and two biological replicates for GRO-seq. b The fraction of cells expressing Il12b by smRNA-FISH (top left) and the number of transcripts per cell (bottom left) are plotted for cells expressing at least five transcripts at the indicated time points after LPS activation of wild-type and Rad21−/− macrophages. Numbers represent adjusted P values. Wilcoxon signed-rank test with Bonferroni correction for multiple testing. Example images are shown on the right. N = 40,941 cells, four independent biological replicates. c Left: heatmap of inducible macrophage genes30 detected by scRNA-seq in at least 20% of wild-type (grey) or Rad21−/− macrophages (red). Middle: Fraction of cells with detectable transcripts for inducible genes. Right: The expression level of LPS-inducible transcripts in cells with detectable transcripts is plotted as ln(CPM +1). Numbers represent adjusted P values. Two-sided Wilcoxon signed-rank test with Bonferroni correction for multiple testing. N = 1362 cells. d Burst fraction and burst intensity measurements for inducible immune genes in wild-type (black) and Rad21−/− macrophages (red). N = 11,8058 cells, 3–7 independent biological replicates. 8 NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 ARTICLE a Wild-type Rad21-/- 6 6 log2 mean mRNA counts log2 mean mRNA counts 4 4 2 2 mRNA copy # versus 0 0 burst fraction −2 −2 −4 −4 r = 0.77 r = 0.71 −6 P = 6.2e-8 −6 P = 1.5e-6 −10 −5 0 −10 −5 0 Ifnb1 log2 burst fraction log2 burst fraction Il12b Wild-type Rad21-/- Ifit1 Cxcl10 6 6 log2 mean mRNA counts log2 mean mRNA counts 4 4 mRNA copy # 2 2 versus 0 0 burst intensity −2 −2 −4 −4 r = 0.035 r = -0.019 −6 P = 0.84 −6 P = 0.91 0 1 2 3 4 5 0 1 2 3 4 5 log2 burst intensity log2 burst intensity b log2 mean mRNA counts Cxcl10 Ifit1 Ifnb1 Il12b 6 6 6 6 4 4 4 4 2 2 2 2 mRNA copy # versus 0 0 0 0 burst fraction −2 −2 −2 −2 r = 0.94 r = 0.93 r = 0.9 −4 r = 0.63 −4 −4 −4 P (correlation) = 3.6e-7 P (correlation) = 3.2e-5 P (correlation) =4.0e-7 P (correlation) =2.2e-4 −6 P (wt vs ko) = 0.36 −6 P (wt vs ko) = 0.45 −6 P (wt vs ko) = 0.18 −6 P (wt vs ko) = 0.46 WT −10 −5 0 −10 −5 0 −10 −5 0 −10 −5 0 Rad21-/- log2 burst fraction log2 burst fraction log2 burst fraction log2 burst fraction log2 mean mRNA counts 6 6 6 r = -0.51 6 P (correlation) = 0.051 4 4 4 P (wt vs ko) = 0.85 4 mRNA copy # 2 2 2 2 versus 0 0 0 0 burst intensity −2 −2 −2 −2 −4 r = 0.046 −4 r = 0.5 −4 −4 r = -0.13 P (correlation) = 0.87 P (correlation) = 0.12 P (correlation) = 0.51 −6 P (wt vs ko) = 0.61 −6 P (wt vs ko) = 0.99 −6 −6 P (wt vs ko) = 0.86 0 1 2 3 4 5 0 1 2 3 4 5 0 1 2 3 4 5 0 1 2 3 4 5 log2 burst intensity log2 burst intensity log2 burst intensity log2 burst intensity Fig. 6 Transcriptional bursting is productive in the absence of cohesin. a Transcript copy number correlates with transcriptional burst fraction in the absence of cohesin. Transcript copy number, a burst fraction (top) and burst intensity (bottom) were determined for wild-type (left) and Rad21−/− macrophages (right) as determined by smRNA-FISH using intronic probes for Il12b and exon probes for Ifit1, Cxcl10, and Ifnb1. Shown are individual replicates for each gene in wild-type and Rad21−/− macrophages after 0, 30, 60, 90, 120 or 180 min of LPS stimulation. The overall correlation was tested. r Pearson correlation, error bands: 95% confidence interval. Two-sided Pearson's product-moment correlation (null hypothesis: correlation = 0). N = 78,747 cells, 3–4 independent biological replicates per transcript. b Transcript copy number, a burst fraction (BF, top) and burst intensity (BI, bottom) for wild-type (black) and Rad21−/− macrophages (red) as in Fig. 4d. For each gene we tested the overall r: Pearson correlation, error bands: 95% confidence interval. Two-sided Pearson's product-moment correlation against the null hypothesis: correlation = 0). Statistical test for differences between wild-type and Rad21−/− macrophages: two-way ANOVA in R with interaction formula l log2.mean.mRNA.counts ∼ log2. Burst fraction * Genotype + log2. Burst intensity * Genotype. N = 78,747 cells, 3–4 independent biological replicates per transcript. NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications 9 ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 a 0.4 Il12b HSS1 enhancer Prdm1 enhancer 0.8 Peli1 enhancer Egr2 enhancer Ifnb1 L2 enhancer P = 3.3e-10 0.4 P = 0.83 P = 8.5e-221 1.0 P = 2.1e-106 0.10 P = 3.4e-11 0.3 0.6 P = 6.3e-17 P = 6.1e-86 ~ P = 1.0 Enhancer 0.2 0.4 0.2 0.5 0.05 burst fraction P = 0.0091 P = 1.7e-128 0.1 0.2 0.0 0.0 0.0 0.0 0.00 0 90 0 60 0 120 0 60 0 90 WT Time after LPS (min) Time after LPS (min) Time after LPS (min) Time after LPS (min) Time after LPS (min) Rad21-/- b Il12b promoter Prdm1 promoter Peli1 promoter 1.00 Egr2 promoter Ifnb1 promoter 0.15 0.075 P = 1.2e-13 P = 0.0 P = 1.0e-16 P = 4.5e-6 1.0 P = 0.0 P = 3.8e-60 0.75 0.10 Promoter 0.10 0.050 P = 0.0 P = 1.0 P = 2.6e-135 0.50 burst fraction 0.5 0.05 0.025 0.05 N/A 0.25 0.00 0.000 0.0 0.00 0.00 WT 0 90 0 60 0 120 0 60 0 90 Rad21-/- Time after LPS (min) Time after LPS (min) Time after LPS (min) Time after LPS (min) Time after LPS (min) c Il12b Prdm1 1.00 Peli1 Egr2 d 0.10 P = 0.0042 P = 0.0017 P = 0.0015 P = 0.0250 Enhancer-promoter coupling Prom. burst fraction 0.06 0.75 0.50 0.05 high Promoter 0.50 burst fraction 0.04 0.25 0.00 0.25 low WT 0.02 0.00 + LPS −0.05 0.00 Enh. burst fraction Rad21-/- 0.1 0.2 0.3 0.1 0.2 0.0 0.2 0.4 0.2 0.4 0.6 + LPS Enh. burst fraction Enh. burst fraction Enh. burst fraction Enh. burst fraction e 6 By TADs 6 By contact domains P = 0.0006 P = 1.99e-05 Gene log2 fold-change Gene log2 fold-change 5 5 p = 5.67e-05 4 4 P = 4.09e-05 3 3 2 2 P = 0.10 P = 0.028 1 1 0 Wild-type Wild-type 0 0 1 2 3 4 5 6 Rad21-/- 0 1 2 3 4 5 6 Rad21-/- Enhancer log2 fold-change Enhancer log2 fold-change Fig. 7 Cohesin couples transcriptional burst fractions of inducible genes and inducible enhancers. a Transcriptional burst fractions of the indicated enhancers in wild-type (grey) and Rad21−/− macrophages (red). N = 135,185 cells, four independent biological replicates. b Transcriptional burst fractions of the indicated promoters in wild-type (grey) and Rad21−/− macrophages (red). N = 135,185 cells, four independent biological replicates. c Transcriptional burst fractions of the indicated promoters and enhancers in wild-type (grey) and Rad21−/− macrophages (red). Statistical test for the interaction: two-way ANOVA, two-sided (design formula: burst fraction gene ∼ burst fraction enhancer*Genotype), P values < 0.05 indicate that the relationship between enhancer and promoter burst fractions is different between Rad21−/− and wild-type macrophages. N = 135,185 cells, four independent biological replicates. Error bands: 95% confidence interval. d Schematic: Relationship between transcriptional burst fractions of inducible genes and inducible enhancers. e Correlation between the LPS-induced transcription of inducible enhancers (GRO-seq enhancer log2 fold-change 1 h after LPS) and the transcription of protein-coding genes located in the same TADs (left) or contact domains (right) as the enhancers (RNA-seq gene log2 fold-change 2 h after LPS) in wild- type (grey) and Rad21−/− macrophages (red). For TADs and contact domains with more than one enhancer, the mean fold-change is shown. TADs: Pearson correlation wild-type = 0.22, P = 4.09e-05. Pearson correlation Rad21−/− = 0.12, P = 0.028. Two-way ANOVA Enhancer log2 fold- change:Genotype P = 0.00061. Contact domains: Pearson correlation wild-type = 0.32, P = 5.67e-05. Pearson correlation Rad21−/− = 0.14, P = 0.10. Two- way ANOVA Enhancer log2 fold-change:Genotype 1.99e-05. Error bands: 95% confidence interval. observed in cohesin-deficient macrophages. To ask whether cohesin pathway linking mutations in cohesin-related genes to a CdLS- and BD2 might converge on a shared pathway we explored the associated Y430C mutation in the BD2 of BRD4 (Ref. 55). In impact of selective BD inhibitors on burst fractions in wild-type and rapidly dividing mouse embryonic stem cells, corresponding BD2 cohesin-deficient macrophages. Selective BD1 inhibition reduced mutations caused a delay in cell cycle progression, susceptibility promoter burst fractions in wild-type and in both cohesin-deficient to DNA damage, and delayed DNA damage repair55. We, macrophages. Interestingly, however, selective inhibition of BD2 therefore, examined cell cycle- and DNA damage-related path- reduced Egr2 and Il12b promoter burst fractions in wild-type, but ways in wild-type and Rad21−/− macrophages. not in cohesin-deficient macrophages (Fig. 9b). Hence, BD2 inhi- Genes related to DNA damage, DNA repair and the cell cycle bition reduced Egr2 and Il12b transcriptional burst fractions in the were downregulated in response to LPS activation in both wild-type presence of cohesin, but had no additive effect to the loss of cohesin and Rad21−/− macrophages (Supplementary Fig. 7a). Other than in Rad21−/− macrophages. This finding suggests the possibility that that, Rad21−/− macrophages showed no significant enrichment of cohesin and BD2 may act in a shared pathway. gene ontology terms related to DNA damage, DNA repair or the cell Cohesin is essential for genome integrity during the cell cycle39, cycle in either resting or LPS-activated conditions (Supplementary and DNA damage and repair have been proposed as one potential Fig. 7a). Gene set enrichment analysis of Hallmarks related to DNA 10 NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 ARTICLE Egr2 enhancer to Egr2 Prdm1 enhancer to Prdm1 Il12b HSS1 enh. to Il12b Peli1 enhancer to Peli1 1.0 1.0 1.0 1.0 Cumulative probability 0.8 0.8 0.8 0.8 0.6 0.6 0.6 0.6 Wild-type 0.4 Wild-type Wild-type Wild-type 0.4 Rad21-/- 0.4 Rad21-/- 0.4 Rad21-/- Rad21-/- P = 1.7e-6 0.2 P = 0.71 P = 0.4 P = 8.4e-18 0.2 0.2 0.2 0.0 0.0 0.0 0.0 0.0 2.0 4.0 6.0 8.0 10.0 0.0 2.0 4.0 6.0 8.0 10.0 0.0 2.0 4.0 6.0 8.0 10.0 0.0 2.0 4.0 6.0 8.0 10.0 Gene-enh. distance (μm) Gene-enh. distance (μm) Gene-enh. distance (μm) Gene-enh. distance (μm) 1μm 1μm 1μm 1μm zoom zoom zoom zoom 1μm 1μm 1μm 1μm 1.0 1.0 1.0 Wild-type 1.0 Wild-type Cumulative probability Rad21-/- Rad21-/- 0.8 0.8 0.8 0.8 0.6 0.6 0.6 0.6 Wild-type Wild-type 0.4 Rad21-/- 0.4 Rad21-/- 0.4 0.4 P = 0.001 P = 0.46 P = 0.45 P = 3.0e-15 0.2 0.2 0.2 0.2 0.0 0.0 0.0 0.0 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 Gene-enh. distance (μm) Gene-enh. distance (μm) Gene-enh. distance (μm) Gene-enh. distance (μm) Fig. 8 Similar enhancer–promoter distances in wild-type and cohesin-deficient macrophages. Cumulative distances between enhancer and promoter smRNA-FISH signals in LPS-activated macrophages at the indicated loci. Number of wild-type cells: 8382 (Egr2), 3985 (Prdm1), 16,604 (Il12b), 19,864 (Peli1), at least three biological replicates per locus. Number of Rad21−/− cells: 7658 (Egr2), 3992 (Prdm1), 7352 (Il12b), 18,600 for Peli1), at least three biological replicates per locus. Top: all observed distances. Number of wild-type gene-enhancer pairs: 1754 (Egr2), 59 (Prdm1), 370 (Il12b), 5004 (Peli1). Number of Rad21−/− gene-enhancer pairs: 128 (Egr2), 35 (Prdm1), 43 (Il12b), 4288 (Peli1). Bottom: Distances up to 1 μm. Number of wild-type gene- enhancer pairs: 534 (Egr2), 23 (Prdm1), 116 (Il12b), 1440 (Peli1). Number of Rad21−/− gene-enhancer pairs: 52 (Egr2), 20 (Prdm1), 19 (Il12b), 956 (Peli1). P values < 0.05 indicate that two cumulative distributions are different by the two-sided Kolmogorov–Smirnov test. damage and the cell cycle56 (MSigDB) revealed no major differences Inducible immune gene expression during the LPS response of between wild-type versus Rad21−/− macrophages (Supplementary macrophages is transient and involves a chain of hierarchical events Fig. 7b). Rad21−/− macrophages showed reduced expression of that alter mRNA distributions dynamically. Therefore, these dis- G2M checkpoint and E2F targets compared to wild-type at 2 h after tributions cannot be assumed to be in steady-state, especially at LPS activation, and marginally higher expression at late time points early time points30. During this transient response, the expression not used in our smRNA-FISH experiments. Genes that are of many inducible immune genes is biphasic38, indicating that only downregulated during the response to UV damage were marginally a subset of cells activates the expression of a given gene during any upregulated in resting Rad21−/− macrophages. Inflammatory one activation cycle. In this setting, cell to cell variability likely pathways are known to be downregulated in Rad21−/− macro- reflects the probability of a gene turning on at least once during the phages and are shown as a control (Supplementary Fig. 7b). The activation cycle, and the smRNA-FISH burst fraction captures the absence of DNA damage or cell cycle effects is consistent with the probability of alleles initiating transcription. Despite the observed quiescent state of mature macrophages (~90% G1) and the lack of increase in enhancer burst fractions in response to activation, loss of p53 target gene activation in Rad21−/− macrophages noted in our cohesin reduced coupling of enhancer and target gene bursting. previous work9. Taken together, these results indicate that the Accordingly, scRNA-seq data and smRNA-FISH counts in cohesin- pathways shared by cohesin loss and BD2 inhibition in macrophages deficient macrophages show a reduction in the frequency of are unrelated to DNA damage. transcript-expressing cells, rather than in transcript copy number per cell. These observations are reminiscent of models where Discussion enhancer function affects the probability rather than the level of Here we dissect the regulation of burst parameters in inducible gene expression57 and support the notion that cohesin couples gene transcription and the cooperation between enhancers and target gene transcription to enhancer activity. their target genes. We show that inducible immune genes and A simple interpretation of these results would be that cohesin enhancers in primary mouse macrophages are regulated pre- mediates enhancer–promoter coupling by increasing the proxi- dominantly by the fraction of transcriptional bursting, while burst mity between enhancers and promoters in 3D nuclear space. intensities remain largely unchanged. We extend previous Consistent with this model, enhancers are major determinants of observations that inducible immune genes show heterogeneous promoter burst frequency15,58,59, DNA looping is known to expression at the single cell level by scRNA-seq (ref. 38) by predominantly affect burst frequencies58–60, and cohesin loss has demonstrating that cohesin loss altered the frequency of cells that previously been linked to altered promoter burst frequencies in express inducible immune gene transcripts, rather than the level HCT-116 cells8. While appealing, the interpretation of cohesin as of expression per cell. Consistent with these results, smRNA- a mere enhancer–promoter glue faces several challenges: First, FISH showed that loss of cohesin predominantly reduced burst while enhancer–promoter looping is sufficient to activate tran- fractions, while burst intensities were largely unaffected. scription of certain genes57,61,62, direct enhancer–promoter NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications 11 ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 a 50 Egr2 enhancer Egr2 P = 5.1e-06 Il12b HSS1 enhancer P = 2.2e-09 Il12b P = 2.4e-09 P = 0.0024 80 35 180 P = 1.0 Gene burst intensity Gene burst intensity P = 0.00027 Enh. burst intensity P = 0.76 P = 0.71 Enh. burst intensity 70 160 40 P = 8.1e-09 P = 5.1e-11 P = 1.0 30 P = 0.00081 P = 0.12 P = 0.019 60 P = 3.1e-05 P = 0.81 140 25 30 50 120 20 100 40 20 15 80 30 60 10 10 20 40 10 5 20 0 0 0 0 LPS - + + + LPS - + + + LPS - + + + LPS - + + + Inhibitor - - BD1 BD2 Inhibitor - - BD1 BD2 Inhibitor - - BD1 BD2 Inhibitor - - BD1 BD2 Egr2 enhancer Egr2 Il12b HSS1 enhancer Il12b 0.7 0.7 P = 5.7e-40 0.4 0.20 P = 4.8e-07 P = 1.0 P = 5.2e-09 Gene burst fraction Gene burst fraction Enh. burst fraction Enh. burst fraction 0.6 0.6 P = 2.5e-106 P = 2.0e-61 P = 4.4e-16 P = 3.9e-84 0.5 0.5 P = 1.6e-14 0.3 0.15 P = 1.3e-63 P = 3.7e-26 P = 7.3e-32 0.4 0.4 P = 3.4e-09 P = 4.3e-23 0.2 0.10 0.3 0.3 0.2 0.2 P = 0.0 0.1 0.05 P = 3.9e-110 0.1 0.1 0 0 0 0 LPS - + + + LPS - + + + LPS - + + + LPS - + + + Inhibitor - - BD1 BD2 Inhibitor - - BD1 BD2 Inhibitor - - BD1 BD2 Inhibitor - - BD1 BD2 b 0.7 Egr2 wild-type 0.14 Egr2 Rad21-/- Il12b wild-type Il12b Rad21-/- P = 6.5e-42 P = 1.0 0.20 0.04 P = 1.0 P = 4.8e-20 Gene burst fraction Gene burst fraction Gene burst fraction Gene burst fraction 0.6 P = 6.3e-99 0.12 P = 0.008 P = 3.3e-113 P = 2.2e-08 0.5 0.10 0.15 0.03 P = 6.2e-11 P = 4.2e-06 P = 0.15 P = 9.9e-37 P = 5.2e-12 0.08 P = 9.1e-05 0.4 0.06 0.10 0.02 0.3 0.2 P = 2.1e-274 0.04 0.05 P = 1.9e-63 0.01 0.1 0.02 0 0 0 0 LPS - + + + LPS - + + + LPS - + + + LPS - + + + Inhibitor - - BD1 BD2 Inhibitor - - BD1 BD2 Inhibitor - - BD1 BD2 Inhibitor - - BD1 BD2 Egr2 Wild-type Rad21-/- Il12b Wild-type Rad21-/- 6.5e-42 1.0 4.8e-20 1.0 2.0 2.0 1.9e-63 2.2e-08 Gene burst fraction 6.3e-99 Gene burst fraction 4.2e-06 0.0081 9.1e-05 2.1e-274 3.3e-113 (norm. to LPS) 6.2e-11 5.2e-12 (norm. to LPS) 1.5 0.15 1.5 9.9e-37 1.0 1.0 0.5 0.5 0.0 0.0 LPS - + + + - + + + LPS - + + + - + + + Inhibitor - - BD1 BD2 - - BD1 BD2 Inhibitor - - BD1 BD2 - - BD1 BD2 Fig. 9 Perturbations of cohesin and the second BET bromodomain BD2 have non-additive effects on burst fractions. a Transcriptional burst fractions of the Egr2 enhancer, Egr2, the Il12b HSS enhancer and Il12b in wild-type macrophages 1 h after LPS activation in the presence of the BD1-selective inhibitor GSK778 (1μM), the BD2-selective inhibitor GSK046 (1 μM), or carrier (DMSO). Egr2 enhancer Egr2: n = 12,032 cells, three independent biological replicates. Il12b HSS enhancer Il12b: n = 17,610 cells, 3–5 independent biological replicates. b Transcriptional burst fractions of Egr2 and Il12b in wild-type (data as in b) and Rad21−/− macrophages 1 h after LPS activation in the presence of the BD1-selective inhibitor GSK778 (1μM), the BD2-selective inhibitor GSK046 (1 μM), or carrier (DMSO). Egr2: n = 12,032 wild-type and 9097 Rad21−/− cells, three independent biological replicates. Il12b: n = 28,863 wild- type and 23,155 Rad21−/− cells, 3–5 independent biological replicates. looping appears to be neither strictly required nor sufficient for contacts in the range of 10–100 kb are often retained, while the activation of target genes40–43. For example, transcriptional contacts in the range of 100 kb and 1 Mb are mostly lost63. In bursting of Sox2 appears uncorrelated with enhancer–promoter agreement with these findings, enhancer and promoter bursts proximity41,44. Second, while active transcription of genes like occurred over a range of distances in our experimental system. Shh, Fos, and Arc is clearly associated with enhancer–promoter Within the resolution limits of our experimental system, cumu- looping10,45,46, contacts between enhancers and promoters at lative enhancer–promoter distances showed little—if any—dif- these10,45 and other genes47–49 occur robustly in the absence of ferences between wild-type and cohesin-deficient cells. cohesin across distances as far as ~100 kb. Third, promoter Selective inhibition of the BD2 bromodomain by GSK046 capture Hi-C in RAD21-depleted HeLa cells indicates that reduced gene transcription, but not enhancer activity at the Il12b 12 NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 ARTICLE locus. This phenocopies the impact of cohesin depletion and Immunofluorescence staining and confocal microscopy. Sterile 13 mm cover- implicates BD2 as well as cohesin in the functional coupling slips were seeded with 2 × 105 cells and kept in a cell culture incubator at 37 ◦C overnight. Coverslips were washed with PBS and fixed for 5 min with 300 μL ice- between enhancer activity and target gene transcription. More- cold 100% methanol at room temperature followed by three washes with PBS and over, BD2 inhibition reduced Egr2 and Il12b transcriptional burst stored at 4 °C overnight. Cells were permeabilised with 200 μL/sample 0.5% Triton fractions in wild-type macrophages but had no additive effect to X-100 in PBS for 10 min, incubated with 10% serum for 30 min, and primary the loss of cohesin in RAD21-deficient macrophages. Taken antibody incubation was performed at 4 °C overnight using RAD21 (Abcam together with the observation that inducible genes appear to be ab154769) 1:500 in 10% serum. Cells were washed three times with PBS and secondary antibody incubation was performed with goat α-rabbit 488 (Invitrogen particularly dependent on both BD2 (ref. 54) and cohesin9,10, A11034 1:500) for 1 h at room temperature. Confocal stacks were acquired in Leica these findings point to similarities between cohesin and BD2. To SP8 microscope with a voxel size of 0.144 × 0.144 × 0.99 μm, ×40 oil objective). For complement the inhibitor data presented here, future experiments the analysis, nuclei and cell outlines were identified from maximum projections by will address the impact of BD2 mutations on transcriptional burst CellProfiler v2.2.0 (ref. 70) as primary and secondary objects. To normalise variable levels of background, RAD21 expression for each cell was defined as the mean parameters. nuclear intensity of the antibody channel minus the cytoplasmic mean. Example Human CdLS is mainly associated with mutations in cohesin- images were filtered to remove speckle background noise using a CellProfiler filter. related genes50. Strikingly, a missense mutation in the BD2 bro- modomain of BRD4 causes CdLS-like developmental smRNA-FISH. Cells were seeded at a density of 2 × 105 per 13 mm coverslips or abnormalities51. Cell cycle delay, susceptibility to DNA damage, 1 × 105 per well in eight well chamber slides (Ibidi, 80826), kept in a cell culture and defective DNA repair have been proposed as shared pathways incubator at 37 ◦C overnight, and activated with LPS where indicated. In experiments using exon probes, cells were stained with CellTracker (Invitrogen C34552 or C2925, triggered by mutations in cohesin-related genes and in BD2 5 nM, 15 min, 37 °C) to facilitate the assignment of cytoplasmic smRNA-FISH signals (Ref. 55). Our analysis shows that in quiescent macrophages, loss to individual cells. Locus-specific exon, intron or enhancer smRNA-FISH probes of cohesin does not trigger cell cycle- or DNA damage-related (Supplementary Table 1) and View RNA ISH cell assay kits (QVC0001) were pur- pathways. Our results do not question the observation that chased from Affymetrix and used as advised by the manufacturer. To increase throughput, we used the QuantiGene View RNA HC Screening Assay (QVP0011, mutations in cohesin-related genes and in BD2 can lead to cell Affymetrix) protocol on chamber slides with the following modifications: First, cycle aberrations and DNA damage in rapidly proliferating volumes were doubled. PBS buffers contained 2 nM of an RNAse inhibitor (Ribo- cells55. They do, however, show that neither cell cycle aberrations nucleoside Vanadyl Complex, NEB S1402S). Cells were fixed in 4% for- nor DNA damage are required for the shared effects of cohesin maldehyde + 2.5% of Glacial Acetic Acid to detect nuclear transcripts. Chamber slides and BD2 on transcription. were stained with DAPI (Thermo Scientific, 62248) 1:3000 in PBS for 5 min at room temperature. Coverslips were mounted using ProLong™ Gold Antifade with DAPI BET proteins can impact transcription by a variety of (Thermo Scientific, P36931). Samples were stored at 4 °C and imaged within 5 days. mechanisms, which include RNA polymerase pause release and transcriptional elongation64–66, as well as modulation of eRNA smRNA-FISH image acquisition and analysis. Coverslips were imaged using a transcription and enhancer activity32,66–68. The availability of Leica TCS SP8 inverted microscope with an oil immersion 40×/1.30 NA. For each selective BD inhibitors and the identification of BD2 mutations in sample, 15 to 30 imaging fields were selected manually with a pixel size of 101.4 nm human disease51 provide a path for future studies to link specific and a z-step size of 300 nm. Image analysis used the MATLAB version of Cell- Profiler (v1.0.0) with specific image-based transcriptomics modules71,72 for nuclei, BET functions to specific protein domains. Such studies will be cytoplasm and dot identification. motivated not only by the desire to dissect the regulation of To increase throughput, chamber slides were imaged with a Leica TCS SP5 transcription but also by the relevance of BET proteins to human inverted microscope with a dry 20 × /0.75 NA with a pixel size of 189 nm and a development51,52 and cancer54,64,65. In this regard, mapping the z-step size of 500 nm. The combination of a dry objective, chamber slides and MatrixScreener software allowed imaging of up to eight samples per night. The mechanisms that functionally link BET proteins and cohesin is Leica TCS SP5 was used in resonant mode with a bidirectional scan, with a line important for understanding the selective advantage of somatic average of 8–16 to reduce noise. Imaging fields were selected automatically, with cohesin mutations in haematopoietic malignancies and solid human supervision only to ensure that cells were present. Fluorescence reflection tumours69. In summary, our results delineate cohesin-dependent was used for automatic focusing. Nuclei and cell outlines were identified from maximum projections using CellProfiler v2.2.0 modules as primary and secondary and -independent steps in enhancer-driven gene activation. objects70. A size filter was applied to avoid scoring cell clusters as individual cells. Cohesin functionally couples the burst frequencies of inducible Outlines were visually inspected, and then loaded in FISHQuant v3a (ref. 73). enhancers and promoters by mechanisms beyond spatial FISHQuant counted mRNA molecules by fitting 3D Gaussian in the 3D stacks. proximity. Identification, quantification, and locus assignment of transcriptional bursts. During a transcriptional burst, active TSSs accumulate signals that are brighter than individual transcripts. Burst intensities were estimated as the ratio between the Methods integrated intensity of pixels at the TSS and the same area around averaged individual Mice and cell culture. Mouse work was performed according to the Animals transcripts. Enhancers and TSSs were assigned to alleles based on distance (1 μm) (Scientific Procedures) Act Mouse work was done under a project licence issued by using a custom Julia script (Julia v1.6, https://github.com/IreneRobles/TSSs). The the UK Home Office, UK following review by the Imperial College London Animal expected co-occurrence of transcriptional bursts was calculated as the probability of Welfare and Ethical Review Body (AWERB). Bone marrow cells from Rosa26- one multiplied by the other. TSS-TSS distances were determined based on xy coor- ERt2Cre Rad21WT/WT or Rad21lox/lox mice9 on a mixed C57BL/6 129 background dinates. Chromatic aberration was measured at centroid distances using 0.1 µm were cultured in complete DMEM medium (10% FCS, 1% Penicillin-Streptomycin, TetraSpeck™ microspheres and was 90.8 nm between 488 and 633 nm (Ifnb1 - Ifnb1 0.05 mM β-mercaptoethanol, 2 mM L-glutamine, 1 mM Na Pyruvate), 20% L929- L2 enhancer), 102.3 nm between 561 and 488 nm (Il12b - Il12b HSS1 enhancer and conditioned media. Cre was induced on day 4 by 200 nM 4-hydroxy tamoxifen Peli1-Peli1 enhancer), and 41.4 nm between 561 and 633 nm (Prdm1 - Prdm1 (Sigma-Aldrich H7904). Macrophages were stimulated on day 7 with 10 ng/ml of enhancer). In each case, the chromatic aberration was smaller than the xy pixel size of LPS from Salmonella typhosa (Sigma-Aldrich L7895). Where indicated, 1 μm 189 nm, and the data was not corrected for chromatic aberration. GSK778 or GSK046 or carrier (DMSO) were added at the same time as LPS. Burst size and burst frequency based on the moments of mRNA distribution. To infer mean μ and variance σ2, 50 cells were sampled from each distribution 1000 times. Burst size was calculated as bm = σ2/μ and burst frequency as fm = μ/(bm−1). Immunoblots. Cell lysates were separated by SDS-PAGE on 7.5% gels (100 V, 90 min) and transferred to nitrocellulose membranes (90 V, 90 min). Membranes were blocked with 5% milk in Tris-buffered saline (TBS) with 0.1% Tween-20 and RNA-seq analysis. 100 bp paired-end RNASeq reads were aligned to mouse incubated with primary antibodies against RAD21 (Abcam 154769, 1/1000) and genome mm9 using Tophat2 (ref. 74) with arguments '–library-type fr-first strand Tubulin (Sigma T9026, 1/1000) in TBS 3% milk overnight at 4 °C, washed 3x in -b2-very-sensitive -b2-L 25' with gene annotation from Ensembl version 67. Read TBS, incubated with Alexa fluor 680 goat anti-rabbit (Life Technologies A21109 1/ counts on genes were summarised using HTSeq-count (ref. 75). The RUVseq 10.000) to detect RAD21 or Alexa fluor 680 goat anti-mouse (Life Technologies Bioconductor package was applied with RUV k = 3 (1.18.0; ref. 76) for ERCC spike- A21057, 1/10.000) to detect Tubulin for 1 h at 4 °C, washed 3x in TBS, imaged (LI- ins and replicate samples as controls. Differentially expressed genes were identified COR Biosciences) and analysed (Image Studio Lite, LI-COR Biosciences). by DESeq2 (ref. 77). NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications 13 ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 GRO-seq analysis. GRO-seq libraries were sequenced as 50 bp single-end reads in 4. Dekker, J. & Mirny, L. The 3D genome as moderator of chromosomal two biological replicates. The ten most 3′ bases were discarded based on fastqc communication. Cell 164, 1110–1121 (2016). quality assessment. Reads were aligned to mouse genome mm9 using bowtie with 5. Merkenschlager, M. & Nora, E. P. CTCF and cohesin in genome folding and arguments '-l 30 -m 10 -n 2 –trim3 10'. Read counts on enhancers were computed transcriptional gene regulation. Annu. Rev. Genomics Hum. Genet. 17, 17–43 using the summarizeOverlaps function from GeomicAlignments R Package. Dif- (2016). ferentially transcribed enhancers were identified using DESeq2. 6. McCord, R., Kaplan, N. & Giorgetti, L. Chromosome conformation capture and beyond: toward an integrative view of chromosome structure and scRNA-seq. WT and Rad21-deficient macrophages were sorted by flow cytometry function. Mol. Cell 77, 688–708 (2020). into 386 well plates, and MARS-Seq libraries78 were prepared and sequenced using 7. Nora, E. P. et al. Targeted degradation of CTCF decouples local insulation of two-level indexing of 192 cells per pool and multiple pools per sequencing lane. chromosome domains from genomic compartmentalization. Cell 169, Sequencing was carried out as paired-end reads, wherein the first read contains the 930–944.e22 (2017). transcript sequence and the second read the cell barcode and UMIs. A quality check 8. Luppino, J. M. et al. Cohesin promotes stochastic domain intermingling to of the generated reads was performed with the FastQC quality control suite. Samples ensure proper regulation of boundary-proximal genes. Nat. Genet. 52, that reached the quality standards were then processed to deconvolute the reads to 840–848 (2020). single-cell level by de-multiplexing according to the cell and pool barcodes. Reads 9. Cuartero, S. et al. Control of inducible gene expression links cohesin to were filtered to remove polyT sequences. Sequencing reads from human, mouse, or hematopoietic progenitor self-renewal and differentiation. Nat. Immunol. 19, canine cells were mapped with the RNA pipeline of the GEMTools 1.7.0 suite79 using 932–941 (2018). default parameters (6% of mismatches, minimum of 80% matched bases, and mini- 10. Calderon, L. et al. Reliance of neuronal gene expression on cohesin scales with mum quality threshold of 26) and the genome references for human (Gencode release chromatin loop length. eLife 11, e76539 (2022). 24, assembly GRCh38.p5), mouse (Gencode release M8, assembly GRCm38.p4), and 11. Oh, S. et al. Enhancer release and retargeting activates disease-susceptibility dog (Ensembl v84, assembly CanFam3.1). The analysis of spike-in control RNA genes. Nature 595, 735–740 (2021). content allowed us to identify empty wells and barcodes with more than 15% of reads 12. Glass, C. K. & Natoli, G. Molecular control of activation and priming in mapping to spike-in artificial transcripts were discarded. Cells with less than 60% of macrophages. Nat. Immunol. 17, 26–33 (2016). reads mapping to the reference genome or more than 2 × 106 total reads were dis- 13. Filippakopoulos, P. et al. Histone recognition and large-scale structural carded. Gene quantification was performed using UMI-corrected transcript infor- analysis of the human bromodomain family. Cell 149, 214–231 (2012). mation to correct for amplification biases, collapsing read counts for reads mapping 14. Larke, M. S. C. et al. Enhancers predominantly regulate gene expression during on a gene with the same UMI (allowing an edit distance up to two nucleotides in UMI differentiation via transcription initiation. Mol. Cell 81, 983–997.e7 (2021). comparisons). Only unambiguously mapped reads were considered. Genes not 15. Larsson, A. J. M. et al. Genomic encoding of transcriptional burst kinetics. expressed in at least 5% of the cells were discarded. Nature 565, 251–254 (2019). 16. Golding, I., Paulsson, J., Zawilski, S. M. & Cox, E. C. Real-time kinetics of gene Statistical analysis. Pearson correlation tests were done with the cor.test function in activity in individual bacteria. Cell 123, 1025–1036 (2005). R. To test significant changes in event frequency across replicates we used a 17. Chubb, J. R., Trcek, T., Shenoy, S. M. & Singer, R. H. Transcriptional pulsing Cochran–Mantel–Haenszel test with replicate as strata variable in R. Bonferroni of a developmental gene. Curr. Biol. 16, 1018–1025 (2006). correction for multiple testing was implemented in MutipleTesting.jl Julia package. To 18. Raj, A., Peskin, C. S., Tranchina, D., Vargas, D. Y. & Tyagi, S. Stochastic test significant changes in distribution values across replicates, we used one-way mRNA synthesis in mammalian cells. PLoS Biol. 4, e309 (2006). ANOVA with Tukey HSD correction in R. To test whether the relationship between 19. Suter, D. M., Molina, N., Naef, F. & Schibler, U. Origins and consequences of two variables differs between conditions we implemented two-way ANOVA in R. transcriptional discontinuity. Curr. Opin. Cell Biol. 23, 657–662 (2011). 20. Falo-Sanjuan, J., Lammers, N. C., Garcia, H. G. & Bray, S. J. Enhancer priming Reporting summary. Further information on research design is available in the Nature enables fast and sustained transcriptional responses to Notch signaling. Dev. Research Reporting Summary linked to this article. Cell 50, 411–425.e8 (2019). 21. Lee, C., Shin, H. & Kimble, J. Dynamics of Notch-dependent transcriptional bursting in its native context. Dev. Cell 50, 426–435.e4 (2019). Data availability 22. Stavreva, D. A. et al. Transcriptional bursting and co-bursting regulation by scRNA-seq data generated for this study have been deposited at the Gene Expression steroid hormone release pattern and transcription factor mobility. Mol. Cell Omnibus (GEO) under accession code GSE190622. smRNA-FISH images generated for this 75, 1161–1177.e11 (2019). study have been deposited at the EMBL-EBI BioImage Archive under accession code 23. Bagnall, J. et al. Gene-specific linear trends constrain transcriptional variability S-BIAD338. The GRO-seq, RNA-seq, H3K27ac ChIP-seq data9 used in this study are of the toll-like receptor signaling. Cell Syst. 11, 300–14.e8 (2020). available in the GEO database under accession code GSE108599. Macrophage Hi-C data 24. Bass, V. L., Wong, V. C., Bullock, M. E., Gaudet, S. & Miller-Jensen, K. TNF used in this study are available in the GEO database under accession code GSE115524. PU.1 stimulation primarily modulates transcriptional burst size of NF-kappaB- ChIP-seq data80 used in this study are available in the GEO database under accession code regulated genes. Mol. Syst. Biol. 17, e10127 (2021). GSE56121. STAT2 ChIP-seq data80 used in this study are available in the GEO database 25. Senecal, A. et al. Transcription factors modulate c-Fos transcriptional bursts. under accession code GSE56123. IRF3 ChIP-seq data81 used in this study are available in the Cell Rep. 8, 75–83 (2014). GEO database under accession code GSE67343. Source data are provided with this paper. 26. Hendy, O., Campbell, J. Jr., Weissman, J. D., Larson, D. R. & Singer, D. S. Differential context-specific impact of individual core promoter elements on Code availability transcriptional dynamics. Mol. Biol. Cell 28, 3360–3370 (2017). Code used in this study is available on Github, including code to make paper figures82 27. Dar, R. D. et al. Transcriptional burst frequency and burst size are equally https://github.com/IreneRobles/Code_Paper, TSS quantification and analysis from modulated across the human genome. Proc. Natl Acad. Sci. USA 109, FISHquant outlines83https://github.com/IreneRobles/TSSs, processing of scRNA-seq 17454–17459 (2012). data84 https://github.com/IreneRobles/SingleCellExperiment, dimensionality reduction 28. Chen, L. F. et al. Enhancer histone acetylation modulates transcriptional on matrix data85 https://github.com/IreneRobles/DimensionalityReduction and scRNA- bursting dynamics of neuronal activity-inducible genes. Cell Rep. 26, seq data86 https://github.com/IreneRobles/SingleCellDimensionalityReduction, and to 1174–1188.e5 (2019). generate input files for FISHquant and CellProfiler from microscope lif files87 https:// 29. So, L. H. et al. General properties of transcriptional time series in Escherichia github.com/IreneRobles/Make2Dand3Dstacks_fromlif. coli. Nat. Genet. 43, 554–560 (2011). 30. Bhatt, D. M. et al. Transcript dynamics of proinflammatory genes revealed by sequence analysis of subcellular RNA fractions. Cell 150, 279–290 (2012). Received: 17 December 2021; Accepted: 6 June 2022; 31. Ostuni, R. et al. Latent enhancers activated by stimulation in differentiated cells. Cell 152, 157–171 (2013). 32. Kaikkonen, M. U. et al. Remodeling of the enhancer landscape during macrophage activation is coupled to enhancer transcription. Mol. Cell 51, 310–325 (2013). 33. Zhou, L. et al. An inducible enhancer required for Il12b promoter activity in References an insulated chromatin environment. Mol. Cell Biol. 27, 2698–2712 1. Fudenberg, G. et al. Formation of chromosomal domains by loop extrusion. (2007). Cell Rep. 15, 2038–2049 (2016). 34. Banerjee, A. R., Kim, Y. J. & Kim, T. H. A novel virus-inducible enhancer of 2. Rao, S. S. P. et al. Cohesin loss eliminates all loop domains. Cell 171, the interferon-β gene with tightly linked promoter and enhancer activities. 305–320.e24 (2017). Nucleic Acids Res. 42, 12537–12554 (2014). 3. Schwarzer, W. et al. Two independent modes of chromatin organization 35. De Santa, F. et al. A large fraction of extragenic RNA pol II transcription sites revealed by cohesin removal. Nature 551, 51–56 (2017). overlap enhancers. PLoS Biol. 8, e1000384 (2010). 14 NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 ARTICLE 36. Kim, T. K. et al. Widespread transcription at neuronal activity-regulated 67. Kanno, T. et al. BRD4 assists elongation of both coding and enhancer RNAs enhancers. Nature 465, 182–187 (2010). by interacting with acetylated histones. Nat. Struct. Mol. Biol. 21, 1047–1057 37. Andersson, R. et al. Nuclear stability and transcriptional directionality (2014). separate functionally distinct RNA species. Nat. Commun. 5, 5336 (2014). 68. Rahnamoun, H. et al. RNAs interact with BRD4 to promote enhanced 38. Shalek, A. K. et al. Single-cell RNA-seq reveals dynamic paracrine control of chromatin engagement and transcription activation. Nat. Struct. Mol. Biol. 25, cellular variation. Nature 510, 363–369 (2014). 687–697 (2018). 39. Nasmyth, K. & Haering, C. H. Cohesin: its roles and mechanisms. Annu Rev. 69. De Koninck, M. & Losada, A. Cohesin mutations in cancer. Cold Spring Harb. Genet. 43, 525–558 (2009). Perspect. Med. 6, a026476 (2016). 40. Chen, H. et al. Dynamic interplay between enhancer–promoter topology and 70. Carpenter, A. E. et al. CellProfiler: image analysis software for identifying and gene activity. Nat. Genet. 50, 1296–1303 (2018). quantifying cell phenotypes. Genome Biol. 7, R100 (2006). 41. Alexander, J. M. et al. Live-cell imaging reveals enhancer-dependent Sox2 71. Battich, N., Stoeger, T. & Pelkmans, L. Image-based transcriptomics in transcription in the absence of enhancer proximity. eLife 8, e41769 (2019). thousands of single human cells at single-molecule resolution. Nat. Methods 42. Benabdallah, N. S. et al. Decreased enhancer-promoter proximity 10, 1127–1133 (2013). accompanying enhancer activation. Mol. Cell 76, 473–484.e7 (2019). 72. Stoeger, T., Battich, N., Herrmann, M. D., Yakimovich, Y. & Pelkmans, L. 43. Heist, T., Fukaya, T. & Levine, M. Large distances separate coregulated genes Computer vision for image-based transcriptomics. Methods 85, 44–53 (2015). in living Drosophila embryos. Proc. Natl Acad. Sci. USA 116, 15062–15067 73. Mueller, F. et al. FISH-quant: automatic counting of transcripts in 3D FISH (2019). images. Nat. Methods 10, 277–278 (2013). 44. Huang, H. et al. CTCF mediates dosage- and sequence-context-dependent 74. Kim, D. et al. TopHat2: accurate alignment of transcriptomes in the presence transcriptional insulation by forming local chromatin domains. Nat. Genet. of insertions, deletions and gene fusions. Genome Biol. 14, R36 (2013). 53, 1064–1074 (2021). 75. Anders, S., Pyl, P. T. & Huber, W. HTSeq-a Python framework to work with 45. Kane, L. et al. Cohesin is required for long-range enhancer action. Preprint at high-throughput sequencing data. Bioinformatics 31, 166–169 (2015). biorxiv https://doi.org/10.1101/2021.06.24.449812 (2021). 76. Risso, D., Ngai, J., Speed, T. P. & Dudoit, S. Normalization of RNA-seq data 46. Beagan, J. A. et al. Three-dimensional genome restructuring across timescales using factor analysis of control genes or samples. Nat. Biotechnol. 32, 896–902 of activity-induced neuronal gene expression. Nat. Neurosci. 23, 707–717 (2014). (2020). 77. Love, M. I., Huber, W. & Anders, S. Moderated estimation of fold change and 47. Lavagnolli, T. et al. Initiation and maintenance of pluripotency gene dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550 (2014). expression in the absence of cohesin. Genes Dev. 29, 23–38 (2015). 78. Jaitin, D. A. et al. Massively parallel single-cell RNA-seq for marker-free 48. Rinzema, N. J. et al. Building regulatory landscapes: enhancer recruits cohesin decomposition of tissues into cell types. Science 343, 776–779 (2014). to create contact domains, engage CTCF sites and activate distant genes. Nat. 79. Marco-Sola, S., Sammeth, M., Guigó, R. & Ribeca, P. The GEM mapper: fast, Struct. Mol. Biol. 29, 563–574 (2022). accurate and versatile alignment by filtration. Nat. Methods 9, 1185–1188 49. Hsieh, T-H. S. et al. Enhancer-promoter interactions and transcription are (2012). maintained upon acute loss of CTCF, cohesin, WAPL, and YY1. Preprint at 80. Mancino, A. et al. A dual cis-regulatory code links IRF8 to constitutive and biorxiv https://doi.org/10.1101/2021.07.14.452365 (2021). inducible gene expression in macrophages. Genes Dev. 29, 394–408 (2015). 50. Deardorff, M. A., Noon, S. E. & Krantz, I. D. In GeneReviews® (eds. Adam, M. 81. Tong, A. J. et al. A stringent systems approach uncovers gene-specific P. et al.) (University of Washington, 2018). mechanisms regulating inflammation. Cell 165, 165–179 (2016). 51. Olley, G. et al. BRD4 interacts with NIPBL and BRD4 is mutated in a Cornelia 82. Robles-Rebollo, I. Cohesin couples transcriptional bursting probabilities of de Lange–like syndrome. Nat. Genet. https://doi.org/10.1038/s41588-018- inducible enhancers and promoters. Github:IreneRobles/Code_Paper. https:// 0042-y (2018). doi.org/10.5281/zenodo.6489770 (2022). 52. Alesi, V. et al. Confirmation of BRD4 haploinsufficiency role in Cornelia de 83. Robles-Rebollo, I. Cohesin couples transcriptional bursting probabilities of Lange-like phenotype and delineation of a 19p13.12p13.11 gene contiguous inducible enhancers and promoters. Github:IreneRobles/TSSs. https://doi.org/ syndrome. Ann. Hum. Genet. 83, 100–109 (2019). 10.5281/zenodo.6490112 (2022). 53. Crump, N. T. et al. BET inhibition disrupts transcription but retains 84. Robles-Rebollo, I. Cohesin couples transcriptional bursting probabilities of enhancer-promoter contact. Nat. Commun. 12, 223 (2021). inducible enhancers and promoters. Github:IreneRobles/SingleCellExperiment. 54. Gilan, et al. Selective targeting of BD1 and BD2 of the BET proteins in cancer https://doi.org/10.5281/zenodo.6490124 (2022). and immuno-inflammation. Science https://doi.org/10.1126/science.aaz8455 85. Robles-Rebollo, I. Cohesin couples transcriptional bursting probabilities of (2020) inducible enhancers and promoters. Github:IreneRobles/Dimensionality 55. Olley, G. et al. Cornelia de Lange syndrome-associated mutations cause a Reduction. https://doi.org/10.5281/zenodo.6490176 (2022). DNA damage signalling and repair defect. Nat. Commun. 12, 3127 (2021). 86. Robles-Rebollo, I. Cohesin couples transcriptional bursting probabilities 56. Liberzon, A. et al. The molecular signatures database (MSigDB) hallmark gene of inducible enhancers and promoters. Github:IreneRobles/Single set collection. Cell Syst. 1, 417–425 (2015). CellDimensionalityReduction. https://doi.org/10.5281/zenodo.6490213 57. Walters, M. C. et al. Enhancers increase the probability but not the level of (2022). gene expression. Proc. Natl Acad. Sci. USA 92, 7125–7129 (1995). 87. Robles-Rebollo, I. Cohesin couples transcriptional bursting probabilities of 58. Bartman, C. R., Hsu, S. C., Hsiung, C. C.-S., Raj, A. & Blobel, G. A. Enhancer inducible enhancers and promoters. Github:IreneRobles/Make2Dand3Dstacks_ regulation of transcriptional bursting parameters revealed by forced fromlif. https://doi.org/10.5281/zenodo.6490130 (2022). chromatin looping. Mol. Cell 62, 237–247 (2016). 59. Fukaya, T., Lim, B. & Levine, M. Enhancer control of transcriptional bursting. Cell 166, 358–368 (2016). 60. Mermet, J. et al. Clock-dependent chromatin topology modulates circadian Acknowledgements transcription and behavior. Genes Dev. 32, 347–358 (2018). We thank Dr. Thomas A Milne (MRC Molecular Haematology Unit, Oxford) for the 61. Deng, W. et al. Controlling long-range genomic interactions at a native locus discussion and Dr. Vahid S Shahrezaei and Ilinca Patrascan for critical reading of the by targeted tethering of a looping factor. Cell 149, 1233–1244 (2012). manuscript. This work was supported by the Medical Research Council UK and The 62. Morgan, S. L. et al. Manipulation of nuclear architecture through CRISPR- Wellcome Trust (Investigator Award 099276/Z/12/Z to M.M.). mediated chromosomal looping. Nat. Commun. 8, 15993 (2017). 63. Thiecke, M. J. et al. Cohesin-dependent and -independent mechanisms mediate chromosomal contacts between promoters and enhancers. Cell Rep. Author contributions 32, 107929 (2020). I.R.-R., A.C.-S., S.W. and A.G.C. did experiments, I.R.-R., A.C.-S., E.M., S.W. and A.G.C. 64. Jang, M. K. et al. The bromodomain protein Brd4 is a positive regulatory analysed data, I.R.-R. curated data, I.R.-R., M.M.K., C.W. and D.D. wrote scripts, H.H. component of P-TEFb and stimulates RNA polymerase II-dependent devised the scRNA-seq pipeline, C.W. and D.D. devised the microscopy pipeline, I.R.-R., transcription. Mol. Cell 19, 523–534 (2005). S.C., S.M., I.R., R.K.P., M.P.H.S., A.G.F. and M.M. conceptualised the study, IR-R and 65. Dawson, M. A. et al. Inhibition of BET recruitment to chromatin as an MM made figures and wrote the manuscript. effective treatment for MLL-fusion leukaemia. Nature 478, 529–533 (2011). 66. Winter, G. E. et al. BET bromodomain proteins function as master transcription elongation factors independent of CDK9 recruitment. Mol. Cell Competing interests 67, 5–18.e19 (2017). The authors declare no competing interests. NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications 15 ARTICLE NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-022-31192-9 Additional information Open Access This article is licensed under a Creative Commons Supplementary information The online version contains supplementary material Attribution 4.0 International License, which permits use, sharing, available at https://doi.org/10.1038/s41467-022-31192-9. adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Correspondence and requests for materials should be addressed to Matthias Commons license, and indicate if changes were made. The images or other third party Merkenschlager. material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the Peer review information Nature Communications thanks Pradeepa Madapura, James article’s Creative Commons license and your intended use is not permitted by statutory Davies, Pawel Paszek and the other, anonymous, reviewer(s) for their contribution to the regulation or exceeds the permitted use, you will need to obtain permission directly from peer review of this work. Peer reviewer reports are available the copyright holder. To view a copy of this license, visit http://creativecommons.org/ licenses/by/4.0/. Reprints and permission information is available at http://www.nature.com/reprints Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in © The Author(s) 2022 published maps and institutional affiliations. 16 NATURE COMMUNICATIONS | (2022)13:4342 | https://doi.org/10.1038/s41467-022-31192-9 | www.nature.com/naturecommunications