Introduction

Growth hormone (GH) regulates hepatic expression of enzymes and transporters that play critical roles in lipid metabolism and in the detoxification of many drugs and other lipophilic foreign chemicals [1, 2]. Dysregulation of hepatic GH signaling can lead to liver metabolic disorders, including the development of fatty liver disease and non-alcoholic steatohepatitis, with males more susceptible than females, as seen in both mice and humans [3–5]. Correspondingly, GH regulates many liver-expressed genes in a sex-dependent manner, enabling each sex to meet its specific metabolic and hormonal requirements [6]. This program of sex-dependent gene expression is controlled by the sex-dependent temporal patterns of pituitary GH secretion [7, 8], which emerge at puberty but are programmed earlier in life by neonatal exposure to androgen [9–11]. Pituitary GH secretion, and plasma GH profiles, are intermittent (pulsatile) in pubertal and adult males, whereas they are near-continuous (persistent) in pubertal and adult females, as seen in rats, mice and humans [12]. These plasma GH profiles, in turn, regulate sex-specific gene transcription in the liver through both positive regulatory mechanisms (class I sex-biased genes) and negative regulatory mechanisms (class II sex-biased genes) [6].

The sex-dependent hepatic actions of GH require the transcription factor STAT5 [13], which is activated by phosphorylation on a single tyrosine residue catalyzed by the GH receptor-associated tyrosine kinase JAK2 [14]. Each successive male plasma GH pulse induces a cycle of STAT5 tyrosine phosphorylation, dimerization and nuclear translocation, followed by STAT5 tyrosine dephosphorylation and return back to the cytosol in time to reset the overall signaling pathway for the subsequent plasma GH pulse [1]. In male mouse liver, GH pulse-activated STAT5 binds to the DNA motif TTCNNNGGA at open chromatin regions strongly enriched for proximity to male-biased genes, but the relationship between STAT5 binding, chromatin opening and sex differences in chromatin accessibility is unknown. In female liver, persistent activation of STAT5 by the near-continuous presence of circulating GH is associated with a significant enrichment of STAT5 binding nearby female-biased genes [15], however, the underlying mechanisms, including chromatin features that distinguish these STAT5 binding sites from male-biased STAT5 binding sites are poorly understood.

Liver STAT5 is activated by male plasma GH pulses within minutes, enabling STAT5 to rapidly induce the transcription of several male-biased genes [16]; however, a majority of male-biased genes respond slowly to changes in GH secretory patterns [17], likely reflecting the time required for secondary changes, including changes in histone modifications and the underlying chromatin state [18]. Continuous infusion of GH in male mice overrides the endogenous plasma GH pulses and substantially feminizes liver gene expression over a period of days, with individual female-biased genes already in an active chromatin state in male liver generally responding earlier than genes in an inactive chromatin state [17]. These chromatin state changes include changes in chromatin accessibility at genomic regions discovered as DNase-I hypersensitive sites (DHS), a hallmark of epigenetic regulation. DHS can be assayed by DNase-seq, which has identified ∼70,000 DHS in male and female mouse liver [19], including thousands of enhancers, promoters, and insulator regions [20]. More than 4,000 of the 70,000 liver DHS show sex differences in chromatin accessibility, as well as sex-biased binding of STAT5 and other GH-regulated transcription factors, which reinforce sex differences in liver gene expression [18, 19] and regulate downstream sex differences in disease susceptibility [21–23]. Importantly, male-biased transcription factor binding is strongly enriched at male-biased DHS nearby male-biased genes, and female-biased transcription factor binding is strongly enriched at female-biased DHS nearby female-biased genes. Sex differences in chromatin structure and accessibility are thus key features of sex-differential liver gene expression [17, 18]. However, little is understood about the mechanisms linking the sex-dependent temporal GH secretion patterns to the robust sex differences in chromatin accessibility and transcription factor binding that regulate liver gene expression.

Here, we elucidate the relationship between the sex-dependent patterns of plasma GH stimulation of hepatocytes and sex differences in liver chromatin accessibility. We identify more than 800 male-biased enhancer DHS regions, where the direct binding of plasma GH pulse-activated STAT5 induces a dynamic cycle of male liver chromatin opening and closing at sites that comprise 31% of all male-biased DHS. Thus, the pulsatility of plasma GH stimulation per se confers significant male bias in chromatin accessibility and STAT5 binding at a substantial fraction of the genomic sites linked to sex-biased liver gene expression and liver disease.

Furthermore, we establish that a single physiological replacement dose of GH given as a pulse to hypophysectomized (hypox) mice recapitulates, within 30 min, the pulsatile re-opening of chromatin seen in pituitary-intact male mouse liver. We elucidate key epigenetic features distinguishing this dynamic, STAT5-driven mechanism of male-biased chromatin opening from that operative at a second, distinct class comprised of static male-biased DHS, which are constitutively open in male liver but closed in female liver. Finally, our analysis of repressive histone marks enriched at each class of sex-biased DHS indicates that distinct epigenetic mechanisms mediate sex-specific gene repression in each sex. Together, these studies identify pulsatile chromatin opening as a novel mechanism controlling sex differences in chromatin accessibility and transcription factor binding closely linked to sex differences in gene expression and liver disease.

Results

Sex-biased DHS are primarily distal enhancers that target sex-specific genes

– Open (accessible) chromatin regions identified in mouse liver by DNase-seq analysis have been classified as enhancers, insulators and promoters based on their chromatin mark patterns and CTCF binding activities (n = 70,211 standard reference DHS set) [20] (Table S2A). Several thousand of these DHS show greater chromatin accessibility in male than female liver (n = 2,729 male-biased DHS) or greater accessibility in female than male liver (n = 1,366 female-biased DHS) [19]. These sex-biased DHS regions were enriched for enhancer marks, with >85% classified as enhancers or weak enhancers in mouse liver, as compared to 66% for all sex-independent DHS (Fig. 1A).

Fig. 1.

Assignment of the sex-biased liver DHS to their putative gene targets (closest RefSeq gene or multi-exonic lncRNA gene transcription start site in the same topologically associating domain (TAD) [20] revealed that sex-biased DHS were highly enriched for genes of the corresponding sex bias, but not for genes of the opposite sex bias (Table S2B). Cumulative plots of the distance from each DHS to its target gene revealed that a majority of sex-biased DHS classified as enhancers or insulators are distal to their target genes. Sex-biased DHS with H3K4me3 marks, which comprise only 1-2% of all sex-biased DHS, are proximal to their target genes, validating their classification as promoter DHS (Fig. 1B). Thus, a majority of sex-biased DHS are marked as positively acting distal enhancers, consistent with our finding that sex-dependent DNA looping at an intra-TAD scale is common in mouse liver [24].

No major differences in chromatin state distributions between male-biased, female-biased and sex-independent enhancer DHS were seen based on separate chromatin state maps developed for male and female liver [18]. Thus, all three sets of enhancer DHS showed a high frequency (50-66%) of chromatin state E6, whose emission probabilities (Fig. 1C) indicate a high frequency of DHS in combination with the activating chromatin marks H3K27ac and H3K4me1, and a lower frequency (12-16%) of chromatin state E5 (high frequency of DHS but low frequency of H3K27ac and H3K4me1) (Fig. 1D). Sex-biased and sex-independent insulator DHS also showed similar chromatin state distributions (state E5 ∼ state E6), as did sex-biased and sex-independent promoter DHS, which were primarily in state E7, characterized by a high frequency of H3K4me3 and H3K4me1 marks (Fig. 1D).

The absence of major differences in overall chromatin state distributions between male-biased and female-biased liver DHS prompted us to investigate other factors that may provide insight into underlying mechanisms regulating sex-differences in hepatic chromatin accessibility and their link to sex differences in gene expression.

Impact of endogenous pulses of GH and STAT5 activity in male liver

- Sex differences in pituitary GH secretion – pulsatile in males vs. near continuous (persistent) in females – regulate the sex-dependent expression of hundreds of genes in adult mouse liver. This regulation requires the GH-activated transcription factor STAT5 [25, 26]. We hypothesize that the pulsatile activation of STAT5 seen in male liver in direct response to plasma GH stimulation [15, 16, 27] dynamically alters the chromatin accessibility landscape of male mouse liver. More specifically, we propose that the repeated activation of STAT5 in male mouse liver by plasma GH pulses induces dynamic cycles of chromatin opening and closing at a subset of liver DHS, and that this response can be discovered by comparing chromatin accessibility profiles in livers of individual male mice euthanized at a peak vs. at a trough of hepatic STAT5 activity (Fig. 2A, top). To test this hypothesis, we analyzed DNase-seq libraries prepared from liver nuclei purified from 21 individual adult male mice. In parallel, we determined the STAT5 DNA-binding activity of each liver by EMSA (electrophoretic mobility shift analysis) of whole liver extracts.

Fig. 2.

Ten of the 21 livers had high STAT5 EMSA activity (STAT5-high activity livers) and 11 livers had very low or no detectable STAT5 EMSA activity (STAT5-low activity livers) (Fig. 2B, Fig. S1). Next, we performed DNase-seq analysis on genomic DNA fragments released by light DNase-I digestion of nuclei purified from each liver, followed by diffReps analysis [28] comparing the DNase-seq-released DNA fragments from each group of livers. We thus discovered genomic regions showing significant differential chromatin accessibility between STAT5-high and STAT5-low male livers. Principal component analysis using either the top 200 or the top 600 most significant diffReps-identified differentially accessible regions revealed that 18 of the 21 livers yielded consistent, STAT5 activity-dependent patterns of DNase-released fragments. Two of the STAT5-high livers and one of the STAT5-low livers were outliers (Fig. 2C). These outliers were also identified by their discordant DNase-seq read count distributions (Fig. 2D) and were excluded from all downstream analyses.

We re-analyzed the remaining 18 STAT5-high and STAT5-low livers using diffReps then filtered the output list to retain regions identified as DHS peak regions by MACS2 analysis of the same 18 DNase-seq datasets. We thus identified n=2,832 genomic sites where chromatin opening is associated with high liver STAT5 activity and n=123 other sites where chromatin opening is associated with low liver STAT5 activity (Table S4A). As each liver represents a time point of peak (STAT5-high livers) or trough (STAT5-low livers) levels of GH pulse-activated liver STAT5 activity (Fig. 2A), the STAT5-high/STAT5-low differential sites identify genomic regions where chromatin dynamically opens or closes in male mouse liver in response to GH pulse activation of STAT5. The greater chromatin opening in STAT5-high livers as compared to STAT5-low male livers, and compared to female livers, was visualized in aggregate plots of normalized DNase-seq cuts across the 2,832 genomic regions (Fig. 2E, left). In contrast, chromatin accessibility was greater in the STAT5-low livers at the 123 STAT5-low diffReps-identified genomic sites, where STAT5 activation is associated with chromatin closing (Fig. 2E, right). This binary pattern, where individual male mouse livers show either high or low DNase-seq read count distributions at the top differential genomic sites, was also seen in an independent set of nine male liver DNase-seq samples from a second mouse strain generated by the ENCODE consortium: 4 of the 9 livers showed read distributions very similar to the STAT5-high livers, while 5 of the 9 livers were similar to the STAT5-low livers (Fig. 2D, green bars).

Dynamic vs static liver DHS

– Comparison of the set of 2,832 STAT5-high genomic sites with our standard reference set of 70,211 liver DHS revealed that 31% (n=834) of the 2,729 male-biased liver DHS described above are STAT5-high sites, i.e., they respond dynamically to endogenous male plasma GH pulses of liver STAT5, whereas only 0.5% of female-biased liver DHS and 2% of sex-independent liver DHS showed this dynamic STAT5 response (Fig. 2F, green). This strong enrichment of GH/STAT5 pulse-induced chromatin opening at male-biased DHS as compared to sex-independent DHS (ES = 12.9; p < 1E-05) identifies intermittent chromatin opening induced by male plasma GH pulses as a mechanism that can explain the male bias in chromatin accessibility for a significant subset (31%) of male-biased DHS.

To determine whether STAT5 binding per se is an important feature of the observed pulsatile changes in chromatin accessibility, we examined normalized DNase-I cut site aggregate plots for the subset comprised of n=1,307 male-biased DHS that bind STAT5 in ChIP-seq analyses. STAT5-high livers showed the highest level of chromatin opening, followed by up to a 2.8-fold lower level of chromatin opening at those same genomic regions in STAT5-low male livers and in female livers (Fig. 3A, plot 1). In contrast, male-biased DHS that did not bind STAT5 (n=1,422 DHS) showed nearly equal chromatin accessibility in STAT5-high vs. STAT5-low male livers but ∼2-fold lower accessibility in female livers (Fig. 3A, plot 2). Thus, the male bias in accessibility at the STAT5-bound but not at the non-STAT5 bound male-biased DHS can be explained by STAT5-induced pulsatile chromatin opening in male liver. Moreover, it is also apparent that the male bias of the latter DHS set is associated with chromatin closing in female liver (Fig. 3A, plot 2, black). The latter conclusion is supported by the low mean normalized DNase cutting frequency at the 1,422 male-biased DHS in female liver when compared to that of the full genome-wide set of 53,404 STAT5-unbound, sex-independent DHS (Fig. 3A, plot 2, black versus plot 4). Chromatin accessibility was much higher at the subset of sex-independent DHS that bound STAT5 (n=12,712 DHS, plot 3). Strikingly, the high chromatin accessibility at these sex-independent DHS was seen in livers of both sexes and was largely invariant between STAT5-high and STAT5-low male livers. We conclude that STAT5 binding is associated with increased chromatin opening, and that dynamic chromatin opening and closing associated with STAT5 binding is a defining characteristic of a distinct subset of male-biased DHS but occurs infrequently at female-biased and sex-independent DHS.

Fig. 3.

The above findings allow us to identify two distinct subsets of male-biased DHS, which differ in their responsiveness to liver STAT5 activation: 1) dynamic male-biased DHS are characterized by a male bias in chromatin accessibility linked to an increase in chromatin opening following the binding of STAT5 when activated by a plasma GH pulse in male liver every 3-4 h [1, 16] (n=834 dynamic male-biased DHS, 31%); and 2) static male-biased DHS are open in male liver constitutively, i.e., throughout the pulsatile, on/off cycles of GH-induced STAT5 activity, and are comparatively closed in female liver (n=1,895 static male-biased DHS, 69%). Consistent with this proposal, mean chromatin opening at the set of 834 dynamic male-biased DHS was 3.8-fold higher in STAT5-high male livers than in STAT5-low male livers or in female livers (Fig. 3A, plot 5). In contrast, the higher chromatin accessibility in male than female liver at the set of 1,895 static male-biased DHS was largely independent of the male liver’s STAT5 activity status. Of note, the mean level of chromatin opening at those 1,895 sites was 2.6-fold lower than the peak level of accessibility seen at the dynamic male-biased DHS regions in STAT5-high male livers (Fig. 3A, plot 6 vs. plot 5) (Table 1A).

Table 1:

Chromatin accessibility at the set of static female-biased DHS (1,359 sites; Fig. 2F) was 2.7 to 2.9-fold higher in female liver than in male liver, where chromatin accessibility was independent of male plasma GH/STAT5 pulses (Fig. 3A, plot 7). The aggregate female liver accessibility of these static female-biased DHS was very similar to that of the genome-wide set of 64,584 static sex-independent DHS (Fig. 3A, plot 8 vs. plot 7; peak normalized DNase-I activity value of 699 vs 683, Table S6).

STAT5 binding is closely associated with dynamic male-biased DHS

- The presence of a canonical STAT5 motif, TTC-NNN-GAA, is a distinguishing feature of dynamic male-biased DHS: 81% of 834 dynamic male-biased DHS contained a STAT5 motif versus only 38% of 1,895 static male-biased DHS (c.f., 25% background STAT5 motif frequency for genome-wide set of sex-independent static DHS) (Fig. 3B). In addition, multiple STAT5 motifs are more frequently found at dynamic DHS than at static DHS (Fig. 3C). Consistently, STAT5 binding, determined experimentally by ChIP-seq, occurs at a greater fraction of dynamic than static male-biased DHS (85% versus 32%; Fig. S2A), and the level of STAT5 binding (normalized STAT5 ChIP-seq read counts) was significantly higher at the STAT5-bound subsets of the dynamic vs. static DHS sets (Fig. S2B). De novo motif discovery supported these findings, with top scoring motifs matching the STAT5B motif found in dynamic DHS (Fig. S3A, Fig. S3B) but not in static DHS sets (Fig. S3C-S3E). A close association between STAT5 binding and chromatin opening is also indicated by the higher chromatin accessibility at the DHS subsets associated with STAT5 binding (Fig. S2C, top row vs. bottom row).

These findings establish a close link between STAT5 binding, which is pulsatile in male liver [15], and the repeated opening and closing of chromatin at dynamic DHS. Nevertheless, pulsatile chromatin opening was also found at 124 male-biased DHS (15% of all dynamic male-biased DHS) that did not bind STAT5 (Fig. S2C, plot 1B). Moreover, STAT5 binding alone is not sufficient to insure dynamic, male-biased chromatin opening and closing, insofar as 90% of all sex-independent DHS that bind STAT5 do not undergo dynamic chromatin opening and closing (Fig. S2C, plot 5A (11,473 sites, 90.3%) vs plot 4A (1,239 sites, 9.7%)). Furthermore, only 54% of male-biased DHS that bind STAT5 are dynamic DHS (710 sites), the other 46% (597 sites) being static male-biased DHS, even though they bind STAT5 (Fig. S2C, plot 1A vs 2A), albeit at a lower level than the dynamic male-biased DHS and similar to the STAT5-bound static sex-independent DHS (Fig. S2B). Of note, we observed a 5.6-fold higher frequency of dynamic chromatin opening and closing at the male-biased DHS with STAT5 bound (54%) than at the sex-independent DHS with STAT5 bound (9.7%). Thus, genomic features other than pulsatile STAT5 binding per se are required for dynamic chromatin opening and closing to occur.

Enrichment of other GH-regulated liver transcription factors at dynamic vs static DHS

– We investigated whether dynamic and static male-biased DHS differ with respect to the binding of other, sex-biased transcription factors (Table S1B). We considered two GH/STAT5-regulated transcriptional repressors that reinforce STAT5 regulation of liver sex differences. One factor, BCL6, is a male-biased protein that can compete for STAT5 binding to chromatin and preferentially represses female-biased genes in male liver [15, 18, 22, 29]. The second factor, CUX2, is a female-specific repressor whose binding sites in female liver are enriched nearby male-biased genes, which enables CUX2 to repress those genes in female liver [30, 31]. Both repressors showed significant enrichment for binding to genomic regions defined by male-biased DHS as compared to a background set of static sex-independent DHS. Thus, male-biased BCL6 binding was significantly enriched at the dynamic male-biased DHS (ES=2.0, p=E-14) but showed minimal enrichment at static male-biased DHS (ES=1.3, p=1.6E-04). BCL6 binding was also enriched at dynamic sex-independent DHS (Table S7B). In contrast, CUX2 was most highly enriched for binding in female liver at genomic sites that correspond to static male-biased DHS (ES=4.6, p=E-40). This binding of CUX2 may contribute to the greater closure of those DHS in female as compared to male liver (c.f., Fig. 3A, plot 6). Finally, male-biased STAT5 binding sites showed 2.5-fold greater enrichment at dynamic than at static male-biased DHS (ES=58, p=E-301 vs ES=24, p=E-190) (Table S7B), consistent with our findings, above, implicating pulsatile STAT5 binding in dynamic chromatin opening at those sites.

We reanalyzed published mouse liver ChIP-seq data for FOXA1 and FOXA2 [32], pioneer factors implicated in chromatin opening, [33, 34], to identify sex-dependent binding sites for each factor. Male-biased binding sites discovered for each FOXA factor showed strong enrichment for binding at male-biased DHS, with ∼2-fold higher enrichment seen at the static male-biased DHS (for FOXA1: ES=6.1 (static DHS) vs. ES=3.7 (dynamic DHS); for FOXA2: ES=44 (static DHS) vs. ES=21 (dynamic DHS) (Table S1B). Consistent with this, de novo motif discovery identified a Fox family factor as a close match for one of the top enriched motifs in the set of static but not in the set of dynamic male-biased DHS (Fig. S3C vs Fig. S3A, p-value 9.2E-06 vs. 2.0E-03).

Finally, female-biased FOXA2 binding sites, but not female-biased FOXA1 binding sites, showed strong enrichment for static female-biased DHS (ES=18, p=E-75; Table S7B). Taken together, these findings support the proposal that FOXA2, and to a lesser extent FOXA1, contribute to sex-dependent chromatin opening, in particular at static DHS.

Impact of hypophysectomy on liver chromatin accessibility

– Surgical removal of the pituitary gland (hypophysectomy) ablates pituitary GH secretion and thereby abolishes liver STAT5 activation [16], leading to widespread loss of sex-specific liver gene expression [6]. We hypothesized that by ablating GH-induced STAT5 activation and DNA binding, hypophysectomy will lead to closure of many of the open chromatin regions that regulate sex-specific gene expression. Furthermore, we proposed that restoration of a pulsatile GH signal, in the form of a single exogenous GH pulse given by i.p. injection, will reopen chromatin at many of those sites (Fig. 4A). To test this hypothesis, we used DNase-seq to compare the chromatin accessibility profiles of liver nuclei from hypox male and female mice to those of pituitary-intact control liver nuclei. Hypophysectomy induced changes in chromatin accessibility at several thousand sites, including large numbers of sex-independent DHS, with many more genomic regions undergoing chromatin closing than chromatin opening in male liver, but not in female liver (Fig. 4B, Table S4). Importantly, male-biased DHS that responded to hypophysectomy were almost exclusively closed in male liver following hypophysectomy (Fig. 4C, first two rows). A much smaller percentage of static female-biased DHS responded to hypophysectomy, and the responses observed generally involved chromatin opening in male liver and chromatin closing in female liver (Fig. 4C, Table S4G). Thus, sex-biased DHS are subject to both positive and negative pituitary hormone regulation and respond differently between the sexes.

Fig. 4.

DHS responses to pituitary ablation link to class I and class II sex-biased gene regulatory responses

- Hypophysectomy abolishes GH-regulated liver sex differences via two distinct mechanisms, which identify two classes of sex-biased genes [6]: hypophysectomy leads to down regulation of class I sex-biased genes in the sex where the gene is more highly expressed, and it leads to up regulation of class II sex-biased genes in the sex where the gene shows lower expression (Fig. S4, model). Building on this classification, we hypothesized that DHS that close following hypophysectomy are enhancers mapping to sex-specific genes subject to either positive regulation by GH (class I sex-biased genes) or negative regulation by GH (class II sex-biased genes). Supporting this proposal, the DHS that close in male liver following hypophysectomy showed strong, specific enrichment (ES=5.7, p=1.3E-69) for mapping to class I male-specific genes, which are down regulated in hypox male liver; whereas the DHS that close in hypox female liver showed strong, specific enrichment (ES=7.8, p=6.6E-46) for mapping to class I female-specific genes, which are down regulated in hypox female liver (Fig. 4D). In contrast, the DHS that open in hypox male liver were specifically enriched for mapping to class II female-specific genes (ES=3.6, p=1.8E-15), which are up regulated (de-repressed) in male liver following hypophysectomy; whereas the DHS that open in hypox female liver were specifically enriched albeit only moderately (ES=1.9, p=4.8E-04) for mapping to class II male-specific genes, which are up regulated in female liver following hypophysectomy (Fig. S4). The enrichment patterns exhibited by these four sets of hypophysectomy-responsive DHS mirror the corresponding sex-specific gene response patterns (Fig. 4D, bottom). Thus, all four DHS sets respond to hypophysectomy in a manner consistent with positively-acting regulatory elements linked to sex-biased genes. DHS that are linked to class I sex-biased genes, and whose chromatin closes following hypophysectomy, require pituitary hormone to maintain open chromatin; when pituitary hormones are ablated by hypophysectomy, chromatin closes and their class I sex-biased target genes are repressed. In contrast, DHS that are linked to class II sex-biased genes of the opposite sex bias, and whose chromatin opens following hypophysectomy, are kept in a closed chromatin state by pituitary hormone; when pituitary hormones are ablated, chromatin opens locally at those DHS and their class II sex-biased target genes are de-repressed: expression of female-biased class II genes increases in hypox male liver and expression of male-biased class II genes increases in hypox female liver.

A single exogenous GH pulse rapidly reopens chromatin at dynamic DHS in hypox male liver

– Next, we investigated whether the loss of pituitary GH accounts for the widespread chromatin closing seen in hypox male liver. Fig. 4B shows that exogenous GH treatment induces chromatin opening within 30 min at more than 3,000 DHS, including 71% of all dynamic DHS that closed in male liver following hypophysectomy (Fig. 4E). The fraction of reopening chromatin regions increased to 83% when considering the set of dynamic male-biased DHS that closed following hypophysectomy (Table S2A, column AK). Chromatin reopening was sustained for at least 90 min and then decreased substantially back toward baseline by 240 min (Fig. 4B; Fig. 5A, plot 1), at which time liver STAT5 signaling has terminated [16]. Importantly, the mean level of chromatin re-opening induced by an exogenous GH pulse at dynamic male-biased DHS was very similar to that of the same DHS set in STAT5-high male liver (Fig. 5A, plot 1 vs. Fig. 3A, plot 5; peak DNase-I activity value 1886 vs. 1848). Thus, the rapid chromatin opening induced by a single exogenous pulse of GH recapitulates the effects of an endogenous GH pulse in opening dynamic male-biased DHS.

Fig. 5.

Hypophysectomy led to closure of many fewer static male-biased DHS than dynamic male-biased DHS (Fig. 4C). Exogenous GH treatment partially reversed chromatin closing at static male-biased DHS within 30 min, with the effect persisting, even after 240 min (Fig. 5A, plot 2). This persistence is consistent with static male-biased DHS remaining constitutively open in male liver when liver STAT5 activity dissipates between endogenous plasma GH pulses. Static female-biased DHS also showed a decrease in mean chromatin opening following hypophysectomy, with further decreases in accessibility seen following GH pulse treatment (Fig. 5A, plot 3). Finally, in an important control, the mean chromatin accessibility at the set of 64,584 static sex-independent DHS was unchanged following hypophysectomy or GH pulse stimulation (Fig. 5A, plot 4).

Stratification of the full set of 2,729 male-biased DHS by the presence or absence of STAT5 binding, as determined by ChIP-seq, highlighted the STAT5 dependence of the rapid increases in chromatin accessibility stimulated by GH treatment (Fig. 5B, plot 1 vs plot 2). GH stimulated chromatin reopening at the static male-biased DHS was also primarily associated with the STAT5-bound DHS subset (Fig. S5, plot 2A vs 2B). In contrast, GH induced a much smaller increase in chromatin opening at the STAT5-bound subset of sex-independent DHS (Fig. 5B, plot 3 vs plot 1), where chromatin is already open and largely independent of plasma GH pulses (Fig. 3A, plot 3). GH decreased chromatin accessibility within 30 min at the 123 genomic sites with lower chromatin accessibility in STAT5-high liver as compared to STAT5-low liver (c.f., Fig. 2E), followed by a return to baseline by 240 min (Fig. 5C); this response pattern is consistent with the repression of chromatin accessibility at these sites by endogenous GH/STAT5 pulses (Fig. 2E). Finally, sex-independent DHS that do not bind STAT5 were unresponsive to both hypophysectomy and GH pulse treatment (Fig. 5B, plot 4).

A subset comprised of 354 dynamic, STAT5-high responsive DHS that close following hypophysectomy did not reopen even 240 min after GH pulse treatment (Fig. 4E, Fig. S6). These 354 DHS showed significantly lower differential chromatin accessibility between STAT5-high and STAT5-low livers than did the dynamic STAT5-high DHS subset that responded to an exogenous GH pulse (Fig. S6C, set 4 vs set 5). These 354 dynamic DHS may require multiple GH pulses to reopen or may be co-dependent on other pituitary-regulated hormones ablated by hypophysectomy.

Enrichment of sex-biased histone marks at static and dynamic male-biased DHS

– Our initial analyses revealed no major differences between dynamic and static male-biased DHS regarding enhancer vs insulator vs promoter classifications (Fig. 1A) or overall chromatin state distributions (Fig. 1D, Fig. S7). We therefore examined the two classes of male-biased DHS for differences in sex-biased histone marks (Table 1C). Male-biased enhancer marks (H3K27ac and H3K4me1) showed strong enrichment at both dynamic and static male-biased DHS when compared to a background set of 64,584 static, sex-independent DHS (Fig. 6A, Table S7B). In contrast, male-biased H3K36me3 marks, which are characteristic of transcribed regions, were enriched at static but not at dynamic male-biased DHS. In addition, static but not at dynamic male-biased DHS were significantly depleted of female-biased enhancer marks (H3K27ac and H3K4me1) (Fig. 6B). This finding is consistent with the above finding that static male-biased sites are in a comparatively closed state in female liver (Fig. 3, plot 2; model in Fig. 6C). Static female-biased DHS showed strong enrichments for female-biased enhancer marks (H3K27ac and H3K4me1), for female-biased H3K4me3 promoter marks, and for female-biased H3K36me3 transcribed region marks (Fig. 6B), a pattern that was very similar to the enrichments seen at static male-biased DHS (Fig. 6A).

Fig. 6.

Examination of two repressive histone marks, H3K27me3 and H3K9me3, revealed their use in a unique way in each sex to enforce sex differences in chromatin states at sex-biased DHS (Fig. 6C). Male-biased H3K27me3 marks were specifically associated with static female-biased DHS (Fig. 6A), consistent with our prior work showing deposition of these marks by the Ezh1/Ezh2 enzymatic component of PRC2 as a specific mechanism to repress many female-biased genes in male liver [35]. In contrast, female-biased H3K9me3 repressive marks were significantly enriched at both dynamic and static male-biased DHS (Fig. 6B); however, we did not find any corresponding enrichment of male-biased H3K9me3 marks at female-biased DHS. This novel finding suggests that H3K9me3 marks are specifically used to repress male-biased genes in female liver.

Chromatin state analysis at sex-biased DHS

– To further investigate chromatin state differences between dynamic and static DHS, we computed enrichment scores and depletion scores at each DHS set for each of 14 chromatin states. These 14 states were defined by a panel of activating and repressive chromatin marks and by the presence of open chromatin (Fig. 1C) and were determined separately for male and female mouse liver [18]. Enrichments were calculated using a genome-wide background set comprised of all 64,584 static sex-independent DHS (Fig. 7). In male liver, enhancer state E6 (characterized by the presence of DHS, H3K27ac and H3K4me1; Fig. 1C) was enriched at both dynamic and static male-biased DHS, and at dynamic sex-independent DHS. Promoter states E7 and E8 were significantly depleted, consistent with the relative paucity of promoter states in the full set of sex-biased DHS (Fig. 1A). We also observed strong depletion of the inactive chromatin states E1 (H3K27me3 marks) and E2 (marks not known) and of the transcribed state E14 (H3K36me3 marks), specifically from the dynamic but not static male-biased DHS in male liver (Fig. 7A). Finally, dynamic and static male-biased DHS were both strongly enriched for the inactive state E2 in female liver (Fig. 7B).

Fig. 7.

Dynamic male-biased DHS were enriched for the enhancer state E9 in female liver (ES=2.42, p=1.6E-69; Fig. 7B) but not male liver (Fig. 7A). The emission parameters of state E9 (Fig. 1C) are very similar to those of state E6, namely, it shows a high frequency of the activating chromatin marks H3K27ac and H3K4me1 but lacks the open chromatin (DHS) feature characteristic of state E6. Thus, in female liver, dynamic male-biased DHS contain histone marks that typify an active enhancer but are inactive due to the comparatively closed state of their chromatin. Indeed, the extent of chromatin opening at these DHS in female liver is equivalent to the level of chromatin opening seen at the same set of sites in male liver between plasma GH/liver STAT5 pulses (Fig. 3A, plot 5; model in Fig. 6C). Dynamic male-biased DHS were also enriched for chromatin state E10 in female liver (ES=2.89, p=1.9E-18) but not male liver. State E10 shows a high frequency of H3K27ac marks, but not H3K4me1 marks, and lacks DHS. This pattern of chromatin state E9 and E10 enrichment at dynamic male-biased DHS in female liver, combined with the strong enrichment of H3K9me3 repressive marks (Fig. 6B), indicates that the genomic regions encompassing dynamic male-biased DHS are in a bivalent chromatin state in female liver (Fig. 6C). Conceivably, such a bivalent state may facilitate chromatin opening under certain biological conditions, e.g., chemical exposure or biological stress.

Static female-biased DHS showed little or no enrichment of specific chromatin states in female liver when compared to the genome-wide background set of static sex-independent DHS (Fig. 7B). This supports a model whereby their female bias in accessibility is largely due to active suppression in male liver (Fig. 6C). Indeed, the genomic regions encompassing static female-biased DHS were strongly enriched in male liver for inactive chromatin states E1 (characterized by repressive mark H3K27me3), E2 and E4 (Fig. 7A). Static female-biased DHS showed weak enrichment for enhancer states E10 and E11 in male liver, which are respectively characterized by the active enhancer marks H3K27ac and H3K4me1 but devoid of DHS. This is analogous to the finding, above, that male-biased DHS are enriched for chromatin states replete with active histone marks but deficient in DHS in female liver. Static female-biased DHS also showed modest enrichment in male liver for bivalent state E12, which is characterized by a mixture of activating histone marks and the presence of repressive H3K27me3 marks, consistent with the strong enrichment of the latter marks seen in Fig. 6A, above. Finally, promoter-like states E7 and E8 were significantly depleted from static female-biased DHS (Fig. 7), consistent with the female-biased DHS largely being gene distal sex-biased enhancers.

Distinct gene targets and enriched biological processes of dynamic and static sex-biased DHS

– Dynamic and static male-biased DHS both showed strong enrichment for mapping to male-biased gene targets (8-fold and 11.7-fold enrichments, respectively), as did static female-biased DHS for female-biased gene targets (12.5-fold enrichment) when DHS were mapped to the single nearest transcription start site in the same TAD (Table S2B). Further, mapping of DHS to putative target genes using GREAT, which typically maps each DHS to two genes, revealed how many male-biased genes may be regulated by multiple male-biased DHS. In some cases, all of the associated male-biased DHS were static male-biased DHS (e.g., Cyp4a12a, with 7 static male-biased DHS), while in other cases they were a mixture of dynamic and static male-biased DHS (e.g., Cyp7b1, with 7 dynamic and 8 static male-biased DHS) (Table S9A). Top female-biased genes enriched for nearby static female-biased DHS included Cux2 and three Cyp3a genes, with 10-15 female-biased DHS each (Table S10C). Sex-independent genes that are well-established direct targets of STAT5 [36] include Igf1, a target of 12 dynamic sex-independent DHS and 3 dynamic male-biased DHS, and Socs2, a target of 8 dynamic sex-independent DHS (Table S10D, Table S10A). Other analyses revealed that all three sets of sex-biased DHS were significantly enriched for mapping to hypophysectomy class I-responsive sex-biased genes as compared to hypophysectomy class II-responsive sex-biased genes (Table S2B). This finding is consistent with the sex-biased DHS sets being enriched for positively-acting enhancers that close following hypophysectomy, which leads to decreased expression of their class-I (i.e., hypophysectomy repressed) gene targets.

The gene targets of each sex-biased DHS set were enriched for distinct but partially overlapping Gene Ontology (GO) Biological Processes, as determined by GREAT analysis. Top enriched GO terms common to both dynamic and static male-biased DHS included lipid metabolic process and steroid metabolic process, consistent with the major role that GH plays in the male-prevalence of fatty liver development and liver metabolic disease [3–5]. Unique enriched terms for dynamic male-biased DHS included cell adhesion, gland development and hepaticobiliary development, whereas gene targets of static male-biased DHS were uniquely enriched for cellular response to glucocorticoid stimulus and EGF receptor signaling, among others (Fig. S8, Table S10, bottom).

Static female-biased DHS gene targets were uniquely enriched for long chain fatty acid metabolic pathway and related terms. Finally, the set of dynamic sex-independent DHS mapped to gene targets enriched for terms related to glucose transport and metabolism, IGF receptor signaling and fibroblast proliferation (Fig. S8). This finding is consistent with the widespread metabolic effects that GH has in both sexes, including complex effects on glucose uptake and glucose oxidation, and on hepatic gluconeogenesis and glycogenolysis [37].

Discussion

Sex differences in chromatin accessibility are a central epigenetic feature that enables regulatory proteins, such as the GH-activated transcription factor STAT5, to bind chromatin and regulate hepatocyte gene transcription in a sex-specific manner. However, the underlying mechanisms controlling sex differences in liver chromatin accessibility, which occur at more than 4,000 distinct sites across the mouse genome, are poorly understood. GH, working through its sex-dependent temporal patterns of secretion by the pituitary gland, is the major hormonal factor controlling STAT5-dependent sex differences in liver gene transcription, but little is known about the potential of GH to directly regulate sex differences in the epigenetic landscape required for sex-dependent transcriptional outputs. To address this question, we elucidated the impact of male plasma GH pulses on global patterns of liver chromatin accessibility by analyzing a population of 21 individual male mouse livers euthanized at either a peak or a trough of plasma GH pulse-stimulated hepatic STAT5 DNA-binding activity. Our findings establish that the naturally occurring, endogenous plasma GH pulses characteristic of males induce dynamic cycles of chromatin opening and closing at several thousand DNase-I hypersensitive sites (DHS) in male mouse liver and that these events comprise one of two major mechanisms regulating the male bias in liver chromatin accessibility. Analysis of sex-dependent transcription factor binding patterns and chromatin states elucidated key features distinguishing this dynamic mechanism of male-biased enhancer activation from that of static, GH pulse-unresponsive male-biased DHS and from female-biased DHS (see model, Fig. 6C). GH thus acts at three distinct steps to regulate sex-dependent hepatocyte gene expression, all three involving the GH-activated transcription factor STAT5 (Fig. 8): GH pulse-induced chromatin opening at dynamic male-biased DHS driven by pulsatile GH activation of STAT5, as discussed further, below; direct transcriptional activation of sex-specific genes by GH-activated STAT5; and GH/STAT5-dependent transcriptional regulation of downstream repressors, such as the female-specific CUX2, which binds to male-biased enhancers in female liver to reinforce sex differences in gene transcription.

Fig. 8.

GH pulse activation of dynamic male-biased DHS

- GH-induced chromatin opening in male mouse liver is shown to be directly controlled by the endogenous male rhythm of plasma GH stimulation of hepatocytes, the major liver cell type contributing to sex-biased liver gene expression [38]. The GH pulse-regulated chromatin regions identified here responded to changes in plasma GH levels in a rapid, dynamic manner, with extensive GH-induced chromatin opening occurring within 30 min, as seen when hypophysectomized male mice were given a physiological replacement dose of GH pulse by intraperitoneal injection (Fig. 4A). This time course is consistent with the rapid activation of GH receptor signaling to STAT5, which occurs within 5 min in cell culture [39] and within 15 min in vivo in a hypophysectomized rat model [40]. Importantly, 85% of the dynamic male-biased DHS identified here (710 of 834 dynamic male-biased DHS) were bound by STAT5, which appears to be a key driver of these GH pulse-induced chromatin opening events. Chromatin closing followed the termination of nuclear STAT5 signaling, which is complete within 4 h [16] and occurs in sufficient time to reset the GH receptor-JAK2 signaling complex and resensitize hepatocytes before the next plasma GH pulse [39] (Fig. 2A). STAT5 binding and chromatin opening at dynamic male-biased DHS were both significantly lower in livers of female mice, where plasma GH levels and liver STAT5 activity are persistent (near-continuous) yet ineffective at maintaining an open chromatin state. Indeed, the extent of chromatin opening at these male-biased DHS in female liver is very similar to the level in livers of male mice euthanized between plasma GH pulses (STAT5-low male livers; Fig. 3A, panel 5). Pulsatile chromatin opening stimulated by endogenous plasma GH pulses is thus a unique mechanism for establishing and maintaining male-biased chromatin accessibility and transcription factor binding at 834 male-biased DHS, corresponding to 31% of all male-biased DHS in the liver.

GH pulse-responsive DHS were discovered by comparing global chromatin accessibility patterns in a set of 8 livers collected from mice euthanized at a peak of GH-activated STAT5 DNA-binding activity (STAT5-high livers) to those of 10 other livers from mice euthanized between pulses of STAT5 DNA-binding activity, i.e., when liver STAT5 activity is very low or undetectable (STAT5-low livers). Three other male livers gave DNase-seq profiles inconsistent with their EMSA-determined STAT5 DNA-binding activity. These outlier livers may have come from mice euthanized just after the onset of a STAT5 pulse [27], when more time is needed to induce chromatin opening (2 outliers), or shortly after STAT5 is deactivated by tyrosine dephosphorylation [41] but prior to the reversal of chromatin opening (1 outlier). Importantly, the top 200 genomic regions showing differential chromatin opening between STAT5-high and STAT5-low livers also separated a group of nine C57Bl/6 male mouse liver DNase-seq samples from the ENCODE consortium into two distinct classes, corresponding to STAT5-high activity livers (n=4) and STAT5-low activity livers (n=5), respectively. Thus, the dynamic pattern of chromatin opening and closing at these genomic regions is robust across studies and mouse strains.

STAT5-high male livers showed the highest levels of chromatin accessibility among the male-biased DHS that bind STAT5. Similarly, a higher level of chromatin opening was seen at sex-independent DHS that bind STAT5 compared to those that do not bind STAT5 (Fig. 3A, panel 3 vs 4). There is thus a close linkage between STAT5 binding and chromatin opening. This finding is consistent with the proposal that STAT5 acts as a pioneer factor to enable chromatin opening, as was reported for STAT5 action at the Il9 gene locus in Th9 T-cells [42]. It should be noted, however, that pulsatile chromatin opening also occurred at the small subset (15%) of dynamic male-biased DHS that did not bind STAT5, suggesting a role for other GH receptor signaling pathways [43] in male-biased chromatin opening at those sites. Furthermore, pulsatile STAT5 activation and pulsatile STAT5 DNA binding alone are not sufficient to insure pulsatile chromatin opening, insofar as many sex-independent DHS do bind STAT5, yet are constitutively open in both male and female liver (Fig. S2C, plot 5A). Finally, plasma GH pulse stimulation decreased chromatin accessibility at a small number of liver DHS, most of which were sex-independent. The mechanism for this decrease and its physiological significance are unknown.

Regulation of static male-biased DHS

- We identified a second, but less well-defined mechanism that controls the male bias in chromatin accessibility at a distinct DHS set, comprised of 1,895 static male-biased DHS. These DHS represent 69% of all male-biased DHS and mapped to a set of target genes and enriched biological processes distinct from but overlapping those of the dynamic male-biased DHS. Static male-biased DHS are constitutively open in male liver, with chromatin accessibility largely unchanged across the peaks and valleys of GH-induced liver STAT5 activity. Moreover, in contrast to dynamic male-biased DHS, the male bias in chromatin accessibility at static male-biased DHS is primarily a reflection of their being in a closed state in female liver (Fig. 6C). Static male-biased DHS are relatively deficient in STAT5 binding but showed significant enrichment for binding CUX2, a female-specific repressor protein [30], which we infer contributes to the closure of static male-biased DHS in female liver through its established transcriptional repressor activity [44]. Female-biased DHS showed even greater depletion of STAT5 binding (comparable to genome-wide background of sex-independent DHS), and correspondingly, the vast majority (99%) are static i.e., unresponsive to plasma pulse stimulation in male liver.

Our findings support the following proposal to help explain the widespread closure of male-biased liver DHS, most notably dynamic male-biased DHS (95% closing; Table S2C) when male mice are given GH as a continuous infusion, which mimics the female plasma GH pattern and substantially feminizes liver gene expression [17, 45] (Fig. 6C). First, dynamic male-biased DHS close due to loss of the pulsatile plasma GH pattern and its intrinsic ability to stimulate chromatin opening. Second, static male-biased DHS close due to the strong induction of CUX2 [17], whose binding is enriched at those sites and is associated with chromatin closing, as noted above. The pioneer factors FOXA1 and FOXA2 showed strongest enrichment of their male-biased binding sites at static male-biased DHS and may help maintain these DHS in their constitutively open, male-biased chromatin state. Correspondingly, static female-biased DHS showed strongest enrichment for female-biased FOXA2 binding sites, which likely contributes to their greater accessibility in female compared to male liver (Fig. 6C).

Novel insights from chromatin state analysis

- Dynamic and static male-biased DHS are both largely (∼95%) comprised of promoter-distal enhancer DHS, as indicated by their flanking histone modifications. This supports our prior finding that sex-dependent intra-TAD DNA looping mechanisms are common in mouse liver [24] and indicates such looping mechanisms likely bring both subsets of male-biased DHS into closer proximity with their sex-biased promoter targets. Dynamic male-biased DHS were enriched for the active enhancer state E6 in male liver but were enriched for enhancer states E9 and E10 in female liver. States E9 and E10 are characterized by a high frequency of the activating chromatin marks H3K27ac and H3K4me1 (E9) or H3K27ac alone (E10), similar to state E6, but unlike E6, they are both deficient in open chromatin (DHS) (Fig. 1C). Thus, in female liver, dynamic male-biased DHS regions contain active histone marks but are in a comparatively closed chromatin state. This may in part be due to their enrichment for H3K9me3 (repressive) histone marks in female liver, which gives them bivalent character (Fig. 6C). This bivalent chromatin state may give dynamic male-biased DHS the potential for chromatin opening and activation of their latent enhancer activity in female liver, e.g., in response to chemical exposure or biological stress [46]. In contrast, static male-biased DHS are depleted of female-biased enhancer marks (H3K27ac and H3K4me1) and are in an inactive chromatin state in female liver (Fig. 6C). Of note, the bivalent chromatin state of dynamic male-biased DHS in female liver is distinct from the poised state these genomic regions apparently adopt in male liver between plasma GH pulses (i.e., in STAT5-low male liver) and in hypophysectomized male liver, where a single GH pulse is all that is required to induce rapid chromatin opening, even after several weeks of pituitary hormone ablation.

Static female-biased DHS are deficient in DHS and are strongly enriched for the repressive histone mark H3K27me3 in male liver, which is thought to contribute directly to their closed, inactive chromatin state in male liver (Fig. 6C). Prior studies have revealed heterogeneity in these female-biased DHS regarding the time-dependent loss of H3K27me3 marks and chromatin opening in livers of male mice exposed to GH as a continuous infusion [17], and in the extent of female-biased target gene de-repression following the loss of H3K27me3 marks in livers of male mice deficient in the H3K27-trimethylase complex PRC2 [35]. Novel insight into the fundamental underlying epigenetic mechanisms of sex-biased gene regulation comes from our striking finding that H3K27me3 and H3K9me3, two widely used repressive histone marks, are used in a unique way in each sex to enforce sex differences in chromatin states at sex-biased DHS. Thus, whereas male-biased H3K27me3 repressive marks were highly enriched at, and specifically associated with, static female-biased DHS, consistent with our prior work [18, 35], we show here that female-biased H3K9me3 repressive marks are highly enriched at both dynamic and static male-biased DHS, but without a corresponding enrichment of male-biased H3K9me3 marks at static female-biased DHS. This novel discovery supports the proposal that sex-specific gene repression is mediated by two distinct mechanisms, with H3K9me3 marks specifically used to repress male-biased genes in female liver, and H3K27me3 marks specifically used to repress female-biased genes in male liver. Further study will be required to elucidate the underlying mechanisms whereby these repressive histone marks are respectively deposited in a sex-dependent and site-specific manner and the role of GH in regulating these molecular events.

Regulatory elements associated with class I and class II sex-biased genes

- DNase-seq profiling identified more than 5,000 liver DHS that open or close following hypophysectomy. Strikingly, a single physiological replacement dose of GH restored, within 30 min, liver STAT5 activity and chromatin accessibility at 83% of dynamic male-biased DHS that closed following hypophysectomy. Remarkably, these chromatin sites are maintained in a primed (poised) state in hypox male liver over a period of several weeks, despite the prolonged deficiency of GH and other pituitary hormones and the dysregulation of several thousand liver-expressed genes [6, 16]. DHS that close in male liver following hypophysectomy were enriched for proximity to class I male-specific genes, which are down regulated in hypox male liver, whereas DHS that close in hypox female liver showed the strongest enrichment for class I female-specific gene targets, which are down regulated in hypox female liver (model, Fig. S4). Further, DHS that open in hypox male liver were enriched for proximity to class II female-specific genes and DHS that open in female liver were enriched for proximity to class II male-specific genes. These enrichments are consistent with the definition of class II sex-biased genes, i.e., genes that are up regulated by hypophysectomy in the sex where they show lower expression in intact mice. Together, these findings lend strong support for the functional role of each class of sex-biased and pituitary hormone-dependent DHS as regulatory elements with intrinsic, positively acting enhancer potential, keeping gene expression on in intact liver by maintaining open chromatin, e.g., in the case of class I male-biased DHS and class I male-biased genes in male liver, or keeping gene expression off by maintaining a closed chromatin state, e.g., in the case of class II male-biased DHS and class II female-biased genes in male liver. Further study will be required to elucidate the underlying regulatory factors and molecular mechanisms governing these gene regulatory circuits.

Methods

Animal treatments and DNase-seq analysis

- All mouse work was carried out in accordance with ARRIVE Essential guidelines 2.0 [47] for study design, sample size, randomization, experimental animals and procedures, and statistical methods, and with approval of the Boston University Institutional Animal Care and Use Committee. Male and female CD1 mice (ICR strain, strain code # 022), 7-weeks old and purchased from Charles River Laboratories (Wilmington, MA), were kept on a 12-hour light/dark cycle with food and water without restriction. Livers were collected from individual untreated mice between 8 and 10 weeks of age. Where indicated, male and female mice were hypophysectomized (hypox) by the supplier at ∼7-8 weeks of age.

Completeness of hypophysectomy was verified by the absence of weight gain over a 2-3 week period and by the lack of Mup protein in urine samples (SDS-PAGE analysis) [16]. Hypox male mice were treated with recombinant rat GH (purchased from Prof. Arieh Gertler, Protein Laboratories Rehovot Ltd, Rehovot, Israel), given as a single intraperitoneal injection at 125 ng of GH per gram body weight [6], or vehicle (control), and were euthanized by cervical dislocation 30, 90 or 240 min later [16]. A protein extract prepared from a small piece of each liver to be used for DNase-seq analysis (see below) was assayed for liver STAT5 DNA-binding activity by EMSA, as described [16]. Results of this assay allowed us to classify each individual male mouse as STAT5-high activity (n=10) or STAT5-low activity (n=11) at the time of euthanasia and liver collection. EMSA was also used to verify the ablation of liver STAT5 activity following hypophysectomy and the effectiveness of exogenous GH administration at restoring endogenous STAT5-high liver levels of STAT5 EMSA activity within 30 min.

DNase-seq libraries, Illumina sequencing, and data processing

- Livers used for these analyses were obtained from intact male and female mice, from hypox male and hypox female mice, and from hypox male mice given a single injection of GH and euthanized 30, 90 or 240 min later (n=6-12 livers per group). Nuclei were purified from individual fresh livers as described [19] and stored at -80 °C. Nuclei were subsequently treated with DNase I to release DNA fragments from hypersensitive genomic regions to identify liver DNase hypersensitive sites (DHS) [16]. The DNA fragments released from each liver were combined (3 livers per pool) to give n=2 to n=4 independent pooled samples, which were used to prepare DNase-seq libraries for each mouse group (Table S1). In the case of STAT5-high and STAT5-low liver DHS analysis, however, DNase-seq libraries were prepared directly from DNase-I digested fragments released from each individual liver without pooling, as described below. DNase-seq libraries were prepared using the NEBNext Ultra Library Prep Kit for Illumina (New England Biolabs). Illumina sequencing was performed on an Illumina HiSeq instrument to a depth of 12 to 33 million mapped 40-50 bp single-end sequence reads per sample. Detailed sequencing statistics are shown in Table S1. Raw and processed data are available at www.ncbi.nlm.nih.gov/gds under accession numbers GSE131848 and GSE131852 (SuperSeries GSE131853).

Sequencing data was analyzed using a custom in-house DNase-seq pipeline [46]. Briefly, the pipeline processes raw FASTQ files and outputs various control metrics, including FASTQC reports (FastQC v0.11.7), confirmation of read length, verification of the absence of read strand bias, and identification of contaminating adaptor sequences (Trim_galore v0.4.2). Reads were mapped to mouse genome mm9 using Bowtie2 (v2.2.6) [48]. Regions of DNase hypersensitivity (DHS) were discovered as peaks identified by MACS2 (v2.1.0.20150731) [49] using the option (-nomodel –shift -100 –extsize 200) to inhibit read shifting, and the option (-keep-dup) to retain all reads that contribute to the peak signal. Peaks were discovered for each individual DNase-seq library then filtered to remove peaks that overlap ENCODE blacklisted regions [50] as well as peaks comprised of 5 or more identical reads that do not overlap any other read (“straight peaks”). Nine additional individual male mouse liver DNase-seq datasets generated by the ENCODE consortium (FASTQ files downloaded at http://genome.ucsc.edu/cgi-bin/hgFileUi?db=mm9&g=wgEncodeUwDnase; GEO accession: GSM1014195, liver replicates #5 to #13) were analyzed using the same DNase-seq pipeline (cell line, liver; strain, C57BL/6; age, adult 8 weeks; sex, male). These DNase-seq datasets were single end sequencing reads, 36 bp in length, with a sequencing depth of 19 to 37 million total mapped reads per sample.

Liver DHS classification and DHS gene target expression data

- A set of 72,862 mouse liver DHS regions previously identified by DNase-seq analysis of male and female mouse liver [19] was filtered to remove 515 DHS regions that overlapped ENCODE blacklisted regions. 2,136 DHS regions that could not be classified according to a 5-class DHS model defined previously [20] were also removed to obtain a standard reference set of 70,211 mouse liver DHS. Each DHS was designated as a promoter, weak promoter, enhancer, weak enhancer, or insulator DHS based on the 5-class model, which is primarily based on ChIP-seq signals for the H3 histone marks H3K4me1and H3K4me3, combined with CTCF ChIP-seq binding in adult male mouse liver [20]. Weak enhancers were identified by their low levels of H3K27ac ChIP-seq signal compared to DHS classified as enhancers, combined with a distance from RefSeq gene transcription start sites inconsistent with a weak promoter designation [20]. Further, individual DHS were designated male-biased (n=2,729), female-biased (n=1,366), or sex-independent (n= 66,116) based on significant differences in chromatin accessibility between male and female mouse liver [19]. Except where noted, each DHS was assigned a single putative gene target (RefSeq gene or multi-exonic lncRNA gene) corresponding to the closest transcription start site within the same TAD [20]. Table S2A lists DHS target genes, their overlap with STAT5 binding sites determined by ChIP-seq [15], their DNase-seq activity (chromatin accessibility) in male and female mouse liver, and their responses to hypophysectomy and to hypophysectomy + GH treatment. Where indicated, DHS target genes were defined using GREAT (v.4.0.4) [51, 52] (http://great.stanford.edu/public/html/index.php) and the default settings for the Basal + extension gene association rule, which typically assigns two genes to each DHS (see output in Table S9 and Table S10, below).

Differential DHS between STAT5-high and STAT5-low activity livers

– Nuclei purified from each individual STAT5-high activity (n=10) and STAT5 low-activity mouse liver (n=11) (Table S1), classified based on EMSA as described above, were analyzed by DNase-seq to discover genomic regions (i.e., DHS) that responded dynamically to endogenous plasma GH pulses and the associated changes in liver STAT5 DNA-binding activity, as follows. diffReps analysis [28] was used to discover genomic sites that were more open or were more closed (DNase-seq normalized intensity |fold-change| > 2 and FDR < 0.05 [Benjamini-Hochberg adjusted p-value]) for the set of EMSA-identified STAT5-high activity livers as compared to the set of EMSA-identified STAT5-low activity livers. The diffReps nucleosome option (200 bp window size) was used and the option (-frag) was set to zero for all comparisons. The differential sites that diffReps identified as significant were further analyzed by two methods, principal component analysis and boxplot analysis, to determine the distributions of normalized sequence read counts (reads per kilobase per million mapped reads) for each individual DNase-seq library. Three outlier DNase-seq samples were thus identified: two from livers designated STAT5-high activity based on EMSA (samples G74A_M1 and G74A_M2), which gave DNase-seq read count distributions (top 200 and top 600 diffReps-identified differential sites) more similar to STAT5-low activity livers; and one from a STAT5-low EMSA activity liver (sample G92_M5), which gave a DNase-seq read count distribution more similar to STAT5-high activity livers (Fig. 2C and Fig. 2D, below). These three outlier liver DNase-seq libraries were excluded from all downstream analysis. Next, we implemented diffReps analysis using the same cutoffs described above to compare the DNase-seq activity profiles of the remaining EMSA-identified STAT5-high activity (n=8) and STAT5-low activity (n=10) DNase-seq samples. The resultant set of diffReps differential sites was overlapped with a MACS2 DHS peak union list, which was obtained by merging the MACS2 DHS peak calls from all 18 male mouse liver DNase-seq samples using the BEDTools “merge” command. This overlap analysis yielded a final set of 70,767 MACS2 DHS peak union sites (see Table S4A), of which 2,832 MACS2-identified DHS that were more open (more accessible chromatin state) and 123 DHS were in a more closed chromatin state based on their overlap with the diffReps-identified STAT5-high vs. STAT5-low livers differential regions (see Fig. 2E, below). The other 67,812 DHS regions were designated static DHS, as they did not overlap a diffReps-identified STAT5-high vs. STAT5-low differential region.

Comparison of these 70,767 MACS2-identified DHS peak union sites with our standard reference set of 70,211 mouse liver DHS identified n=2,373 reference set DHS that overlapped a STAT5-high differential site (see Fig. 2F, below), n=98 reference set DHS that overlapped a STAT5-low differential site, and n=45,754 liver DHS that overlapped a static DHS site (Table S2A, columns H and I). The remaining 21,986 standard reference set DHS did not overlap the 70,767 DHS peak union list and were also labeled static DHS. We applied the designation of sex bias determined previously [19], namely, male-biased, female-biased, or sex-independent, to each DHS that overlapped the reference set 70,211 liver DHS (Table S2A, column F). The standard reference set liver DHS were further designated dynamic DHS if they overlapped a STAT5-high > STAT5-low differential DHS (i.e., a GH/STAT5 pulse-opened DHS). Liver DHS not identified as dynamic (including the 98 STAT5-low > STAT5-high DHS peak union sites that overlap a standard reference set DHS) were designated static with respect to chromatin accessibility changes induced by endogenous GH-induced STAT5 pulses. Sex-specific DHS were thus designated dynamic male-biased DHS, static male-biased DHS, dynamic female-biased DHS, or static female-biased DHS with respect to GH/STAT5-induced chromatin opening (Table S2A, column I). Sex-independent DHS were similarly designated dynamic or static sex-independent DHS (Fig. 2F, below).

DNase-I cut site aggregate plots

- Normalized DNase-I cut site aggregate plots were generated using an input DNase-seq dataset and the set of input genomic regions (DHS sequences) to be analyzed by sequence read counting. First, FASTQ files from DNase-seq biological replicates were concatenated to obtain a single combined replicates file for each treatment group. For example, for male liver, we generated a single combined STAT5-high DNase-seq FASTQ file by merging the n=8 biological replicate STAT5-high liver FASTQ files, and separately, we generated a single combined STAT5-low sample by merging the n=10 biological replicate STAT5-low DNase-seq FASTQ files (Table S1). Combined replicate FASTQ files were generated for intact female, hypox female, hypox male, and each of the hypox male + GH treatment time point DNase-seq samples in the same manner. Second, for each combined replicate file, a BED file comprised of positive and negative strand reads was processed to determine the DNase-I cut site, which corresponds to the 5’-end of each sequence read. Third, the set of input genomic regions was processed to generate a list of 2 kb regions centered at the midpoint of each DHS. The BEDTools “coverage” command using the (-d) option was used to count the number of DNase-I cut sites at each nt position of the 2 kb midpoint-centered regions, thus producing a read count matrix composed of 2,000 read counts for each input genomic region. Fourth, a custom R script was used to load the read count matrix, calculate the sum of read counts at each nt position, normalize the raw read counts by the number of reads in the subset of 70,211 standard reference DHS to be analyzed, normalize by the number of input genomic regions, and then generate a plot of the DNase-I signal across the input genomic regions. The resultant normalized DNase-I cutting profiles were then smoothed using a LOWESS smoother as implemented in R (package: gplots v3.0.1.1). An offset was applied to the profile by subtracting an average read count (calculated from the first 200 nt positions) from the normalized read count intensity to standardize the baseline of each profile. Profiles for other, related DNase-seq datasets were processed in the same manner and were plotted on the same set of axes to enable direct comparisons across all such plots in Figs. 2, 3, 5 and Figs. S1, S4 and S5.

Impact of hypophysectomy and GH injection on chromatin opening and closing

- MACS2 was used to discover differential DHS peaks in DNase-seq samples prepared from intact male and intact female mouse liver, hypox male and hypox female mouse liver, and from livers of hypox male mice given a single injection of GH and euthanized either 30, 90, or 240 min later (Table S1). The effects of hypophysectomy on chromatin accessibility were determined by comparing hypox liver DNase-seq samples to the corresponding samples prepared from intact liver (control) samples in male liver and separately in female liver. The effects of a single injection of GH were determined by comparing DNase-seq samples from livers of hypox male mice treated with GH to the untreated hypox male liver controls at each time point. For each comparison, a DHS peak union list was generated by merging the MACS2 DHS peak calls from all individual biological replicates for the corresponding treatment group. For example, a single peak union list for the hypox male compared to intact male liver samples was generated by merging all the MACS2 peaks from each of the n=8 individual male mouse liver DNase-seq samples (n=6 intact male control livers and n=2 hypox male livers; Table S1). Genomic regions that showed significantly differential DNase-seq signal between hypox and intact control livers, or between hypox male + GH and hypox male livers were discovered separately for each comparison using diffReps, using the nucleosome option (200 bp window). Significance was based on |fold-change| > 2 and FDR < 0.05 (Benjamini-Hochberg adjusted p-value) for diffReps-normalized signal intensity values. The diffReps-identified differential sites were then filtered by their overlap with the peak union list for all 8 samples to obtain the sets of differential DHS (e.g., 2,142 DHS that open and 4,856 DHS that close in male liver following hypophysectomy; Table S4, Summary sheet). MACS2 peak union DHS that did not overlap a diffReps region were annotated as static DHS peaks (50,055 sites). Concatenation of the diffReps-identified differential DHS with this set of static DHS yielded the full set of 57,053 DHS peak union sites for this dataset. Each of the 70,211 standard reference set liver DHS (see above) was then labeled based on whether it overlapped a diffReps-identified differential DHS that opened following hypophysectomy, a diffReps differential DHS that closed following hypophysectomy, a static DHS, or “none”, for those reference set DHS that did not overlap any of the set of 57,053 hypophysectomy study DHS peak union sites. Corresponding analyses were performed for the intact and hypox female liver samples, and for the hypox male + GH time-course comparisons mentioned above. Table S4 summarizes the numbers of hypophysectomy and hypophysectomy + GH responsive differential DHS for each of these comparisons, and the subsets that overlap the reference set of 70,211 liver DHS.

Sex-specific genes and hypophysectomy response classification

- RNA-seq data from intact male and intact female mouse liver was previously used to identify sex-specific genes based on a gene list comprised of 24,197 RefSeq genes [16]. This list was expanded to include both RefSeq and multi-exonic sex-specific lncRNA genes, as follows. Sex-specific RefSeq genes were identified from polyA-selected total liver RNA-seq samples from intact male and intact female liver using the “genebody” method of sequence read counting [16] at a threshold of |fold-change| > 1.5, adjusted p-value < 0.05 and FPKM > 0.25 for the sex with a greater signal intensity, which yielded 387 male-specific and 517 female-specific RefSeq genes. Sex-specific multi-exonic lncRNA genes were identified from polyA-selected liver nuclear RNA-seq samples in male and female mouse liver using the “ExonCollapsed” read counting method [53] at a threshold of |fold-change| > 2, adjusted p-value < 0.05 and FPKM > 0.25 for the sex with a higher signal intensity, which yielded 121 male-specific and 102 female-specific multi-exonic lncRNA genes, for a total of 508 male-specific and 619 female-specific genes (Table S3). A stringent set of sex-independent genes was defined by FPKM > 1 in both male and female liver, |fold-change| < 1.2, and adjusted p-value > 0.1, which yielded a total of 7,253 stringently sex-independent genes, including lncRNA genes (Table S3, column N). Sex-biased genes responsive to hypophysectomy were defined by |fold-change| > 2 and adjusted p-value < 0.05 for hypox mouse liver versus intact mouse liver, determined separately for male liver and for female liver. Using these thresholds, RNA-seq data from intact and hypox male and female liver samples were used to group sex-specific genes into classes I and II and their corresponding subclasses (IA, IB, IC, IIA, and IIB) [6, 16]. Class I sex-biased genes are those sex-biased that are down regulated by hypophysectomy in the sex where they show higher expression in intact mice. Class II sex-biased genes are those that are up regulated by hypophysectomy in the sex where they show lower expression in intact mice. Subclasses A, B, and C indicate the response to hypophysectomy in the dominant sex (class II genes) or in the opposite sex (class I genes), as defined in Table 3 of [6]. The number of sex-biased genes in each hypophysectomy response class is shown in Table S3.

STAT5 binding and motif analysis

- ChIP-seq analysis of STAT5 binding in male and female mouse liver identified 15,094 merged peaks comprised of male-enriched and female-enriched, and male-female common STAT5 binding sites [15]. A small subset of these STAT5 binding sites, 75 peaks, overlapped ENCODE blacklisted regions and was excluded from downstream analysis. BEDTools overlap analysis of these STAT5 binding sites allowed us to designate each of the 70,211 reference DHS as “STAT5-bound”, if a STAT5 binding site was present, or “not bound”. Table S2A provides full details, including the STAT5 binding site, the sex-specificity of STAT5 binding (male-enriched, female-enriched, or common to both sexes), the normalized ChIP-seq read counts for STAT5 binding in STAT5-high activity male livers (MH) and in STAT5-high activity female livers (FH), and the average of these two sets of sequence read counts. The STAT5B motif M00459 from the TRANSFAC motif database (Release 2011.1) [54] was used to determine the frequency of STAT5 motif occurrence in each set of DHS sequences. Motifs found in DHS sequences were identified using FIMO (v4.12.0) [55] using the option (--thres 0.0005) to improve detection of short length motifs. The number of STAT5B motif occurrences in each of the 70,211 reference set DHS sequences is shown in Table S2A.

Enrichment analysis

- For all enrichment calculations described below, the significance of the enrichment was determined by a Benjamini and Hochberg adjusted Fisher’s Exact test p-value as implemented in R. Enrichments with adjusted p-value < 1E-03 were considered statistically significant. Enrichments of sex-biased DHS for being enhancer, promoter, or insulator regions (e.g. male-biased DHS for being enhancers compared to the background set of sex-independent DHS) were calculated, as follows: Enrichment score = (ratio A) / (ratio B), where: ratio A = the number of sex-biased DHS that are enhancers, divided by the number of sex-biased DHS that are not enhancers; and ratio B = the number of sex-independent DHS that are enhancers, divided by the number of sex-independent DHS that are not enhancers. For example, 2,551 male-biased DHS were classified as enhancers, and 178 other male-biased DHS were not enhancers (2,551/178=14.3), whereas 43,591 sex-independent DHS were classified as enhancers, and 22,525 sex-independent DHS were not enhancers (43,591/22,525=1.93), which gives an enrichment score A/B=14.3/1.93=7.41. The set of 66,116 sex-independent DHS was used as the background for these enrichment calculations.

Enrichments of sex-biased DHS subsets (enhancer, insulator, or promoter) for mapping to sex-specific genes were calculated (e.g. male-biased enhancer DHS mapping to male-specific genes) as follows (Table S2B): Enrichment score = (ratio A) / (ratio B), where: ratio A = the number of male-biased enhancer DHS that map to male-specific genes, divided by the number of male-biased enhancer DHS that map to sex-independent genes; and ratio B = the number sex-independent enhancer DHS that map to male-specific genes, divided by the number of sex-independent enhancer DHS that map to sex-independent genes. For example, among the male-biased DHS classified as enhancers, 404 male-biased enhancer DHS map to male-specific genes, and 525 male-biased enhancer DHS map to sex-independent genes (404/525 = 0.77), whereas 1,495 sex-independent enhancer DHS map to male-specific genes, and 15,457 sex-independent enhancer DHS map to sex-independent genes (1,495/15,457 = 0.097), which gives an enrichment score A/B = 0.77/0.097 = 8.0. The set of 66,116 sex-independent DHS was used as the background for these enrichment calculations.

Enrichments of hypophysectomy-responsive DHS for mapping to class I or class II sex-specific genes (e.g. DHS that close in response to hypophysectomy in male liver and that map to male class I genes) were calculated for male and female liver as follows: Enrichment score = ratio A/ ratio B, where: ratio A = the number of DHS that open (or that close) following hypophysectomy and that map to class I (or to class II) sex-specific genes, divided by the number of DHS that open (or that close) following hypophysectomy and that map to sex-independent genes; and ratio B = the number of hypophysectomy-unresponsive DHS (static DHS) that map to class I or class II sex-specific genes, divided by the number of static DHS that map to sex-independent genes. For example, in male liver, 217 DHS that close following hypophysectomy map to a male class I sex-specific gene, and 1,203 other DHS that close map to sex-independent genes (217/1,203=0.1803), whereas 444 static DHS map to a male class I sex-specific gene, and 14,053 static DHS map to a sex-independent gene (444/14,053 = 0.0316), which gives an enrichment score A/B = 0.1803/0.0316 = 5.7. The static DHS used for these enrichment calculations correspond to the set of 35,562 static DHS in male liver and 30,394 static DHS in female liver.

Enrichments for overlap with genomic regions containing biologically relevant sets of transcription factor binding sites, chromatin marks and combinations of epigenetic features (chromatin states) were calculated for each of the following 4 sets of DHS: (1) dynamic male-biased DHS (834 sites), (2) static male-biased DHS (1,895 sites), (3) dynamic sex-independent DHS (1,532 sites), and (4) static female-biased DHS (1,359 sites) relative to a background set of static sex-independent DHS (64,584 sites) (see Fig. 2F). Briefly, these sets of DHS were identified based on their sex-specificity [19] and were classified as “dynamic” if they overlapped a STAT5-high vs. STAT5-low differential DHS, otherwise they were classified as “static”. Characterization of the 70,211 reference set DHS members as static or dynamic and their sex-specificity designations are listed in Table S2A. Biologically relevant regions annotated include sex-biased transcription factor binding sites for STAT5 and BCL6 [15], CUX2 [30], and FOXA1/2 [32]. We also examined DHS enrichment for sex-biased chromatin marks and chromatin states defined for male and female liver [18]. FOXA1 and FOXA2 ChIP-seq data [32] was processed and reanalyzed with MAnorm [56] in analyses performed by Gracia Bonilla of this laboratory, which identified sex-biased FOXA1 and sex-biased FOXA2 binding sites in mouse liver. Data from the MAnorm analysis of FOXA1 and FOXA2 ChIP-seq peaks (Table S8) was filtered by log2(fold-change) value (i.e., M-value) to define male-biased and female-biased binding sites. Genomic coordinates for each of the sex-biased transcription factor binding sites, chromatin marks, and chromatin states are provided (Table S5, Table S7).

Enrichments of DHS for overlapping biologically relevant regions were calculated (e.g., dynamic male-biased DHS for containing STAT5-High (Male) binding sites compared to the background set of static sex-independent DHS) as follows: Enrichment score = ratio A/ ratio B, where: ratio A = the percentage of dynamic male-biased DHS with a STAT5 binding site; and ratio B = the percentage of static sex-independent DHS with the same type of binding site. For example, 235 dynamic male-biased DHS contain a STAT5-High (Male) binding site (235/834 = 0.2818), whereas 311 static sex-independent DHS contain a binding site (311/64,584 = 0.0048), which gives an enrichment score A/B = 0.2818/0.0048 = 58.7. The set of 64,584 static sex-independent DHS (Fig. 2F) was used as the background for these enrichment calculations.

Chromatin state map analysis

- Chromatin state maps (14 state model), developed for male mouse liver, and separately for female mouse liver, are based on epigenetic data from a panel of six histone marks and DHS data in each sex [18]. These maps were used to identify sex differences in chromatin state and chromatin structure at each DHS region and their relationships to sex-biased gene expression, as follows. BEDTools was used to assign one of the 14 chromatin states to each of the 70,211 reference DHS based on the overlap of the DHS with the male liver chromatin state map, and separately, based on its overlap with the female liver chromatin state map (Table S2A). DHS whose genomic coordinates span two or more different chromatin states were assigned to the state with the largest number of overlapping base pairs, or in case of equal numbers of base pairs, the chromatin state with the smaller genomic coordinates. Enrichments of sets of static and dynamic male-biased DHS for being in one of the 14 chromatin states in male liver (e.g. dynamic male-biased DHS for being in chromatin state E6, which is an enhancer-like state) were calculated as follows: Enrichment score = (ratio A) / (ratio B), where: ratio A = the number of male-biased DHS that are in a particular chromatin state, divided by the number of male-biased DHS not in that chromatin state; and ratio B = the number of static sex-independent DHS in that chromatin state, divided by the number of static sex-independent DHS not in that chromatin state. For example, 604 dynamic male-biased DHS are classified as chromatin state E6, and 230 other dynamic male-biased DHS are not classified as chromatin state E6 (604/230 = 2.626), whereas 24,082 static sex-independent DHS are classified as chromatin state E6, and 40,502 static sex-independent DHS are not classified as chromatin state E6 (24,082/40,502 = 0.5946), which gives an enrichment score A/B = 2.626/0.5946 = 4.4. The set of 64,584 static sex-independent DHS was used as the background set for these enrichment calculations.

DHS peak normalization

- DNase-seq data to be visualized in the UCSC Genome browser (https://genome.ucsc.edu/) was normalized using the number of sequence reads in each DHS peak region per million mapped sequence reads (reads-in-peaks-per-million, RiPPM) as a scaling-factor [46]. Normalization was carried out using a comprehensive list of DHS peak regions merged across each dataset (i.e., the peak union list), obtained by concatenating FASTQ files for biological replicates of each control and treatment group, as described above. DHS peak regions identified in both individual liver samples and by analysis of the combined samples were concatenated into a single list, and then the BEDTools merge command was used to combine overlapping features to generate a single list of non-overlapping DHS peaks. The fraction of reads in peaks for each sample was then calculated to obtain a scaling factor. Raw sequence read counts were divided by this per-million scaling factor to obtain RiPPM normalized read counts.

Acknowledgements

– Supported in part by NIH grant DK121998 (to DJW).

Authors’ contributions

– AR and DJW conceived of the study. All animal studies and wet lab analyses including STAT5 EMSA gels, isolation of liver nuclei, DNase-seq analysis and library preparation were performed by JC. Data processing, primary computational analysis, preparation of analysis scripts, and data visualization were carried out by AR. All other data analysis and preparation of datasets and figures for publication were carried out by AR and DJW. AR and DJW jointly prepared a preliminary draft of the manuscript. The final manuscript was written by DJW and was reviewed and approved by all the authors. DJW supervised the overall project and edited the final manuscript for publication.

Supplemental figures

Fig. S1.

Fig. S2.

Fig. S3.

Fig. S4.

Fig. S5.

Fig. S6.

Fig. S7.

Fig. S8.