Developmental and housekeeping transcriptional programs display distinct modes of enhancer-enhancer cooperativity in Drosophila

Abstract

Genomic enhancers are key transcriptional regulators which, upon the binding of sequence-specific transcription factors, activate their cognate target promoters. Although enhancers have been extensively studied in isolation, a substantial number of genes have more than one simultaneously active enhancer, and it remains unclear how these cooperate to regulate transcription. Using Drosophila melanogaster S2 cells as a model, we assay the activities of more than a thousand individual enhancers and about a million enhancer pairs toward housekeeping and developmental core promoters with STARR-seq. We report that housekeeping and developmental enhancers show distinct modes of enhancer-enhancer cooperativity: while housekeeping enhancers are additive such that their combined activity mirrors the sum of their individual activities, developmental enhancers are super-additive and combine multiplicatively. Super-additivity between developmental enhancers is promiscuous and neither depends on the enhancers' endogenous genomic contexts nor on specific transcription factor motif signatures. However, it can be further boosted by Twist and Trl motifs and saturates for the highest levels of enhancer activity. These results have important implications for our understanding of gene regulation in complex multi-enhancer developmental loci and genomically clustered housekeeping genes, providing a rationale to interpret the transcriptional impact of non-coding mutations at different loci.

Fig. 1. High-throughput assessment of the individual and combined activities of many enhancers.

a Overview of the STARR-seq reporter assay used to simultaneously measure the individual and combined activities of many enhancers. Random control (in gray) and candidate sequences (colors) are fused to the 5’ and the 3’ ends of a transcriptionally inert spacer and cloned downstream of a core promoter, whose transcription mirrors enhancers’ individual and combined activities. b Correlation between 3’ (x-axis) and 5’ (y-axis) individual activities of 953 candidate sequences. The dotted line represents the identity line (y = x), and Pearson’s correlation coefficient is shown on the top left (r). The color code displays sequences’ activity rank inferred from a previously published STARR-seq dataset. c Quantification of the activity of pairs consisting of two random control sequences (Ctl./Ctl., in gray), one control sequence paired with a candidate sequence either in the 5’ (Enh./Ctl., in blue) or the 3’ (Ctl./Enh., in purple) location, or two enhancer sequences (Enh./Enh., in green). n = 16,032; 90,559; 91,817 and 51,7071 pairs, respectively. Two-sided Wilcoxon test P-values are shown; ****P < 2.2e-308. Box plots show the median (line), upper and lower quartiles (box) ±1.5× interquartile range (whiskers), outliers are not shown. d Heatmap of paired activities (see color legend) ranked by individual activities of the 3’ (x-axis) and 5’ (y-axis) candidate sequences. 3’ and 5’ activities are depicted as bar charts on the x and y axes, respectively, with active sequences being highlighted with dashed lines (log2 individual activity >1). e Correlation between candidate sequence pairs (A/B, x-axis) and the reciprocal combinations (B/A, y-axis).

Fig. 2. Developmental enhancers combine multiplicatively.

a Schematic illustration of an additive model (where the number of RNA molecules produced by each enhancer add up) and a multiplicative model (where the number of RNAs combine multiplicatively). b Scatterplots showing predicted activities (x-axis) based on an additive (left) or a multiplicative model (right) versus observed activities (y-axis), with corresponding R-squared (R²) values. Enhancer pairs in which both candidate sequences are active are highlighted using density lines (in orange), and dotted lines correspond to identity lines (y = x). c Fraction of enhancer pairs (in which both candidate sequences are active) for which the additive (in white) or the multiplicative (in blue) predicted values were the most accurate. d Fitted multiplicative model with interaction term using the 5’ and 3’ individual activities to predict the activities of enhancer pairs, with the corresponding adjusted R² value (top left). Enhancer pairs in which both candidate sequences are active are highlighted using density lines (in orange). Dotted lines correspond to the identity line (y = x). Fitted coefficients and resulting equation are shown at the bottom. e Box plots showing, for all enhancer pairs, predicted values using the three different models and observed values (x-axis). Lines connect the values for a representative set of pairs, spanning the dynamic range of observed activities. n = 517,071 pairs per box plot; box plots show the median (line), upper and lower quartiles (box) ±1.5× interquartile range (whiskers), outliers are not shown. f, g Predicted additive (x-axis) versus observed combined activities for ecdysone-inducible enhancers and OSC-specific enhancers (in red) in ecdysone-treated S2 cells (f) or OSC cells (g). As a reference, a subset of pairs containing enhancers that were also active in S2 cells are shown (see color legend). Dotted lines depict the identity lines (where observed combined activities equal expected additive outcomes), and solid lines represent fitted multiplicative models.

Fig. 3. Developmental enhancer super-additivity is promiscuous.

a Distribution of the residuals of the fitted multiplicative model with interaction term, ranging from negative (weaker super-additivity, in blue) to positive (stronger super-additivity, in orange). b For each pairwise combination of TF motifs (120 × 120 = 14,400), the mean impact on the activity (x-axis) and on the residuals (y-axis) of corresponding enhancer pairs are shown. For example, Trl/Twist pairs (in which the 5’ and 3’ enhancers contain at least one instance of the Trl and Twist motifs, respectively) show globally increased residuals (in orange). AP-1 and GATA motifs strongly affected enhancer activity, while the Dref motif negatively impacted activity and super-additivity. c, d Impact of mutating (c, in pink) or adding (d, in pink) Trl (left) or Twist motifs (right) on the predicted additive (x-axis) versus observed combined activities (y-axis) of developmental enhancer pairs. As a reference, corresponding wild-type (WT) developmental enhancer pairs are shown in gray. Dotted lines depict the identity line (where observed combined activities equal expected additive outcomes), and solid lines represent the fitted multiplicative model. For each condition, the residuals of WT versus mutant pairs were quantified (see box plots on the right) using either the additive or the fitted multiplicative model (Fit. mult, see x-axis) and compared using paired, two-sided Wilcoxon tests. Box plots show the median (line), upper and lower quartiles (box) ±1.5× interquartile range (whiskers), outliers are not shown. e Quantitative comparison of the activity of enhancer pairs from the same locus (≤20 kb distance in situ, in pink) versus an activity-matched set of distant enhancer pairs (>20 kb, in gray). n = 567 pairs per box plot. Two-sided Wilcoxon test P-values are shown; box plots show the median (line), upper and lower quartiles (box) ±1.5× interquartile range (whiskers), outliers are not shown.

Fig. 4. Housekeeping enhancers are additive and developmental enhancers are super-additive independently of the CP type.

a Scatterplots showing predicted activities (x-axis) based on an additive (left) or a multiplicative model (right) versus observed activities (y-axis) using the RpS12 housekeeping CP. Corresponding R-squared (R²) values are shown (top left), and enhancer pairs in which both candidate sequences are active are highlighted using density lines (in orange). Identity lines are shown using dotted lines (x = y). b Fraction of enhancer pairs (in which both candidate sequences are active) for which the additive (in white) or the multiplicative (in blue) predicted values were the most accurate. c Selected housekeeping (left) and developmental (right) enhancer pairs with comparable 5’ and 3’ individual activities, either with a housekeeping (hkCP, in red) or a developmental (dCP, in green) Core Promoter. For each pair, individual and combined measured activities are shown (solid gray bars) and compared to predicted activities (striped bars) using either the additive (Pred. add.) or the fitted multiplicative (Pred. fit. mult.) model. Bar heights correspond to the mean activity values and whiskers to the standard deviations. d Expected additive and observed activities of housekeeping versus developmental enhancer pairs (x-axis) using either a housekeeping (hkCP, in red) or a developmental (dCP, in green) CP. n = 3590, 2392, 3513 and 2328 pairs, respectively. Two-sided Wilcoxon test P-values are shown; box plots show the median (line), upper and lower quartiles (box) ±1.5× interquartile range (whiskers), outliers are not shown.