Genome-Wide Analysis of SREBP1 Activity around the Clock Reveals Its Combined Dependency on Nutrient and Circadian Signals

To evaluate the genome-wide dynamics of SREBP1 binding to its target sites in a physiological context, we prepared liver chromatin from C57BL/6 mice, collecting samples each 4 hours during one day (see Material and Methods). ChIP-seq with an antibody that recognizes both SREBP1 isoforms was performed at each time point. The SREBP1 antibody was tested extensively (Figure S1) and has also been used in previous studies [17]. We obtained an average of 38 millions sequence reads by time point by ultra-high-throughput sequencing (Table S1). The mapping allowed the identification of 448 bona fide SREBP1 binding peaks, above the background. As shown in Figure 1A, the binding of SREBP1 is overall oscillatory, with the maximum for most of the sites at Zeitgeber Times (ZT) 14 or 18 (light is on at ZT0 and is off at ZT12). To systematically evaluate the rhythmicity of SREBP1 recruitment to its targets, a cosine function was fitted to the temporal profile of the binding (see Material and Methods). This allowed to calculate, for each peak, the binding phase and the amplitude of the oscillation together with its associated P-value, as exemplified for the two sites found on the Srebp1 gene itself (Figure 1B). 53% of SREBP1 binding sites were found to be rhythmic (P<0.1 for the amplitude). Four clusters of targets were clearly distinguishable based on binding kinetics (Figure 1A), the first with a phase distributed around ZT15–ZT17, whereas the other ones with the phase peaking around ZT11–ZT12, as determined by cosine function (Figure 1C). The observed kinetics of SREBP1 binding, especially for cluster A peaks, was consistent with its gene expression and nuclear localization, (Figure 1D, 1E and 1F). Collectively, these results show that the activity of SREBP1 oscillates with a pronounced circadian rhythm, in agreement with the previously reported daily variations of its RNA and protein levels [16], [18]–[20].

SREBP1 binding sites identified in this study are grouped in four clusters with a slightly shifted phase. To better investigate the features of these sites we calculated the number of nucleotides spanned by each peak and found that sites belonging to cluster A (236 out of 448) were narrow, with a typical length of about 200 nucleotides (Figure 2A). These peaks were also closer to the nearest annotated transcription start site (TSS) than peaks belonging to clusters B, C or D, whose distance to the nearest annotated TSS roughly matches randomly picked genomic locations (Figure 2B). Moreover, the amplitude of the binding oscillation along the 24 hours was greater for cluster A peaks (Figure 2C). These observations suggest that cluster A sites may be more relevant in the regulation of transcription mediated by SREBP1. In agreement with this hypothesis, a MEME [21] motif search analysis clearly identified the canonical SREBP1 consensus motif in more than 60% of the sites belonging to cluster A (Figure 2D), but only in 6% of the sites assigned to clusters B, C and D. Within cluster A, the MEME analysis also showed that motifs for SP1 and NFY, two transcription factors (TFs) known to cooperate with SREBP1 to regulate the transcription of its target genes, were overrepresented [17]. The consensus motif for the Hepatic Nuclear Factor 4 (HNF4) was also identified in 61 out of the 236 cluster A sites. The discovered motifs were enriched in cluster A peaks with an empirical p-value<0.001, as shown in Table S2. In addition, in the regions belonging to cluster B, C and D, we determined only several highly repetitive sequences as top-scoring motifs (for example ACACACACA in 73 sites out of 212) that could not be associated to any known consensus motif for TFs. This result suggests that SP1, NFY and HNF4 may participate to SREBP1-mediated transcriptional regulation and further supports the functional importance of SREBP1 binding sites assigned to cluster A. Thus, we opted to focus the following analyses on these regions, although we cannot exclude that the other sites might contribute to mediate SREBP1 activity in mouse liver potentially through genome loops.

To explore the cellular processes that are regulated by SREBP1 along the day, we annotated each site with the nearest Ensembl transcript. We used DAVID [22], [23] to identify clusters of genes enriched with functional annotations. As expected, we identified lipid biosynthetic processes and fatty acid metabolism as the most prominent pathways controlled by SREBP1 (Table S3). In addition, we found a significant enrichment of genes involved in carbohydrate metabolism, in the response to nutrient levels, in mitochondrial and endoplasmic reticulum functions and in coenzyme metabolism. In a previous genome-wide study performed in human HepG2 cells, it was shown that unique combinations of SREBP1, SP1 and NFY target distinct functional pathways [17]. Since we found a good enrichment within the SREBP1 binding sites of the consensus motifs for NFY and SP1, but also HNF4, we explored whether a network among these three transcription factors could be highlighted in mouse liver. As shown in Table 1, genes involved in lipid biosynthesis and in the regulation of fatty acid and steroid metabolism were highly represented in all categories. In some cases, however, one biological function was targeted by a unique combination of regulators. For example, the biosynthesis of coenzymes was selectively represented within the genes bearing only the site for SP1, whereas both SP1 and HNF4 motifs were present in genes involved in apoptosis. Likewise, a combination of HNF4 and NFY motifs marked most of the genes involved in immunological processes. Finally, pathways related to carbohydrate metabolism and mitochondria were particularly enriched in genes without NFY or SP1 motifs, suggesting that SREBP1 may cooperate with other regulators at the promoter of these genes (the complete functional annotation clustering is reported in Table S4). Our results suggest that in mouse liver, in physiological conditions, the network SREBP1-SP1-NFY-HNF4 may be important in order to determine the functional effect of SREBP1 binding (Figure 2E).

In Figure 1, we showed that SREBP1 binds to target sites belonging to cluster A with a sharp phase between ZT15 and ZT17. To investigate the functional effects of SREBP1 binding on gene transcription, we checked in our previously reported data set [24] the 24 hours profile of RNA polymerase II (Pol II) recruitment in the proximity of SREBP1 target genes. Importantly, most of the SREBP1 binding sites belonging to cluster A were co-occupied by Pol II (Figure 1A), further supporting the functional relevance of these regions. We next evaluated Pol II binding to the promoter and in the gene body of all putative SREBP1 target genes. In parallel, we measured mRNA levels of the same genes by microarray analysis (Table S5). More than 85% of SREBP1 target genes show an expression level above the median expression level of all the transcripts, suggesting that they are transcribed. Our analyses revealed three clusters of target genes, that we called A1, A2 and A3, with distinct temporal profile of transcription and expression (Figure 3A and Figure S2). In cluster A1, the peak of Pol II binding was concomitant, or even slightly earlier than SREBP1 binding. In contrast, for genes belonging to cluster A2, Pol II association to both promoter and gene body strictly followed SREBP1 binding. Lastly, Pol II recruitment to the genes of the A3 group was shifted by about +8 h with respect to SREBP1. For all clusters, the temporal profile of gene expression was consistent with the dynamics of Pol II association. The distribution of all the expression phases obtained for the genes belonging to the three groups confirmed that SREBP1 target genes are expressed in different moments of the day, in spite of the concomitant binding of the transcription factor (Figure 3B). This observation suggests that other factors participate in the regulation of the various SREBP1 target genes in order to assure their appropriate expression timing. Interestingly, genes that were mainly expressed during the fed state (clusters A1 and A2) were functionally enriched in the regulation of lipid and coenzyme biosynthetic processes, as well as in the response to hormones, such as insulin. In contrast, SREBP1 target genes involved in mitochondrial oxidation and apoptosis were enriched during the fasting period (Table 2). Thus, the promoter specific events that determine the different temporal expression profile of SREBP1 target genes contribute to define the set of cellular functions that are active at a given time.

To understand the molecular mechanism underlying the different temporal expression of SREBP1 target genes, we first explored the possible involvement of the network SREBP1-SP1-NFY-HNF4 in determining the functional effect of the binding of SREBP1 to its targets. To check for the presence of a pattern characterizing the three groups of SREBP1 target genes identified earlier (see Figure 3), we evaluated the presence of different combinations of SP1 and NFY motifs and their orientation with respect to the SREBP1 binding sites (data not shown). However, we could not establish any significant correlation. In contrast, we found that HNF4 motifs were significantly overrepresented (P-value<0.02) in the regions under SREBP1 peaks of the genes expressed during the fasting period (cluster A3), compared to the other clusters (Table S6). The actual recruitment of HNF4 to these putative binding sites was assessed by ChIP on randomly selected SREBP1 Responsive Elements (SREs) (Figure 3C). Besides HNF4, other transcription factors, such as the cAMP response element-binding protein (CREB) or Forkhead box proteins O (Foxo), are important players in the hepatic metabolic regulation upon fasting [25], [26]. However, their known consensus motifs were not found in the proximity of cluster A3 SREBP1 peaks. These observations strongly suggest a specific cross-talk between HNF4 and SREBP1 in the regulation of these genes. To further investigate which control processes dictate the distribution of SREBP1 target gene expression along the day, we then considered the possible role of the circadian rhythm in this regulation. To test this hypothesis, it was necessary to uncouple the circadian rhythm from the response to nutrients. Thus, we fed mice lacking BMAL1 (Bmal1^−/−) only during the darkness period for one week before collecting liver samples every four hours. Upon this experimental conditions, the circadian clock was completely disrupted in Bmal1^−/− mice, as demonstrated by the flattened expression of key core and output components of the clock, such as Clock1, Cry1, Cry2, D site albumin promoter binding protein (Dbp), Rev-Erbα and Kruppel-like factor 10 (Klf10) (Figure 4B). Body weight, daily food intake and glycemia were unchanged in Bmal1^−/− mice (Figure S3). Importantly, the imposed rhythmic food intake restored an oscillatory nutrient response, as shown by the levels of circulating insulin that were comparable to the wild type (Figure 4A). Accordingly, in Bmal1^−/− mice SREBP1 translocated to the nucleus from ZT18 onwards (Figure 4D and 4E) and its expression was still cycling, although the phase was delayed by about 6 h compared to the wild type (Figure 4C). We next checked the dynamics of SREBP1 binding to a panel of the targets previously identified in wild type mice and we found that SREBP1 was recruited to all the tested sites in an oscillatory way, but with an average phase shift of about 4 hours (figure 4F). Finally, to assess the impact of the circadian oscillator impairment on SREBP1-driven transcription, we globally evaluated the expression of SREBP1 target genes in Bmal1^−/− mice (Table S5). Most SREBP1 target genes were still scored as oscillating in Bmal1^−/− upon temporal restricted feeding (57% have a P<0.05). However, the heatmap rendering of their expression patterns (Figure 5A) revealed a temporal profile that was perturbed in Bmal1^−/− compared to WT mice, as most of the genes now had a maximum expression at ZT18, coinciding with the binding of the transcription factor (Figure 5C). Accordingly, the expression phases of the genes belonging to the clusters A1 and A3 identified earlier (Figure 3) were now largely concomitant with those of genes belonging to cluster A2, therefore mostly grouped between ZT14 and ZT24 (Figure 5B and Figure S4). This phase shift was not due to a selective decrease of the number of cycling genes in clusters A1 and A3, as the percentage of the significantly oscillating genes was comparable in the three clusters in WT and Bmal1^−/− mice (cluster A1: 84% in WT vs 94% in Bmal1^−/−; cluster A2: 71% in WT vs 66% in Bmal1^−/−; cluster A3: 62% in WT vs 64% in Bmal1^−/− with P<0.05).

To explore how the core clock components participate to this regulation we checked in published data sets whether key transcription factors such as BMAL1, CLOCK1, CRY1, CRY2, PER1, PER2, NPAS2 and REV-ERBs can differentially bind to the promoters of the three clusters of SREBP1 target genes [27]–[29]. Interestingly, we found that most SREBP1 peaks (148 out of 236, ≈63%) have an overlapping REV-ERBα and/or REV-ERBβ peak. In addition, a REV-ERBα/β binding site was detected also in another 17% of the promoters of SREBP1 target genes, but in a non-overlapping position. This strong occupancy of SREBP1 targets by REV-ERBs is consistent with the previously reported involvement of REV-ERBα in the regulation of lipid metabolic genes [30], and suggests the existence, in mouse liver, of a SREBP1-REV-ERBs network in physiological conditions. The frequency of REV-ERBs recruitment was comparable in clusters A1, A2 and A3 (data not shown), thus arguing against the possible role of these nuclear receptors in determining the distinct phase of expression of these genes. However, due to the presence of REV-ERB binding sites in many SREBP1 target genes, the flattened REV-ERB expression observed in Bmal1^−/− may perturb, at least in part, the phase of several SREBP1 target genes. Indeed, in WT mice, the temporal expression profile of SREBP1 and REV-ERBs is very different (the phases of expression are ZT15 and ZT8, respectively [16]), and these factors are not expected to compete for binding at the same time to the same genes. The other transcription factors tested were recruited to a lesser extent on SREBP1 target gene promoters and for none of them we observed a significant enrichment in clusters A1, A2 or A3 (data not shown).

Taken together, our results confirm that SREBP1 activity is strongly dictated by the rhythmicity of nutrient intake. In addition, our observations indicate that a functional circadian core clock is necessary to assure the correct temporal expression profile of SREBP1 target genes and suggest a role for HNF4 in dictating the phase of expression of genes whose mRNA levels peak when SREBP1 binding is low. Further studies will aim at understanding whether and how the lack of circadian rhythm perturbs HNF4 activity.

All animal experiments and procedures were approved by the Swiss Veterinary Office (authorisation VD-1453.4). C57BL/6 male were purchased from Charles River. Bmal1^−/− mice were a kind gift from Dr. Frédéric Gachon and were generated as previously described [63], [64]. 12–14 week old (at time of sacrifice), mice were housed in a 12 h light/12 h dark (LD) regimen for 2 weeks with food and water freely available during night and day. They were then phase-entrained to a 12 hr/12 hr LD regimen with food access between ZT12 and ZT24 for 7 days (ZT = Zeitgeber time; ZT0 is defined as the time when the lights are turned on and ZT12 as the time when lights are turned off). At each ZT2, ZT06, ZT10, ZT14, ZT18 and ZT22 three to five mice were anesthetized with isoflurane and decapitated. Mice were killed under dim red light at ZTs during the dark phase. The livers were perfused with 2 ml of PBS through the spleen and immediately collected. A small piece of liver tissue was snap-frozen in liquid nitrogen. The remaining liver tissue was immediately homogenized in PBS containing 1% formaldehyde for chromatin preparation.

Perfused livers were processed for chromatin preparation as previously described [65]. The chromatin samples from the mice of the same ZT were then pooled, frozen in liquid nitrogen and stored at −80°C. The following antibodies were used: anti-RPB2 (Santa Cruz Biotechnology, H-201), anti-SREBP1 (Santa Cruz Biotechnology, H-160), anti-HNF4 (Santa Cruz Biotechnology, C-19). Chromatin was subjected to immunoprecipitation of Pol II as described [24]. For SREBP1, the samples were diluted ten times in “sonication buffer” containing 50 mM HEPES (pH 7.9), 140 mM NaCl, 1 mM EDTA, 1% Triton X-100, 0.1% Na-deoxycholate, 0.1% SDS and proteinase inhibitors (Roche). 1 ml of diluted chromatin was immunoprecipitated with 10 µg of antibody as described [66]. Briefly, the immune complexes were collected by adsorption to ten µl of protein-A-Sepharose (25% slurry in sonication buffer), pre-blocked with 10 µg/ml of salmon sperm DNA and BSA at 4°C overnight. The beads were washed twice with “sonication buffer”, twice with sonication buffer containing 500 mM NaCl, twice with 20 mM Tris, pH 8.0, 1 mM EDTA, 250 mM LiCl, 0.5% NP-40, 0.5% Na-deoxycholate and twice with TE buffer. The immunocomplexes were eluted with 50 mM Tris, pH 8.0, 1 mM EDTA and 1% SDS at 65°C for 10 min., adjusted to 200 mM NaCl and incubated at 65°C overnight to reverse the cross-links. After successive treatments with 10 µg/ml Rnase A and 20 µg/ml proteinase-K, the samples were extracted with NucleoSpin Kit (Macherey- Nagel). The DNA concentration was determined by fluorometry on the Qubit system (Invitrogen). 10–12 ng DNA were used for the preparation of the library. Libraries for ultra-high throughput sequencing were prepared with the ChIP-Seq DNA sample kit (Illumina) as recommended by the manufacturer.

About 100 mg of snap-frozen liver tissue were used for RNA preparation with the TRIzol reagent (Invitrogen) followed by purification with miRNeasy Mini Kit (Qiagen), according to manufacturer's instructions. For microarray analysis 500 ng of total RNA from each liver sample at the same time point were pooled and analyzed on Mouse Gene 1.0ST arrays according to the manufacturer's instructions (Affymetrix). All statistical analyses were performed with the statistical language R and various Bioconductor packages (http://www.Bioconductor.org↗). Normalized expression signals were calculated from Affymetrix CEL files using RMA normalization method. For quantitative RT-PCR analysis, the retrotranscription has been done using iScript cDNA synthesis kit (Bio-Rad, Laboratories, Hercules, CA) and following the manufacturer's instructions. The primers sequences are shown in Tables S7 and S8. Real-time monitoring of PCR amplification of cDNA was performed using the FastStart Universal SYBR Green Master (Roche Applied Science, Indianapolis, IN) in an ABI Prism 7900 Sequence Detection System (Life Technologies, Carlsbad, CA). The PCR arbitrary units of each gene were defined as the mRNA levels normalized to the 36b4 and the Rps9 expression level in each sample using the qBase Software.

Nuclear extracts were prepared by the NUN procedure as described previously [67], and Western blotting was performed according to standard protocols using the antibody for SREBP1 indicated above. U2AF and Lamin A were used as loading control (anti-U2AF and anti-Lamin A were from Sigma-Aldrich).

At sacrifice, blood was taken for determination of biochemical parameters and circulating hormones. Insulin levels were determined with ELISA kit from Mercodia, Uppsala, Sweden, following manufacturer's instructions.

At each time point, DNA sequenced reads were mapped to the mouse genome (Mus musculus NCBI m37 genome assembly (mm9; July 2007)) using Bowtie [68] with three mismatches and at most five hits allowed on the genome. When computing genomic read densities, each alignment contributed 1/(total number of hits) to the local density. If several reads in the same library mapped at the same genomic position and on the same strand (redundant tags), we kept only one read for the rest of the analysis. The total numbers of reads per time point are given in Table S1.

The mapped reads were shifted to account for the length of the inserts based on the average fragment size, namely 190 184 204 199 202 208 and 176, for each of the seven libraries. The fragment size was divided by two and half of the read length was subtracted resulting in a shift of 55, 52, 62, 60, 61, 64, and 48 nucleotides from ZT02 to ZT26.

We fragmented the genome into 500 nucleotide blocks and collected the counts within each block for SREBP1 as well as for input experiments. We kept all blocks for which we had a signal equivalent to 40 tags in at least one time point of the SREBP1 experiment. We log transformed the data after adding 1 pseudo count and quantile normalized both SREBP1 and input experiments. We selected blocks with a log2 signal SREBP1/input greater than 2, i.e. at least a four-fold enrichment of the SREBP1 signal in comparison to the input signal in at least one time point. We repeated this procedure by shifting the block definition of half the length of the blocks and merged the overlapping blocks that passed these criteria.

To define proper “peaks” in these wide regions, we look for shorter regions accounting for most of the counts. For this, we repeatedly consider the two borders along 50 nucleotides, and discard the one with less read counts if we keep 75% of the total reads in the remaining region. This operation is repeated with shorter borders until no further refinement is possible.

We used the MEME suite [21] to identify enriched motifs in the sequences corresponding to the refined peaks. We first performed several motif discovery, on all peaks, and on the subset associated or not with Pol II. We searched for 15 motifs between 6 and 10 nucleotides long. The discovered motifs have been associated to known transcription factors (from the TRANSFAC database) with STAMP. To retrieve these motifs in the subsets where they have not been discovered, we searched for them in all sequences (using FIMO). To assess the relevance of the number of observed motifs in our dataset, we counted the occurrence of the same motifs in random sequences. These sequences are selected in the proximity of the TSS of genes expressed in our samples (among the top 10% in the microarray data). Their size and distance to TSS are in the same range as that of our SREBP1 peaks.

Before fitting a cosine function to estimate the amplitude and the phase of the oscillation in a 24 hour period, counts in the refined peaks were quantified and normalized according to the total number of non redundant mapped reads for each given library.

We used the function x(t) = b0+b1*cos(b3+2π*t/24) to perform a least squared fitting of temporal profiles. The parameter b0 represents the mean signal, b1 the amplitude of the oscillation, and b3*24/2π the phase. These parameters were estimated by nonlinear least-squares using the Gauss-Newton algorithm.

For the microarray temporal profile analysis we used the function

while when we compared WT and Bmal1^−/− mice we used the function

where b4 and b5 represent the different batch effects, GT is a dummy variable that indicates the mouse genotype (WT or KO) and b0_gt, b1_gt, b3_gt the associated coefficients.

Illumina sequencing data for the ChIP-seq are available at GEO as the GSE48375↗. Additional processed data and visualization tools are provided at http://cyclix.vital-it.ch↗.

Illumina sequencing data for the ChIP-seq are available at GEO as the GSE48375↗. Additional processed data and visualization tools are provided at http://cyclix.vital-it.ch↗.