Transcriptional regulatory logic of the diurnal cycle in the mouse liver

Many organisms exhibit temporal rhythms in gene expression that propel diurnal cycles in physiology. In the liver of mammals, these rhythms are controlled by transcription–translation feedback loops of the core circadian clock and by feeding–fasting cycles. To better understand the regulatory interplay between the circadian clock and feeding rhythms, we mapped DNase I hypersensitive sites (DHSs) in the mouse liver during a diurnal cycle. The intensity of DNase I cleavages cycled at a substantial fraction of all DHSs, suggesting that DHSs harbor regulatory elements that control rhythmic transcription. Using chromatin immunoprecipitation followed by DNA sequencing (ChIP-seq), we found that hypersensitivity cycled in phase with RNA polymerase II (Pol II) loading and H3K27ac histone marks. We then combined the DHSs with temporal Pol II profiles in wild-type (WT) and Bmal1-/- livers to computationally identify transcription factors through which the core clock and feeding–fasting cycles control diurnal rhythms in transcription. While a similar number of mRNAs accumulated rhythmically in Bmal1-/- compared to WT livers, the amplitudes in Bmal1-/- were generally lower. The residual rhythms in Bmal1-/- reflected transcriptional regulators mediating feeding–fasting responses as well as responses to rhythmic systemic signals. Finally, the analysis of DNase I cuts at nucleotide resolution showed dynamically changing footprints consistent with dynamic binding of CLOCK:BMAL1 complexes. Structural modeling suggested that these footprints are driven by a transient heterotetramer binding configuration at peak activity. Together, our temporal DNase I mappings allowed us to decipher the global regulation of diurnal transcription rhythms in the mouse liver.


26
feedback loops of the core circadian clock and by feeding-fasting cycles. To better understand the 27 regulatory interplay between the circadian clock and feeding rhythms, we mapped DNase I 28 hypersensitive sites (DHSs) in mouse liver during a diurnal cycle. The intensity of DNase I 29 cleavages cycled at a substantial fraction of all DHSs, suggesting that DHSs harbor regulatory 30 elements that control rhythmic transcription. Using ChIP-seq, we found that hypersensitivity 31 cycled in phase with RNA polymerase II (Pol II) loading and H3K27ac histone marks. We then 32 combined the DHSs with temporal Pol II profiles in wild-type (WT) and Bmal1 -/livers to 33 computationally identify transcription factors through which the core clock and feeding-fasting Circadian clocks provide mammals with cell-autonomous and organ-based metronomes that relay 45 diurnal environmental cues to temporal gene expression programs [1,2]. In particular, diurnal 46 rhythms in mRNA transcription result from the combined actions of the autonomous circadian 47 oscillator, systemic signals and other temporal cues such as feeding-fasting cycles [3][4][5][6]. While it 48 is commonly assumed that around 10% of genes exhibit cyclic mRNA levels in the liver [7], this 49 number increases to nearly 50% when only considering liver-specifically expressed genes [8].

50
Moreover, these mRNA rhythms cover a continuum of peak times [9,10]. Although mRNAs can 51 also rhythmically accumulate due to post-transcriptional regulation [6,[11][12][13][14], it is of interest to 52 obtain a more comprehensive view on transcriptional regulators and mechanisms underlying 53 time-specific diurnal transcription. In a light-dark (LD) cycle, two main waves of transcription 54 are found, one during the day (at around ZT10) and the other towards the end of the night (around 55 ZT20), accompanied by dynamic chromatin state modifications [6,11,12].

98
To identify DNA regulatory elements controlling diurnal transcriptional rhythms in the mouse 99 liver, we mapped DHSs every 4 hours during a full light dark (LD) cycle. Specifically, C57BL/6 100 male mice were kept in standard 12h light, 12h dark cycles, and four animals were sacrificed 101 every 4 hours for one day followed by liver dissection (Methods). DNase I hypersensitivity 102 libraries were produced, sequenced, and mapped to the mouse genome using standard methods 103 (Methods). To monitor transcription activity in the same conditions, we generated ChIP-seq 104 samples for the histone modification H3K27ac (marking active regulatory elements [33]) and re-105 sequenced previous total Pol II ChiP-seq libraries [12] at increased coverage (Methods and Table   106 S1). Circadian clock outputs result in the rhythmic transcriptional activation of hundreds of genes, 107 notably through binding of CLOCK:BMAL1 heterodimers [16,17]. To validate our assays, we 108 therefore examined the known circadian output gene Dbp (Movie S1), maximally transcribed at 109 ZT8 [30], to determine whether cutting frequency at DHSs exhibited diurnal variation. We 110 detected several DHSs in the vicinity of Dbp, with high intensity and narrow signals surrounded 111 by low noise levels ( Fig 1A). Overall, enriched sites coincided well with regions identified using 112 classical DHS mapping [30], and overlapped with BMAL1 ChIP-seq regions [17] ( Fig S1). As 113 exemplified by a DHS nearby the transcription start site (TSS) of Dbp, we observed that DHSs 114 were located in regions with lower H3K27ac signals in between H3K27ac-enriched islands, 115 suggestive of TF-induced nucleosome displacement [34-36] (Fig 1B). The DNase I 116 hypersensitivity changed diurnally, notably at the TSS (Fig 1C) where the oscillations in DNase I 117 hypersensitivity, Pol II, and H3K27ac peaked in sync at ZT10 (Fig 1B). Moreover, all DHSs 118 within 15 kb of the Dbp TSS displayed oscillations with the same phase as the TSS (Fig 1D), 119 suggesting regulatory relationships between these regions and gene transcription. 120 121 We next analyzed the Npas2 gene (Movie S2), another known clock target [37]. Npas2 is a target 122 of RORs and peaks in the late night-time around ZT22 [38]. We detected several DHSs along the 123 transcribed region of Npas2 (Fig 1E), including proximal (defined as 1-10kb from a TSS) and 124 distal (defined as >10kb from a TSS) sites. The distal sites displayed high amplitude oscillations 125 of DNase I signals and H3K27ac (Fig 1F). Normalized signals at the Npas2 TSS also peaked at 126 the expected phase with maximal signal at ZT22 for all three marks studied ( Fig 1G). Finally, all 127 DHSs associated with Npas2 (those having Npas2 as their closest TSS), including numerous 128 distal regions, likewise cycled with phases around ZT22 (Fig 1H) (Fig S2E). In the promoter region of Rev-159 erbα (Nr1d1), the detected footprints coincided with E-boxes and high BMAL1 ChIP-seq signals 160 (Fig 2A). Overall, the majority (70%) of DHSs within 1 kb of a TSS contained at least one 161 footprint, while this proportion dropped to one half for proximal (defined as DHSs within 1-10kb 162 of a TSS) or distal (>10kb of a TSS) DHSs (Fig 2B) Fig 3B). Moreover, the peak times of the oscillations in DNase 183 I signals were, except for some small deviations, similarly distributed as peak times in gene 184 transcription and H3K27ac [6,11,12], with a weak evening peak around ZT10 and a marked late 185 night peak around ZT22 (Fig 3C). We next considered the relationships of peak times in the 186 DNase I, Pol II and H3K27ac rhythms. It is known that many chromatin marks exhibit diurnal 187 rhythms that are tied to transcription [6,11,12,16], and similarly, enhancer RNAs (eRNAs) were 188 shown to be transcribed in sync with their cognate transcripts [28]. We observed that DNase I 189 cuts, Pol II, and H3K27ac displayed synchronous oscillations at DHSs (Fig 3D). Such To understand how the circadian clock and the feeding-fasting cycle control diurnal gene 204 expression in liver, we studied mRNA expression and Pol II loading at TSSs in WT and Bmal1 -/-205 mice subject to the same, night restricted, feeding regimen ( Fig S4). First, we observed that a 206 similar number of genes oscillated in the WT and Bmal1 -/genotypes (p < 0.05), however, with an 207 overlap of about 30% for Pol II and 50% for mRNA. This indicates that genes with a diurnal 208 expression differ between WT and Bmal1 -/mice ( Fig S4A). While such comparisons are based on 209 cutoffs, stratifying by peak-to-trough amplitudes clearly showed that high amplitude rhythms are 210 more abundant in WT as compared to Bmal1 -/mice (Fig S4B), and that this was more 211 pronounced for mRNA than for Pol II loading at TSSs. For example, we found twelve genes with 212 greater than ten-fold mRNA amplitudes in WT, and only three in Bmal1 -/mice. Genes with Pol II 213 or mRNA rhythms in both genotypes showed highly correlated phases, with a tendency for a DHSs. This enabled us to take into account, in addition to the proximal promoter, a collection of 9 putative regulatory regions that may control the expression of a given gene ( Fig 4A). Specifically, 233 we considered motifs in DHSs located within a certain window from active promoters, and first 234 estimated the optimal window size according to the quality of the fit. We found that the inclusion 235 of DHSs up to 50kbp was improving the fits in both genotypes (Fig 4B)

244
In Bmal1 -/mice ( Fig 4D, Table S4), activities of E-Box, RRE, and D-Box motifs were not 245 detected or greatly reduced, as expected in the absence of a functional circadian oscillator. On the 246 other hand, transcription factors linked with metabolic functions, notably those associated with 247 feeding rhythms (e.g. FOX, CREB, SREBP) were identified among the strongest contributors in 248 the absence of a functional clock. Similarly, transcription factors whose activity depends on 249 systemic signals (e.g. GR and HSF1) were also found with peak activity times that were similar in 250 the WT and Bmal1 -/mice. Interestingly, CREB was found among the most delayed transcription 251 factor activities, with a predicted delay of six hours (Table S4). To test this prediction, we 252 measured nuclear levels of CREB and pCREB using Western blots of nuclear extract from four 253 independent livers every two hours in WT and Bmal1 -/mice ( Fig 4E and Fig S5). On average, we 254 observed a phase delay of approximately two hours in Bmal1 -/mice. Although this was not 255 significant (p=0.5, Chow test), presumably owing to inter-individual variability in the feeding   Fig 6B). Also, the binding dynamics of BMAL1 at E1-E2-sp7 (tandem E-boxes 295 separated by 7bp) was largely similar to that for E1-E2-sp6, though E1-E2-sp7 had both E-boxes 296 predominantly protected only at ZT6, suggesting spacer-specific binding dynamics ( Fig S7). In 297 contrast, the footprints at BMAL1 binding sites with single E-boxes did not show significant 298 changes in time or in the Bmal1 -/mice (Fig S8), again suggesting that other bHLH transcription 299 factors can also bind at BMAL1 sites. In fact, footprints at DNA regions bound by the bHLH 300 11 transcription factor USF1 in ChIP-seq [73] were largely similar to that of BMAL1 sites with 301 single E-boxes, though the fraction of sites with clear footprints was reduced for USF1 compared 302 to BMAL1 (Fig S9).

303
To better understand the time-dependent footprint at BMAL1 sites and to gain insight into how 304 the CLOCK:BMAL1 heterodimer occupies its tandem E-box-containing target sites, we used 305 recently established 3D protein structures of single BMAL1/CLOCK complexes combined with 306 molecular modeling (Methods). Our models strongly support formation of CLOCK:BMAL1 307 heterodimers in a hetero-tetramer configuration at peak activity of these factors, and residual 308 binding of the dimer or other transcription factors during low activity times. Two 3D models of 309 the hetero-tetramer configuration were constructed. In the first model, the spacing between the 310 two E-boxes was 6 bp (sp6) ( Fig Table S4). Radial scale for activities is arbitrary but 725 comparable in C and D. 726 E. Quantification of Western blots for pCREB and CREB in WT and Bmal1 -/genotypes (log 2 727 (pCREB/CREB)). Nuclear extracts from four independent livers were harvested every two hours.