Circadian clock mechanism driving mammalian photoperiodism

The annual photoperiod cycle provides the critical environmental cue synchronizing rhythms of life in seasonal habitats. In 1936, Bünning proposed a circadian-based coincidence timer for photoperiodic synchronization in plants. Formal studies support the universality of this so-called coincidence timer, but we lack understanding of the mechanisms involved. Here we show in mammals that long photoperiods induce the circadian transcription factor BMAL2, in the pars tuberalis of the pituitary, and triggers summer biology through the eyes absent/thyrotrophin (EYA3/TSH) pathway. Conversely, long-duration melatonin signals on short photoperiods induce circadian repressors including DEC1, suppressing BMAL2 and the EYA3/TSH pathway, triggering winter biology. These actions are associated with progressive genome-wide changes in chromatin state, elaborating the effect of the circadian coincidence timer. Hence, circadian clock-pituitary epigenetic pathway interactions form the basis of the mammalian coincidence timer mechanism. Our results constitute a blueprint for circadian-based seasonal timekeeping in vertebrates.

identifies what is actually being signaled at the "photo inducible phase", I don't think there is a clean way to differentiate the models. Of course, I think it is possible that the expression and phase of expression of EYA3 is in fact the photo inducible phase. If this is the authors' conjecture, they should say so.
I would like the authors to clarify this. I would also like the authors to clarify the difference and similarity of the epigenetic and genomic changes associated with short day breeders such as sheep and deer with long day breeders such as hamsters and quail. In all 4 cases, DIO2 is induced in long days via similar mechanisms (although quail do not employ melatonin as its proximal photic transducer; their hypothalamus is directly photoreceptive). Yet, in sheep and deer, breeding is suppressed, while in hamsters and quail, breeding is activated. This is a short paper. I don't expect a long treatise, but I do think the paper would be greatly strengthened by mention of these two different interpretations.

Vincent Cassone
Reviewer #3 (Remarks to the Author): Attached as 'comments for author' Reviewer #4 (Remarks to the Author): I am afraid I am not an expert in circadian biology so I can only comment on the epigenetic part of the paper.
The analysis of H3K4me3 is potentially interesting and could illuminate the story but at present is underdeveloped. I have the following suggestions: 1. H3K4me3 mostly marks CpG island promoters of expressible genes, not necessarily expressed genes. No causality between H3K4me3 and transcription (or the reverse) has been properly established so this should be reflected in the language that's being used. H3K4me3 is not a measure of openness of chromatin; you would have to measure this directly by ATAC-seq or similar methods. So one cannot make any comparisons between H3K4me3 measurements and those of nuclear size or chromatin density. In fact ATAC-seq would have been an excellent thing to do and I'm also surprised the authors haven't done any H3K27me3 measurements, especially as they evoke similarities with plant vernalisation. 2. Please establish the relationship between H3K4me3 and CpG islands in your datasets. For this (and any other) comparison segment your gene classes into non-expressed, constitutively expressed, and seasonally changing (in one direction or other). 3. Please examine all seasonally changing genes (with the controls as above) as to their relationship with H3K4me3 density. Basically as I understand the hypothesis, seasonal genes show an association with ups and downs in K4 methylation, so this needs to be explored thoroughly not just with example genes. It would also be useful to check if K4 methylases and demethylases are seasonally expressed or not. 4. The choice of timepoints for the experiments needs to be explained better so that a nonexpert reader can understand. Fig 1c please explain the choice of timepoints. 5. There doesn't seem to be any SP D28 in Figure 1 g? 6. Please refrain from making parallels with vernalisation in plants as I believe there is no evidence to suggest there are any.

Response to reviewers: NCOMMS-20-01753 -Circadian clock mechanism driving mammalian photoperiodism AUTHOR RESPONSES TO IN RED
We thank the 4 reviewers for their evident care to detail, and comments. Each reviewer offered helpful and constructive suggestions, and as a result we feel that the paper has been significantly improved by their input. For this we are grateful.

Reviewer #1 (Remarks to the Author):
The paper by Wood et al. is a significant and novel piece of work, which provides a plausible mechanistic basis for photoperiodic time measurement in seasonal mammals. Attempting to the investigate the molecular basis of seasonal time measurement in (non-standard laboratory) mammals, such as sheep is no easy task but the authors provide some attractive and interesting data is support of a potential mechanism. The mechanism in particularly interesting as it builds on the Bunning hypothesis of external coincidence timing, originally developed for photoperiodism in plants. Furthermore, it suggests how a similar molecular mechanism may explain photoperiodic time measurement in both plants and mammals, pointing to an evolutionarily conserved process and therefore will have broad interest.
We thank the referee for the generally positive response to the manuscript and detail a response to specific points below: 1. The paper, and the mechanism put forward, hinges on the importance of two proteins -Bmal2 acting as a transcriptional co-activator, to drive long photoperiod expression of Eya3 and hence downstream TSH and DEC1 as a transcriptional repressor, induced by short photoperiod (long duration night time melatonin) that can block Eya3 expression and hence TSH. The transactivation studies showing the respective activities of these two proteins are quite convincing. Nevertheless, the identification of the two proteins (Bmal2 and DEC1) came from RNA seq studies of differentially expressed genes between long and short photoperiod (Bmal2, Fig 2A;  DEC1, Fig 3). While a sound bioinformatic approach was used to narrow in on these two proteins as key players in the photoperiodic timing mechanism, a crucial piece of evidence that is lacking from the study is the demonstration that the Bmal2 and DEC1 are actually found in transcriptional complexes bound to E-boxes in PT tissue. This is important to establish (see below).
We are aware that this is potentially an important issue, and whether the effects we observe are direct or indirect ( (Li et al, J Biol Chem. 2003;278(19):16899-907 2. In Fig 3E, the data to show that Bmal2 changes with both amplitude and phase relative to photoperiod are shown. However, it appears that Bmal1 and Bmal2 are annotated wrongly on the figure, given the description of the data in the text. Should the line marked as Bmal1 be Bmal2 and vice versa?
We thank the referee for this and apologise for this obvious error and have corrected the figure.

3.
A key element of the proposed mechanism is the role of Bmal2 in driving Eya3 expression and the authors note that Bmal2 peaks at ZT4 and ZT20 on LP. Yet in previous work (Dupre et 2010 (Current Biol. 20, 829-835), the same lab. observed that Eya3 peaks at ZT3 and ZT15 on LP. Thus, there seems to be some desynchrony between the timing of the initial peak in Eyas3 and the first peak in Bmal2 (ZT3 v ZT4), which is also evident in this study (in Fig 2b) and a major desynchrony between the second peak in Eya3 (at ZT15) and Bmal2 (ZT20). These observations do not seem to fully match up with the proposed mechanism as Eya3 appears to peak before Bmal2. Also, what is the function of the second peaks in Bmal2 and Eya3 and how do these relate to the proposed mechanism of photoperiodic time measurement. This issue is an issue the authors need to address and hence why it is important to establish the nature of the Additionally, we have collaborated with colleagues at the University of California in some aspects of structural biology studies of BMAL2. These are ongoing, but show that BMAL2 has undergone a substantially different evolutionary path to the "principal circadian regulator BMAL1". There is marked lack of conservation of multiple domains, with the marked exception of the N-terminal "Gregion. Strikingly, the rodents appear as an "out-group" with a marked loss of selection. Impressively, there is greater sequence variation in BMAL2 between rodents and other vertebrate lineages than there is between Coelacanth and Primates. A current (un-tested) hypothesis is that the conserved G-region may significantly contribute to photoperiodism. This protein has clearly taken a markedly different evolutionary path to BMAL1, or indeed other clock genes. Figure 1F, it is shown that the size of the nucleus changes in PT folliculo-stellate (FS) cells as well as in PT thyrotrophs, yet from this study we know that DEC1 is not expressed in PT FS cells. This suggests that photoperiod driven epigenetic events are not necessarily unique to PT thyrotrophs and begs the questions (a). what role do the PT FS cells play in the photoperiodic timing mechanism, if epigenetic events are important as the authors suggest? and (b) does DEC1 explain all?

In the supplementary
FS cells may play an important role in pituitary endocrine regulation, and indeed we reported such an effect in an earlier paper on seasonal prolactin regulation (Dupre et al 2010). Here, we have focused on the EYA3/TSH circuitry, which is exclusively a property of the PT thyrotroph. We do not discount that some of our seasonal changes in RNA and histone modification could be present in the FS cells. Previously we have shown that there are extensive zona adherens, desmo-somes, and gap junctions between FS and thyrotope cells in the PT (Wood et al 2015) and previous studies have shown interactions between FS and other pituitary cell types (Baes et al 1987 Endocrinology 120, 685-691;Allaerts et al 1990. Mol. Cell. Endocrinol. 71, 73-81). Thus changes in protein expression in the TSH cells may be transmitted to FS cells via these junctions, but the factors involved will require further investigation, outwith the scope of the present paper. Therefore with specific reference to the reviewer's comments: ( Minor point: 1. Figure numbers in the text and on the figures need to be made consistent. Figure  1A, B,C etc are used in the text, but Figure 1a,b,c etc are used on the figures and in the legends.

Reviewer #2 (Remarks to the Author):
Wood et al. present extensive and interesting data from Scottish black face castrated male sheep, a short day breeder, showing genomic and epigenetic correlates of seasonal, photoperiodic timing. They show that photoperiodic changes in the duration of melatonin, which is well established to faithfully reflect the length of the scotoperiod, induce BMAL2 expression when the duration of melatonin is short, which in turn triggers the transcription of EYA3. EYA3 has already been shown to instigate TSHB activity, which increases DIO2 deiodinase activity. Conversely, they show that long durations of melatonin, indicative of the short days of winter induces the repressor DEC1 (as well as DEC2 and REVERBa and b), which suppresses BMAL2 activity and shutting off the EYA3/TSHb process.
The authors indicate that their data support the notion that photoperiodism in the sheep is synchronized by a "circadian coincidence" hypothesis in which the photoperiod entrains a "circadian rhythm of photosensitivity , and the expression of summer or winter biology depends on whether or not light coincides with the phase of high sensitivity". This is specifically the "external coincidence model" of Erwin Bunning (1936), who used this model to explain plant photoperiodism. The authors rightfully cite Bunning's work.
However, in course of their writing, the authors lose the "external" part of the coincidence and merely talk about "coincidence". This may seem to be trivial, but there is another coincidence model. The "internal coincidence model" proposed by Nanda and Hamner (1958) suggested that the same phenomenon could be explained by two circadian oscillators that are differentially entrained to dawn and dusk. Then, as photoperiod changes, the phase relationship between the two oscillators change, inducing or suppressing a seasonal event.
We thank the reviewer for this comment, and we are aware of the distinction in the 2 models. We also admit to shamelessly ducking the issue in the submitted version of the paper, because it is complicated, leading us to use the undefined politically neutral term "co-incidence". However, we have now updated the introduction to highlight the differences: Page 3 line 53-57: "An "internal coincidence" model has also been proposed where the role of light is to entrain two circadian oscillators and the phase relationship between the two oscillators determines the response to photoperiod 2,3. In either case the role of the circadian clock is central and formal studies support the universality of a "coincidence timer" in animals 2-7, but we lack understanding of the mechanisms involved".
The authors present an interesting set of data in Figure 3 that shows that the phase and amplitude of expression of several "clock genes" or clock-controlled genes as well as occupancy of selected promoter-motifs differentially entrain to the 2 photoperiods. In this reviewer's mind, this is consistent with an internal coincidence model, not Bunning's hypothesis. Frankly, until someone identifies what is actually being signaled at the "photo inducible phase", I don't think there is a clean way to differentiate the models. Of course, I think it is possible that the expression and phase of expression of EYA3 is in fact the photo inducible phase. If this is the authors' conjecture, they should say so. I would like the authors to clarify this. figure 4). We have also addressed the complexities internal vs external coincidence timing in the discussion but have not been expansive because this is a nuanced and conceptual view that may not be useful for a wider audience:

Differentiating between internal and external coincidence is a challenge, which is why we have stuck to the neutral "coincidence model". With the data set from figure 3 we did attempt to assess the phase relationships and protein-protein interactions between the circadian regulated genes in LP and SP to see if we could fit an internal coincidence model. However, these analyses were somewhat inconclusive, and therefore were not included in the manuscript. Our analysis in figure 3 does show a clear shift in peak phase expression (figure 3c & d) but also a very clear peak at ZT4 LP (12 hours after dark) which could be attributed to the presence of a photoinducible phase (external) or a change in phasing (internal)but there are a number of unique genes expressed at this time point, coherent with the concept of a release from repression on short photoperiods in the presence of light (external). And indeed we identify SP repressors that target these LP ZT4 induced genes (supplementary
Page 12 line 333-334: "While we cannot differentiate between an internal or external coincidence model of photoperiodic time measurement in this study, generally the…" I would also like the authors to clarify the difference and similarity of the epigenetic and genomic changes associated with short day breeders such as sheep and deer with long day breeders such as hamsters and quail. In all 4 cases, DIO2 is induced in long days via similar mechanisms (although quail do not employ melatonin as its proximal photic transducer; their hypothalamus is directly photoreceptive). Yet, in sheep and deer, breeding is suppressed, while in hamsters and quail, breeding is activated.
This is a short paper. I don't expect a long treatise, but I do think the paper would be greatly strengthened by mention of these two different interpretations. While these are interesting points and we do not think a discussion of long day and short day breeders will add to our PT-focused manuscript, and believe other researchers cover this more comprehensively.

Reviewer #3 (Remarks to the Author):
Importance: While some key molecular factors in the regulation of photoperidicity have been uncovered in the last ~10 years, the mechanism is not well understood. Elucidating this mechanism is important since photoperiodicity drives much of reproductive behaviour in economically important mammalian species. Moreover, uncovering the basic mechanism may well yield insights into how environmental cues interact with molecular signalling pathways that modify a range of physiological processes including reproduction and feeding behaviour. The coincidence model predicts that photoperiod entrains a circadian rhythm of photosensitivity, and that the expression of summer or winter phenotype depends on whether the external light period matches this photosensitive period. Previous work by this group and others have identified the pars tuberalis (PT) of the pituitary gland as a key target organ of pineal melatonin production. The onset of darkness should trigger a circadian oscillator that reaches maximum amplitude during daylight period in summer and night period in winter. Simultaneously, the pineal gland produces melatonin during the hours of darkness. The integration of these steps should output a signal that reaches maximum amplitude in the photosensitive phase in summer, but not in winter. Previous work has demonstrated that the transcription factor EYA3 obeys this behaviour, and is necessary to activate downstream targets of the seasonal response such as bTSH. The outstanding question is: How do circadian oscillator pathways interact with melatonin signalling to modulate EYA3 expression. In order to answer this question the investigators performed a series of temporal transcriptional profiling experiments in order to identify genes encoding transcriptional regulators of EYA3 in the target organ, the ovine PT. Overall, this manuscript represents a considerable step forward in the description of the process of photoperiod regulation. The authors have generated high quality datasets from a well phenotyped model system, and these data will be a considerable resource for the field for years to come. In addition, the authors have identified some new molecular players. However, the paper suffers from a lack of clarity about many of the experimental protocols used, and would be significantly improved by extending the methods sections and showing more supporting evidence of their experimental protocols. Finally, some conclusions are not supported by the experimental data. My detailed critique is below.
We thank the reviewer for their positive comments and hope we have addressed the issues of clarity in the experimental protocols and that our additional analysis and textual clarifications ensure our conclusions are supported by the data. See detailed response below.
Part 1 (mostly Fig2 and associated datasets). The investigators clearly justify the timepoints chosen and demonstrate that seasonal endocrine pathways and EYA3 transcription are appropriately dynamically differentially regulated. Next, the transcriptional datasets at these timepoints are interrogated to discover differentially expressed genes that have high expression in the LP but are not induced in the SP. Several genes are identified, but the investigators further refine their investigation to a single target, BMAL2. This is a known Ebox gene that had not been previously implicated in seasonality. In in-vitro assays the team shows that overexpression of BMAL2 can activate the EYA3 promoter in concert with BMAL1, CLOCK and TEF, via the PAS-B domain. Finally, using archival histological material, the authors demonstrate PTspecific induction of BMAL2 gene expression by light in the LP but not the SP.
 The in-vivo methodology and transcriptomics experiment is not well describedcrucially, how many replicates were performed per time-point? This should be clearly stated throughout.
We apologise for the lack of clarity and have made the following changes to clarify: In the results section: Page 5 line 102: " (Fig. 1b, Supplementary Fig. 1a- Fig. 3a (diurnal comparison at day 28), and, 2). The experiment presented in Fig. 4a (melatonin implant study). Animals were blood sampled throughout the study and terminally sampled at the indicated time-points (Figure 1b, 3a, 4a) (Fig. 1b)". Page 16 line 422-424:" Ovine melatonin was measured by radioimmunoassay as previously described 44 for animals in the melatonin implant study (Fig. 4a, n=6  We have also updated the ChIP-seq, ISO-seq and CAGE-seq, and transcription factor binding site analysis. The methods state that the time-course data was fitted to a 5th order polynomial (photoperiod x time), and significant genes generated from this. What is the justification for this method (should be referenced)?
We have added further extensive clarification of our model selection approach within the methods, see above response and the manuscript. To give specifics in response to the reviewer; the polynomial regression approach we used is similar to that within the maSigPro pipeline for time-series analysis (Nueda, et al 2014). We used Akaike information criterion (AIC) to investigate the optimal model selection for expressed genes, balancing model overfitting and underfitting (using the Oshlack and Gordon selectModel implementation in limma). It is not possible to select a single model that is optimal across all genes, however for genes > 0 log2CPM and with an amplitude > 1.5 we found that including orthogonal polynomials up to 5 th order was optimal for the most genes in both SP and LP time-series. Panel F shows the number of optimal genes for each model investigated, and is divided into bins for expression ranges to demonstrate how including higher order polynomials becomes more important for the more highly expressed genes.
 Other genes (apart from BMAL2) display the predicted pattern of gene expression and could presumably amplify/modulate the EYA3 response. Is it justified for the authors to conclude that BMAL2 is the only important player?
Our seasonal analysis focuses on ZT4 and identified 48 genes that were upregulated at the first long day (fig 2a, sup table 3 Figure 2c,  d, e, supplementary figure 2g and h, figure 4h. Note that supplementary figure 5b, c and d have not altered because these are representative plots, we apologise for the oversight and have updated all the figure legends accordingly. Part 1 overallthe investigators demonstrate that BMAL2 is an important new player in the induction of EYA3 activity, with potential role as part of the key PT molecular driver of photoperiodism. However, experiments should be explained much more clearly, controls (e.g. protein overexpression) shown and stats recalculated with appropriate n.
We have tried to address all the reviewers points above and hope they suffice.
Part 2 (mostly Figure 2 and associated data). Next, the investigators perform experiments to identify the factor that prevents high amplitude expression of BMAL2 in the SP. They reason that such a factor should have high expression during the period that is dark in the winter but light in the summer (ZT8-16). By performing high density temporal transcriptomic profiling the authors generate a high quality, detailed molecular description of the PT in the diurnal cycle in both SP and LP. In ZT20 of the SP functional enrichment analysis identified negative regulators of transcription, notably including E-box regulators. Several enriched genes at ZT4 in the LP contained an E-box in their promoters.
 In Figure 2f the authors show that BMAL2 expression is high in the light period of LP but not SP pituitary. Moreover, they state that BMAL1 expression is not altered by photoperiod. The data in Figure 3e appears to contradict these points. Is the graph miss-labelled? This was a mistake and the plot was mislabelled.
 The String analysis linking ZT20 repressors and ZT4 targets is poorly described. What is the statistical basis for this analysis? Would any 2 randomly chosen timepoints show a similar level of connectivity, for example? This method should be clearly described and justified. Is there a more rigorous way of showing this interaction? We have clarified the approach in the methods section and in the results section. Furthermore we have added an additional plot to supplementary figure 4, now 4b. In brief we used the known protein-protein interactions to test whether genes, which are increased in expression at LP ZT4, were more likely to be targets by the products of genes expressed at SP ZT20 (our repressor genes). See the updated text for the results section below: Page 10 line 258 -267: "This led us to ask; are SP ZT20 repressors more likely to interact, and therefore potentially repress, LP ZT4 up-regulated genes. We used curated experimental protein-protein interaction (PPI) observations from the STRING database, which contain known protein-protein interactions and functional associations 47. We found that LP ZT4 up-regulated genes are more enriched within the SP ZT20 repressor network (P-value = 0.001) than down-regulated genes ( Supplementary Fig 4a, b). We did not find significant enrichment when considering genes that are differentially up or down regulated across the whole day (P-value = 0.25) (Supplementary Fig. 4b)".
The methods were clarified, Page 24 line 688, 695-697: ". The significance of enrichment of PPI repressor connected genes within upregulated vs down-regulated genes was evaluated using fishers two-way exact tests".
The figure legend was updated for figure 4b: "b. Fig. 4a). The significance of enrichment (fishers two-way exact test) for PPI connectivity within the up-regulated vs downregulated genes is shown."

Number of differentially expressed genes (DEGs), up-regulated (yellow) and down-regulated (blue) for a daily mean (all 24 hour timepoints) between SP vs LP on day 28, compared to LP ZT4 vs SP ZT20 contrast. Connectivity via protein-protein intereactions (PPI), defined by STRING to transcriptional repressors is indicated by the checkered shading (also represented in
 Part 2 overall. The experiments described generate a valuable dataset, but the analysis is poorly described and superficial, without adequate justification. We have made extensive additions to the explanations of the methodologies used throughout the manuscript, see response to reviewer in point 1 and below: Page 22-23 line 625 -651: "For gene annotation, five tissue samples were sequenced over two experimental runs using PacBio Iso-Seq. In the first run PT and PD samples were sequenced from an RNA pool of SP and LP Scottish blackface sheep (N=1) and a pineal from a commercial mule sheep from Manchester, UK. This RNA was sent to GATC Biotech (Konstanz, Germany) for cDNA library preparation using their in-house method with mRNA 5' cap and poly(A) tail selections and sequencing on a PacBio RSII system. GATC made full length normalized RNA libraries.size selected for <2kb, 2kb-4kb, >4kb. sequenced across 75 PacBio RS II SMRT cells (SRX7688275). In a second run, PT from a pool of sheep in LP, and SP (N=3) were sequenced. RNA was extracted using RNeasy Mini Kit (Qiagen) with on-column Dnase digestion. A fulllength cDNA library was constructed for each sample using the TeloPrime Full-Length cDNA Amplification Kit V1 (Lexogen) and amplified using PrimeSTAR GXL DNA Polymerase (Takara Bio) with 22 PCR cycles of 98 °C denaturation for 10 seconds, 60 °C annealing for 15 seconds, and 68 °C extension for 10 minutes. PacBio SMRTbell libraries were prepared using SMRTbell Template Prep Kit 1.0 and each library was sequenced on two SMRT Cells v2 LR using 20-hour movies on a Sequel platform at the IMB Sequencing Facility (University of Queensland, SRX7688271). All Iso-Seq data was first processed using software IsoSeq v3.1 to obtain full-length non-concatemer reads with at least 3 full sequencing passes, which were then mapped to the sheep reference genome GCA_002742125.1 using GMAP version 2018-05-30. TAMA Collapse from the TAMA tool kit 81 was used to generate unique gene and transcript models, which were further merged with RNAseq-based annotation data using TAMA Merge to incorporate any transcript models that were identified by RNAseq but not Iso-Seq. Functional annotation of transcripts was carried out using Trinotate (v3.1.1)".
Page 23-24 line 657 661 & 664-666: "We applied cap analysis gene expression (CAGE) to identify the location and relative expression of TSS regions of the PT across both LP and SP. When combined with IsoSeq and RNASeq derived transcript annotation this provided a comprehensive identification of TSS in the genome which allowed us to more accurately apply DNA binding motif analysis to promoter regions" ……." We sequenced archived RNA samples from the PT in both SP and LP (ZT4, week 12) 29. We also sequenced RNA from PD (both SP and LP), and Pineal for comparison as outgroups. Reads were trimmed using fastx toolkit 0.0.14 and cutadapt 1.4. Reads were mapped using BWA 0.7.17 to the 5th release of the sheep genome (Oar_rambouillet_v1.0; assembly GCA_002742125.1). CAGEr 1.26.0 was used for processing and cluster analysis of TSS (Supplementary Table 4). We filtered reads for a mapping quality > 30 and sequencing quality > 20. Tag counts were normalised using the power law method with an alpha of 1.12 and T of 106 (deterimined by plotting the reverse cumulatives of PT samples). We clustered TSS with >1 TPM together using the distclu methods allowing a max distance between TSS of 20 nucleotides.".
Part 3 (Figures 3 and 4). In this part the investigators explore the regulation of SP repressors by melatonin, using animals that are kept constantly in light to inhibit melatonin release, and by replacement. The model is validated by demonstrating the expected expression levels of the previously described melatonin-responsive CRY1 gene. Putative repressors identified in Part 2, REVERB-a , CHRONO and DEC1 show a dynamic expression pattern consistent with melatonin regulation. DEC1 was chosen for further analysis since it is a known E-box suppressor of the circadian clock. DEC1 is expressed in the right place for its expected function (thyrotroph cells of the PT). In vitro, overexpression of DEC1 was able to suppress BMAL2-induced EYA3 activation.  As in part 1, the luciferase experiments are poorly described, expression of overexpressed protein is not shown and the stats is inappropriate. Please see above response regarding the luciferase assays.
Part 3 overall. This part shows convincingly that REVERB-a, CHRONO and DEC1 are regulated by melatonin and may contribute to the negative regulation of BMAL2 in the SP Thank you for your positive comment.
Part 4. Epigenetic analyses. In supplement to the transcriptional profiling performed in Part 1, the authors carry out ChIP-seq experiments to examine the binding of the histone modification H3K4me3. This is a mark that is associated with the promoters of actively transcribed genesindeed, active RNA-PolII recruits this mark to the promoter in concert with transcription. This mark has also been associated with intragenic regions and enhancer sequences. The authors report coincidence of H3K4me3 binding with gene expression levels, and state that this correlation is 'stronger on LP indicating a global activation of gene transcription on LP associated with increased chromatin accessibility'.
 The analysis of the ChIP seq dataset is poorly described in the methods and impossible to follow. How is H3K4me3 density called? The authors should clearly present their analysis workflow and extend the methods section to include sufficient information to allow this analysis to be replicated. The methodologies have been updated to include a better explanation of the analysis completed, we have also included a workflow diagram, see supplementary figure 7 (and Panel H) Fig. 6a). Furthermore, H3K4me3 peaks were well-associated with CpG islands (CGIs) on the sheep genome as described in the previous study ( Supplementary Fig. 6b)  In addition we included and additional workflow figuresupplementary figure 7, the ChIP-seq workflow is shown in Panel H.
 The correlation of gene expression with 'H3K4me3 density' is a very unusual way of presenting this kind of data. Generally ChIP seq peaks are called using a sliding window approach to differentiate enriched vs background levels. The peak is subsequently treated as digitalie, present or absent (rather than a continuous variable). Peak calls are then mapped to reference genome Panel H -ChIP-seq workflow datasets to identify, e.g, association with transcriptional start sites, introns, repetitive DNA, etc. The authors should show this analysis for their data to confirm i) sequence distribution of their H3K4me3 peaks to different genomic features, ii) overlap of peaks with expressed genes at TSS/elsewhere with a clearly defined window.
We thank the reviewer for their suggestions of additional analyses and have completed these and present the peak call distributions in supplementary figure 6 and Panel I. The distributions of H3K4me3 peaks were closely resembled previous reports and H3K4me3 peaks were well-associated with CpG islands (CGIs) on the sheep genome. For details please see the manuscript at page 22 line 605-623 (also pasted in above response point 10). We also updated the results Page 5 line 11-121.
We also not that the phrase "H3K4me3 density" was a terminology issue on our part and we have updated the text to either use "peaks" or "marks" throughout.
For the correlation analysis we improved our approach and updated figure 1d (see Panel J). Note: We changed the labels of figures to "H3K4me3 peak read counts around TSSs" from "H3K4me3 density". The figure below shows  TSSs. Seasonal genes were identified from differentially expressed genes between SPday84 and LPday1,day7,day28 and between LPday112 and SPday1,day7,day28 (log Table 3), and observed a strong correlation between seasonal gene expression and H3K4me3 peaks around the transcription start sites (TSS's) (Fig. 1d). Importantly, this correlation was absent in non-seasonally regulated genes (Fig. 1d)  What is the evidence that H3K4me3 density drives chromatin accessibility? In any case, the correlation analysis described in Figure 1c does not show this. It shows rather that genes that are down-regulated in the critical interval are associated with reduced relative H3K4me3 deposition, and genes that are upregulated have increased H3K4me3 deposition. This is entirely expected since PolII activity recruits H3K4me3. The gold standard to show that chromatin accessibility is increased globally in LP (as the authors imply by their inclusion of nuclear density measures (Fig1d and suppl)) is to perform chromatin-conformation capture type experiments, which measure the interactivity of distal parts of the genome. We have updated the manuscript throughout to clarify that our findings are related to h3k4me3 marks not chromatin accessibility per se. However, the relationship between H3k4me3 marks and chromatin accessibility has been previously reported with h3k4me3 regulating gene expression through chromatin remodelling by the NURF complex making chromatin more accessible for transcription factors (Wysocka Nature. 2006;442(7098):86-90). Obviously it would be interesting to conduct an ATAC-seq experiment to look specifically at accessibility, as our EM data clearly indicates this is a feature of LP and a whole range of histone and chromatin modifiers could be responsible for this. The manuscript has been updated and improved, especially at pages 5 -6 line 103-150.  In the discussion the authors state that the progressive increase in H3K4me3 binding at EYA3 over the weeks leading into the LP is at odds with previous data where they show that at the individual level there is a binary switch in gene expression. I disagreethe quantitative increase in H3K4me3 that the authors observe at the TSS of EYA3 could represent a binary H3K4me3 TSS deposition in an increasing number of cells across the time interval. Again, the histone mark is likely to be secondary to the transcriptional activation. We entirely agree with the reviewer and apologise because we must not have written these statements clearly. We have updated the text as follows, page 13 line 364-370: "Our earlier work shows that at the individual cell level the transition between winter and summer physiology is a binary, all-or-nothing phenomenon 29,56. Integrating these two findings, we suggest that individual thyrotroph cells of the PT exhibit a distribution of critical day length requirements/sensitivity for circadian triggering of the summer physiology leading to a binary switch in cell phenotype, which in a whole tissue assay would appear as a progressive change in epigenetic status".
 The authors note the interesting finding that EYA3 appears to have a seasonal TSS. Is this a more generalised phenomenon? The authors should interrogate their dataset to assign H3K4me3 peaks to 'canonical' and 'noncanonical' TSSs. Should there be multiple instances of this, motif analysis of these sites could potentially uncover seasonality-specific regulatory mechanisms. We thank the reviewer for pointing this out and have completed the suggested analysis, which is now included in supplementary figure 1g,h (Panel K) and updated the manuscript accordingly: Page 6 line 134-144: "We noted that approximately 70% of seasonally DEGs (Supplementary Table 3) had more than one TSS compared to only ~20% of the PT genomic background (Supplementary Fig. 1g, Supplementary Table 4). Next we took genes that were up-regulated in either SP or LP and plotted the proportion of genes with multiple TSS's, and repeated this for H3K4me3 marked TSS's, this revealed that H3K4me3 marks are more likely to occur on genes with multiple TSS's ( Supplementary Fig. 1h) and highly expressed seasonal DEGs have a greater prevalence of multiple TSSs than non-seasonal genes expressed at the same level. Furthermore, this phenomenon was more pronounced in LP than SP. This indicates an enrichment of multiple TSS's in LP up-regulated genes which is associated with the H3K4me3 mark".
Page 42 line 1217 -1220: "g. Percentage of genes with a given number of transcription start sites in the genomic background (all >0 log2CPM expressed genes) of the pars tuberalis (grey bars) as compared to all seasonally differentially expressed genes (white bars)".  Part4 overall. This part of the work is the weakest. It does not really increase the informative value of the paper in the current state. Moreover, I feel that the conclusions drawn are not supported by the data. These data could be removed with minimal impact on the overall quality of the work. Should the analysis of seasonal-specific TSS be extended and prove a generalised phenomenon, then this is a significant ant interesting new finding. We believe that this is the first demonstration of genome wide seasonal changes in epigenetic status and therefore is of value to report. We have taken the reviewers useful comments on board and strengthened the data regarding the correlation of seasonal gene expression and h3k4me3 counts at TSSs and demonstrated that the multiple TSS phenomenon is a more generalized feature. We believe these are significant and interesting findings that should remain in the manuscript.

Reviewer #4 (Remarks to the Author):
We thank the reviewer for their helpful comments, much of this reviewers suggested changes and analyses overlapped with reviewer 3, so to avoid repeating information we endeavour to reference to those responses where relevant and note figure numbers and page/line references in the manuscript rather than full copying of text.
Panel K -g. Percentage of genes with a given number of transcription start sites in the genomic background (all >0 log2CPM expressed genes) of the pars tuberalis (grey bars) as compared to all seasonally differentially expressed genes (white bars). H. Comparison of the prevalence of multiple (>1) TSS across different gene cohorts. The cohorts are LP 28 days up-regulated DEGs (solid red), SP 28 days up-regulated DEGs (solid blue) and all PT expressed genes as the background (solid black). Prevalence of multiple TSS on genes is shown for all thresholds (%) of the uppermost expressed genes (i.e. increasing thresholds for the upper quantile of gene expression). The equivalent gene expression (log2CPM) values for upper quantiles for a lower threshold cutoff upper x-axis. Dashed lines indicate the proportion of gene in the cohorts with multiple H3K4me3 (>1) marked TSS. hypothesis, seasonal genes show an association with ups and downs in K4 methylation, so this needs to be explored thoroughly not just with example genes. It would also be useful to check if K4 methylases and demethylases are seasonally expressed or not.
Reviewer 3 (point 11) also encouraged us to improve our correlation analysis and we have update figure 1d, see panel J in this document and our response to reviewer 3. In the manuscript Page 6 line 122 -133 and Page 36-37 line 1037-1045.
The reviewer's second point relates to the regulators of the h3k4me3 marker. Histone modifications are precisely balanced by methyltransferases ("writers"), demethylatases ("erasers") and effector proteins ("readers"), so we took the known writers, readers and erasers (Hyun 2017) and checked our rnaseq data. We found no seasonal changes in expression. These data are added to Supplementary Table 5 and the manuscript is updated: Page 6 line 127-133: "Histone modifications are precisely balanced by methyltransferases ("writers"), demethylases ("erasers") and effector proteins ("readers"), therefore we checked the RNA expression of H3K4me3 readers, writers and erasers but found no seasonal changes (Supplementary Table 5). This suggests that changes in protein activity of H3K4me3 modulators may be key in the observed seasonal alterations in H3K4me3 marks.)".
4. The choice of timepoints for the experiments needs to be explained better so that a nonexpert reader can understand. Fig 1c please explain the choice of timepoints. We have now updated Fig1c (now fig 1d) to include all comparisons to remove a time-point bias, please see response to point 4 above. Furthermore the manuscript has been updated extensively to improve the description of experiments and timepoints, see response to reviewer 3 points 1,4,8,9, and  We have also made a number of changes in the results to clarify and explain the analysis better, especially pages 5 -7 lines 103 -173.
5. There doesn't seem to be any SP D28 in Figure 1 g? Figure legend 1g has been updated to clarify the exclusion of SP day 28 and 84. "Note: SP day 28 and 84 are not included because the H3k4me3 peaks are nondetectable".
6. Please refrain from making parallels with vernalisation in plants as I believe there is no evidence to suggest there are any. We acknowledge that our references to vernalisation "muddy the waters" and have undertaken to remove direct parallels and only use it in the introduction as a means to set up the idea of seasonal epigenetic regulation.