This Article Abstract Full Text (PDF) All Versions of this Article: biophysj.104.053975v1 88/4/2369    most recent Alert me when this article is cited Alert me if a correction is posted Services Similar articles in this journal Similar articles in PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Google Scholar Google Scholar Articles by François, P. Search for Related Content PubMed PubMed Citation Articles by François, P.

A Model for the Neurospora Circadian Clock

Laboratoire de Physique Statistique, Centre National de la Recherche Scientifique, UMR 8550, Ecole Normale Supérieure, Paris, France

Correspondence: Address reprint requests to Paul François, Ecole Normale Supérieure, Physics, 24 rue Lhomond, Paris 75231, France. Tel.: 33-1-44-32-3763; E-mail: francois{at}tournesol.lps.ens.fr.

	ABSTRACT

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION ACKNOWLEDGEMENTS REFERENCES

 
Circadian clocks are important biological oscillators that generally involve two feedback loops. Here, we propose a new model for the Neurospora crassa circadian clock. First, we model its main negative feedback loop, including only experimentally well-documented reactions, the transcriptional activation of frequency (frq) by the white-collar complex (WCC), and the post-transcriptional dimerization of FRQ with WCC. This main loop is sufficient for oscillations and a similar one lies at the core of almost all known circadian clocks. Second, the model is refined to include the less characterized enhancement of white-collar 1 (WC-1) protein synthesis by FRQ, the positive second feedback loop. Numerical testing of different hypotheses led us to propose that the synthesis of WC-1 is enhanced by FRQ monomers and repressed by FRQ dimers. We demonstrate that this second loop contributes significantly to the robustness of the oscillator period against parameter variation. A phase response curve to light pulses is also computed and agrees well with experiments. On a general level, our results show that explicit time delays are not required for sustained oscillations but that it is crucial to take into account mRNA dynamics and protein-protein interactions.

	INTRODUCTION

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION ACKNOWLEDGEMENTS REFERENCES

 
Circadian clocks are important examples of genetic oscillators used to synchronize organisms to the daily cycle of light and dark. Circadian rhythms have been widely studied for many years (Daan and Pittendrigh, 1976), and recent works have unveiled the detailed mechanisms of this internal timing in several organisms (Young, 2002; Reppert and Weaver, 2001, 2002). Clocks from different organisms appear to share common features. Their core component relies on one feedback loop including at least two genetic interactions, a positive and a negative one. At least two proteins or groups of proteins are involved in these genetic interactions. The first group of proteins, the activating proteins, interacts with the DNA and activates the transcription of genes corresponding to the second group of proteins. In coordination with some post-transcriptional modifications, this second group of proteins usually interacts in the cytoplasm and in the nucleus with the activating proteins, forming multimers unable to activate transcription. These proteins are consequently repressing proteins. The aim of this article is to describe and to study these interactions in a simple system where the core components of this main feedback loop have been well described, the Neurospora crassa circadian clock (Loros and Dunlap, 2001), and to compare it to the experiments. For this fungus, circadian rhythmic growth patterns were described 50 years ago (Pittendrigh et al., 1959). With advances in molecular biology, understanding of its circadian clock has improved, and main components of this clock have been determined and in some cases described in a quantitative way (Garceau et al., 1997; Ballario et al., 1998; Merrow et al., 1997; Lee et al., 2000; Froehlich et al., 2003).

In the following, a model of the Neurospora circadian clock main loop is first proposed and compared to available experimental data. Biological interactions are modeled with mass-action laws so that the necessary delays in the clock are the consequence of the well-described chemical reactions. The model of the core loop appears to correctly describe oscillations of frq transcripts and FRQ proteins, but does not account for the observed WC-1 oscillations. To describe them, a second positive loop involving the enhancement of WC-1 synthesis by FRQ (Lee et al., 2000) needs to be taken into account. Refined models are proposed and tested. This leads to the specific proposal that WC-1 translation is enhanced by FRQ monomers and suppressed by FRQ homodimers. The positive feedback loop is found to enhance robustness of the clock to parameter variations. Light response of this model is also computed, and is found to be in good agreement with the experiments.

Some experimental results about the Neurospora circadian clock
The Neurospora circadian clock is based on an autoregulatory negative-feedback loop with three proteins: the FREQUENCY protein, FRQ, the repressing protein; and white-collar proteins WC-1 and WC-2, the activating proteins. Here we summarize the main experimental findings (for detailed reviews, see Loros and Dunlap, 2001; Dunlap et al., 2004).

The gene frq is historically one of the first to have been identified as a part of the core Neurospora's circadian clock (Feldman and Hoyle, 1973). In constant darkness, frq RNA and FRQ proteins concentrations oscillate. The peak of frq transcript is followed after 4–6 h by a somewhat larger peak of FRQ proteins (Fig. 1 A redrawn from Garceau et al., 1997). The circadian cycle can be divided in two precisely defined phases (Merrow et al., 1997). The first phase is the negative feedback itself (repression) in which FRQ represses its own transcription. This phase is 14–18 h long. The second phase (de-repression) is simply the recovery from this repression when frq transcript level returns to high concentrations. This step is 4 h long.

View larger version (18K): [in this window] [in a new window]   FIGURE 1  (A) Redrawn from experimental data from Garceau et al. (1997). The symbol x is the relative concentration of frq RNA; and the symbol + is the relative concentration of FRQ protein. (B) Redrawn from experimental data from Froehlich et al. (2003) showing rhythmic binding of WCC to the FRQ promoter. The symbol + is the densitometry of the C-box-bound complex; and the symbol x is the densitometry of FRQ protein. (C) Redrawn from experimental curves from Lee et al. (2000), showing antiphase oscillations of FRQ and WC-1. The symbol x is the WC-1 protein level; and the symbol + is the FRQ protein level.

The coiled-coil-domain-mediated FRQ-FRQ interaction is also necessary to circadian oscillations (Cheng et al., 2001a). It seems necessary for the interaction between FRQ and WCC, but its precise role in the circadian oscillation is still unknown.

Recent work showed that the level of WC-1 proteins also oscillates, in phase opposition to FRQ proteins, as shown in Fig. 1 C redrawn from Lee et al. (2000). However, the wc-1 transcript level does not vary throughout the day. Dunlap and co-workers established that this mechanism is triggered by the presence of FRQ proteins (Lee et al., 2000). Also, when the frq gene is knocked-out, WC-1 level is very low compared to the level in wild-type cells. These experimental results all suggest an enhancement by FRQ of the production of new WC-1 proteins from the existing transcripts. However, the detailed mechanism of this enhancement is still unknown.

It was finally shown that transcription of WC-2 is also activated by FRQ, but in a nonrhythmic way (Cheng et al., 2001b). In this work, Cheng and co-workers engineered quinic acid (QA)-controlled strains of Neurospora and also observed that despite considerable changes in the levels of WC and FRQ proteins due to induction by QA, the period of the clock changed only slightly. The amplitude of the clock can therefore vary whereas its period remains constant.

	METHODS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION ACKNOWLEDGEMENTS REFERENCES

 
Model characteristics
For simplicity and because of lack of precise experimental data, different cellular compartments and separate concentrations for the nucleus and cytoplasm are not considered. All concentrations are referenced to the cell volume, so that concentrations represent the effective number of molecules present in the cells. As WC-2 is expressed in large excess compared to WC-1, we assume that WC-2 quickly interacts with WC-1 to form WCC.

The one-loop model
The first two equations model the transcriptional regulation of frq transcripts by WCC:

(1)

(2)

Note that frq represents the frequency gene without WCC bound to its promoter; frq:WCC denotes the gene bound to WCC; and RNA stands for frq transcripts. As there is only one copy of each gene per cell during vegetative growth, [frq] + [frq:WCC] = 1 gene per cell, so that only one equation is necessary to describe regulation of the promotion. In details, WCC proteins can bind to the frq gene with a rate and the bound protein is released with a rate (Eq. 1). When WCC is bound to the frequency gene, transcription is initiated with a rate _FRQ. We assume a first-order degradation for RNA with a constant degradation rate _RNA.

The following differential equations stand for protein productions and regulations:

(3)

(4)

(5)

In these equations, [FRQ], [WCC], and [T] denote, respectively, the concentration of FRQ protein, of WCC, and of the multimer formed by FRQ and WCC.

FRQ is translated from the transcripts with a rate ß. The complex WCC is assumed to be produced with a fixed rate _WCC. FRQ and WCC proteins can form a multimer T with a rate . It is supposed that FRQ only binds to the free form of WCC and does not bind to the WCC protein bound to the frq promoter. Additionally, the complex formed by WCC and FRQ is not able to promote transcription. Finally, the dissociation of the complex is neglected. A schematic representation of the one-loop model is presented in Fig. 2 A.

View larger version (12K): [in this window] [in a new window]   FIGURE 2  (A) Schematic representation of the one-loop model; (B) schematic representation of the second two-loop model.

Models with another positive feedback loop
A positive feedback loop has been observed in Neurospora: thanks to a post-transcriptional mechanism, FRQ enhances WC-1 production. However, the precise mechanism mediating this enhancement is still unknown. To provide some indication on the type of possible biological mechanism, we study two additional models. Both of these models rely on the fact that FRQ might interact with the wc-1 transcripts and enhance their transcription. In the following, WC-2 is still supposed to be in large excess and to quickly form a dimer with WC-1.

First two-loop model
In this model, we suppose that FRQ proteins interact with wc-1 transcripts, forming a complex. The complex between FRQ and wc-1 transcripts is translated with a delay after its formation. So additional equations describe the dynamics of wc-1 transcripts, a species not explicitly taken into account in the previous one-loop model:

(6)

(7)

The value RNA_W stands for the normal form of the wc-1 transcripts, and RNA_W+ stands for its enhanced form. The wc-1 transcripts are supposed to be produced with a fixed rate _WCC, and their normal form to interact with FRQ to form a complex with a second-order reaction rate . The reverse reaction is possible with a rate µ. Note that after some time, the total concentration of the wc-1 transcripts is constant. The equation for FRQ must be modified as

(8)

Finally, the enhancement of translation is simply supposed to appear after a delay (which phenomenologically accounts for undescribed biological processes),

(9)

There are two different translation rates in this last equation: ß_– stands for the translation rate of the normal form and ß₊ stands for the translation rate of the complexed form. [RNA_W+] stands for [RNA_W+](t – ), where is the delay in translation of the second type of RNA.

Second two-loop model
As will be seen, the delay needs to be quite long to reproduce experimental data. This delay should result from well-defined biochemical interactions and several hypothetical interactions were tested to see which model could agree with the experimental data, as there is no precise description of the activation of WCC by FRQ yet.

This led us to propose a second two-loop model without any explicit delay. In this second two-loop model, FRQ proteins are still supposed to directly interact with wc-1 transcripts, and form a complex. The complex between FRQ and wc-1 transcript is translated after its formation without any delay. In this model we take into account the homodimerization of FRQ. It is hypothesized that another FRQ protein can interact with the FRQ protein bound to the wc-1 transcript, and that the complex formed by this FRQ dimer cannot be translated. The new equations for the wc-1 transcript are

(10)

(11)

(12)

The value RNA_W stands for the normal form of wc-1 transcripts, RNA_W+ stands for the enhanced form complexed with FRQ, and RNA_W– for the form complexed with two FRQ proteins. We suppose that wc-1 transcripts are produced with a fixed rate _WCC, and that the wc-1 transcripts' normal form can interact with FRQ to form a RNA-protein complex with a second-order reaction rate . The reverse reaction is possible with a rate µ. The enhanced form can interact once again with FRQ.

For protein-protein interactions, we suppose that FRQ homodimerizes with second-order rate and that homodimers dissociate into two FRQ proteins with rate . Then FRQ dimer (FRQ₂) can interact with WCC to form the multimer FRQ₂:WCC. Equations for FRQ, WCC, FRQ₂, and FRQ₂:WCC consequently are

(13)

(14)

(15)

(16)

A schematic representation of the second two-loop model is presented in Fig. 2 B.

Parameters choice
A logarithmic plot of the frq transcript concentration from Garceau et al. (1997; not shown) shows that its decay is exponential. The first-order degradation rate is of the order of 0.2 h^–1. This value was already proposed by Ruoff et al. (1999) with a fit of the behavior predicted by the Goodwin model. However, it is not possible in the present model to know if this exponential decay is due to the real degradation rate of RNA or simply to the detachment constant of WCC from the frq promoter. Mathematical analysis of the one-loop model (P. François, unpublished) shows that the equation for RNA decay is r(t) = A exp(–t) + B exp(–_RNAt), where A and B are two constants. If > _RNA, the RNA decay is mainly directed by its own degradation (parameter _RNA), since, after a short time, exp(–_RNAt) >> exp(–t). If < _RNA, it is directed by the dynamics of the detachment of WCC from the FRQ promoter (parameter ). We therefore simulated both behaviors and saw no major qualitative difference between the models.

Experiments from Lee et al. (2000) provide the WC-1 degradation rate, and show that in the presence of FRQ, the WC-1 degradation rate is not affected by FRQ. Therefore, the degradation rate of the multimer WCC-FRQ _T is almost the same as the WC-1 degradation rate _WCC. Examination of the Western blots provides an approximate value of this rate. The WC-1 concentration is divided by 3 within 4 h. This gives

For the FRQ degradation rate, several parameters were tested, with no qualitative differences. Actually, it seems difficult to find its precise value from experimental data, as in the present models, the FRQ concentration decrease is directed by at least three parameters: , _RNA, and _FRQ. We therefore tested two set of parameters for the one-loop model: one with _FRQ = 0.05 h^–1, the other with _FRQ = 0.25 h^–1. The behaviors of the transcripts and the proteins for these two sets of parameters are similar. For the first two-loop model, _FRQ = 0.05 h^–1 gave the best results. For the second two-loop model, we chose _FRQ = 0.3 h^–1, a value similar to the WCC degradation rate and close to the value proposed by Ruoff et al. (1999) from fits of the Goodwin model.

It is not possible to deduce the values of the other parameters from curves with linear rise and decay without knowing absolute concentrations. Actually, in the models, it is possible to rescale parameters to obtain almost any absolute concentrations. For instance, in the one-loop model, if we multiply both _WCC and ß by a same constant c and divide and by the same c, the qualitative dynamics of the system will not be changed, despite the change of absolute concentrations. We therefore chose parameters so that both kinetic constants and proteins concentration seem in a physiological range, and fit the oscillator period and the experimental curves. Similar oscillations occur for a very large set of parameters, so that the qualitative behavior of the oscillator is mostly independent of the choice of parameters.

Numerical methods
Integration of differential equations was performed with a Runge-Kutta algorithm. The time step was reduced until no significant difference in simulations appeared after further reduction. To ensure that the real asymptotic limit-cycle was observed, simulations were performed until no difference appeared between successive oscillations.

All programs were written in C++.

	RESULTS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION ACKNOWLEDGEMENTS REFERENCES

 
The one-loop model is considered and studied (Eqs. 1–5) in the sections "A simple model with only mass-action laws can simulate the Neurospora crassa circadian oscillator" and "RNA control and mechanism of repression". The section "Models with a positive feedback loop give different qualitative behaviors for WCC" is devoted to the analysis of the two-loop models. The section "Parameters variation, period dependence, and compensation" is devoted to an analysis of parameter dependence and to a comparison between models. The section "Phase response curve" studies the light response of the second two-loop model.

A simple model with only mass-action laws can simulate the Neurospora crassa circadian oscillator
The Neurospora one-loop model simulates oscillations of the levels of frq transcripts and FRQ proteins (Fig. 3 A and Fig. 4 A). The delay between the frq transcripts and the FRQ proteins peaks is 6 h, in agreement with experimental observation. The decay of frq transcripts is exponential and requires 18 h in agreement with the experiments. The de-repression phase, when frq transcripts rise to their peak level, is 4 h long. The behavior of the concentration of FRQ proteins is clearly not a simple exponential. This seems also to be the case for the experimental curves. Finally, WCC is observed to rhythmically bind to frq promoter (Fig. 4 B).

View larger version (17K): [in this window] [in a new window]   FIGURE 3  Simulation of the one-loop model of Neurospora circadian clock, with RNA < . (A) Absolute concentration of proteins and transcripts. (B, solid line) Free WCC, + symbol: plot of (ß[RNA] – FRQ[FRQ])/([FRQ]), confirming that we can make a quasistatic assumption to relate [FRQ] to [WCC]. (C, solid line) Free FRQ (absolute concentration), + symbol: plot of (WCC – WCC[WCC])/([WCC]), confirming that we can make a quasistatic assumption to relate [WCC] to [FRQ]. Constants for this first model are WCC = 3.75 mol h–1, FRQ = 7.5 mol h–1, ß = 0.7 h–1, = 0.35 h–1, = 10 mol–1 h–1, WCC = 0.3 h–1, RNA = 0.2 h–1, FRQ = 0.05 h–1, and = 2000 mol–1 h–1, where mol stands for molecules and h for hours. Degradation constant of the complex is the same as WC-1. Parameters were chosen to have oscillations qualitatively similar to the experimental curves presented in Fig. 1 A.

View larger version (20K): [in this window] [in a new window]   FIGURE 4  Simulation of the one-loop model of Neurospora circadian clock, with parameters of Fig. 3. Concentrations on A and B are rescaled by their maximum values for comparisons with experimental curves. (A, dotted line) Total FRQ (free + complexed). (A, solid line) FRQ RNA. (B, solid line) Total FRQ. (B, dotted line) Binding activity on frq promoter.

De-repression phase
During this phase, free (not complexed with WCC) FRQ concentration is low and free WCC concentration is high. The newly formed FRQ quickly interacts with the free WCC, which is in large excess. A consequence of this interaction is that the concentration of free FRQ is proportional to the concentration of frq transcripts, and inversely proportional to free WCC concentration, as shown by the comparison between WCC concentration and this quasistatic assumption provided in Fig. 3 B.

At the same time, the free WCC interacts with the frq gene promoter. The concentration of frq transcripts consequently rises exponentially. When frq transcripts reach high concentration, almost all free WCC disappears and FRQ concentration can rise again. At the end of this phase, frq transcripts are near their maximum level.

Repression phase
FRQ free concentration is now high whereas WCC free concentration is low. This time, the newly formed WCC immediately interacts with the free FRQ in excess, produced with a high rate because of the high concentration of frq transcripts. A consequence of this interaction is that concentration of free WCC is inversely proportional to free FRQ concentration, as shown by the comparison between FRQ concentration and the quasistatic assumption provided in Fig. 3 C.

At the same time, bound WCC is released in an exponential way from the frq promoter. The frq transcription rate decays in the same way, and transcripts are degraded in an exponential way. As the FRQ production rate is proportional to transcript concentration, when the transcript level is low, the production rate becomes too low, and finally all the free FRQ is degraded. WCC can once again accumulate, and a new cycle begins.

For comparison with the experimental curves, the total (free+complexed) concentration of proteins should be taken into account. In the one-loop model (Eqs. 1–5), this concentration is almost constant for total WCC, because the degradation constant of free WCC is the same as the degradation constant of the complex (Fig. 3 A). The FRQ curves in the model are similar to the experimental ones (compare Fig. 1 A with Fig. 4 A).

RNA control and mechanism of repression
If the dynamics of WCC production and of the dimerization is fast enough compared to the characteristic time constants of frq transcripts, there is, effectively, a dynamical switch between WCC and FRQ. The low concentration species is in quasiequilibrium, and its dynamics is slaved to the high concentration species.

In simple terms, when WCC concentration is higher than FRQ concentration, WCC proteins titrate all free FRQ and after a very short time, only the free WCC, with dimers of WCC and FRQ, remains. Inversely, when FRQ concentration is higher than WCC concentration, FRQ proteins titrate all the WCC and after a very short time, the free FRQ, with dimers of WCC and FRQ, remains. Consequently, both proteins cannot be simultaneously present in uncomplexed form in the cell with comparable concentrations, and part of the protein in excess and all the low concentration proteins are stored in the complex. A first consequence of the dimerization is, therefore, that the dynamics of both free proteins is controlled by the concentration of frq transcripts: When FRQ is in excess, free FRQ concentration is controlled by the concentration of the frq transcripts, and controls WCC free concentration thanks to the dimerization; and when WCC concentration is high, the FRQ production rate is proportional to the concentration of the frq transcripts. The produced FRQ proteins quickly dimerize with free WCC so that the free WCC sequestration rate actually is proportional to the transcript concentration.

Finally dimerization explains the repression mechanism: when FRQ is present at high concentration, it titrates WCC and therefore prevents its binding to the frq promoter.

To sum up, the core mechanism of the clock can be reduced to two coupled mechanisms—a slow dynamical process composed by all the transcription machinery, mathematically described by Eqs. 1 and 2, coupled to a rapid switch at the protein level, mathematically described by Eqs. 3 and 4. This switch, in turn, controls the slow process through the transcriptional activation.

Models with a positive feedback loop give different qualitative behaviors for WCC
The one-loop model does not take into account the regulation of WC-1 by FRQ, so that the WC-1 level is almost constant. Introduction of the experimentally observed second feedback loop is needed to explain WC-1 oscillations. It has been established that the activation of WC-1 by FRQ is mediated through a post-transcriptional mechanism (Lee et al., 2000). Here, we supposed that the translation of wc-1 transcripts is enhanced by FRQ.

In experiments, the WC-1 concentration peak is out of phase with the FRQ concentration peak, so that this activation seems delayed. One simple way of modeling such a phase shift is to introduce a phenomenological delay in the equations, which effectively describes some biological processes such as cellular transport, for instance.

This strategy is followed in our first two-loop model in Eqs. 6– 9 (see also Fig. 5). It gives realistic amplitudes for both FRQ and WCC oscillations. However, the delay length is critical, to qualitatively match the WCC behavior seen in experiments and to quantitatively account for the delays between peaks and for the period length. In the present model, to obtain a realistic behavior, the delay needed to be set to 7 h. As this delay is quite long, we therefore tried to find possible mechanisms explaining it.

View larger version (23K): [in this window] [in a new window]   FIGURE 5  (A) Simulation of the first two-loop model with a post-transcriptional activation of WC-1 by FRQ. Oscillations of RNA (solid line), total WCC (dashed line), and total FRQ (dotted line). (B) Rescaled concentrations. Constants for this model are WCC = 0.3 mol h–1, FRQ = 10 mol h–1, ß = 0.6 h–1, = 0.6 h–1, = 1 mol–1 h–1, WCC = 0.3 h–1, RNA = 0.2 h–1, FRQ = 0.05 h–1, = 100 mol–1 h–1, T = 0.3 h–1, = 0.2 mol–1 h–1, µ = 0.01 h–1, ß– = 1 h–1, ß+ = 40 h–1, and = 7 h, where mol stands for molecules and h for hours. Degradation constant of the complex is the same as WC-1. Parameters were chosen to have oscillations qualitatively similar to the experimental curves presented in the top and third panels of Fig. 1.

The second two-loop model, shown in Eqs. 10–16, is an attempt to model such a mechanism. In this model, there are no explicit delays. Instead, we suppose that FRQ monomers bind to wc-1 transcripts to enhance their translation and that, on the contrary, binding of FRQ dimers to wc-1 transcripts represses their translation. Consequently, at high concentration, FRQ homodimerization prevents WCC translation. At low concentrations of FRQ, FRQ proteins essentially exist as monomers and can bind to wc-1 transcripts to strongly enhance their translation, so that, as shown in Fig. 6, WCC peak occurs just after FRQ minimum, which explains the observed out-of-phase relationship between FRQ and WCC. If this homodimerization is switched off, the oscillations disappear and WC-1 is overexpressed (data not shown).

View larger version (20K): [in this window] [in a new window]   FIGURE 6  Simulation of the second two-loop model with a post-transcriptional activation of WC-1 by FRQ. (A) Absolute concentrations of total FRQ (including all complexes and dimers), total WCC, and frq transcripts. (B) Rescaled oscillations of RNA, WCC, and FRQ. (C) Binding of WCC of frequency promoter (solid line): binding activity (dashed line): FRQ concentration. (D) Influence of WCC transcription rate over FRQ oscillations (solid line): WCC = 2.5 mol h–1, and (dashed line): WCC = 5 mol h–1. Constants for this model are WCC = 2.5 mol h–1, FRQ = 75 mol h–1, ß = 1 h–1, = 0.25 h–1, = 0.003 mol–1 h–1, WCC = 0.3 h–1, RNA = 0.2 h–1, FRQ = 0.3 h–1, = 1600 mol–1 h–1, T = 0.3 h–1, = 8000 mol–1 h–1, µ = 0.1 h–1, ß– = 0.001 h–1, ß+ = 10 h–1, = 3000 mol h–1, and = 10 h–1, where mol stands for molecules and h for hours. Degradation constant of the complex is the same as WC-1. Parameters were chosen to have oscillations qualitatively similar to the experimental curves presented in the top panel of Fig. 1.

To gain a better understanding of the model parameter dependence it is first useful to rescale variables. In the following, we set first and and g = [frq].

Second, we define the following rescaled parameters: let be a = ß_FRQ/(_WCCRNA), b = /_RNA, and ; and we rescale time by taking a new time unit t₁ = _RNAt. Taking time t₁ as new time t, the new ODEs for the one-loop model are

(17)

(18)

(19)

(20)

As an example, parameters used for the model in Fig. 3 give = 2.3 x 10^–3, a = 7, b = 1.75, c = 2.2, _F = 5.8 x 10^–4, d = 23, and _W = 3.5 x 10^–3.

For the one-loop model, all parameters can be varied over one order of magnitude without destroying oscillations (data not shown). Noteworthy is that, if we vary the parameters of the one-loop model, there are several subcritical Hopf bifurcations so the oscillations often appear with finite amplitude. This kind of hysteretical transition has previously been proposed as a possible mechanism for noise resistance in genetic oscillators (Barkai and Leibler, 1999). A consequence of these subcritical Hopf bifurcations is that, depending on the initial conditions, a stable limit-cycle can coexist with a stable fixed point. This phenomenon is illustrated in Fig. 7 for the full one-loop model.

View larger version (18K): [in this window] [in a new window]   FIGURE 7  Coexistence of a stable limit cycle and a stable fixed point. (Solid line) Free FRQ for a first set of initial conditions. (Dotted line) Free FRQ for a second set of initial conditions. First set of initial conditions is [WCC] = 0.0017, [RNA] = 10, [FRQ] = 5, and [frq] = 0.54. Second set of initial conditions is [WCC] = 0.0022, [RNA] = 4.44, [FRQ] = 5, and [frq] = 0.9. Units are molecules per cell. Constants for this model are WCC = 9 mol h–1, FRQ = 9 mol h–1, ß = 2.25 h–1, = 0.4 h–1, = 20 mol–1 h–1, WCC = 0.3 h–1, RNA = 0.2 h–1, FRQ = 0.2 h–1, and = 800 mol–1 h–1, where mol stands for molecules and h for hours. Degradation constant of the complex is the same as WC-1.

View larger version (15K): [in this window] [in a new window]   FIGURE 8  Scaling of the amplitude of the oscillator as a function of for the adimensioned system. Amplitude of (W) is plotted as a function of time t (not rescaled). The amplitude of W does not depend on in the limit of small , this shows that the amplitude of W scales as 1/. Parameters are a = 3, b = 5, c = 10, d = 1, and W = F = 0.

1. Period of the oscillation. At the lowest order in the small parameter , five parameters (_RNA, a, b, d_w, and d_f, with d_f = _F/ and d_w = _W/) are crucial for period determination. The period is given by and d_w). The value a corresponds to the ratio of the WCC protein production rate over the FRQ protein production rate, taking into account both transcription and translation (_WCC is an effective rate that can actually be seen as: WCC transcription rate x WCC translation rate ÷ Degradation rates of the wc-1 transcripts). The value b is the ratio between the release constant of WCC by the frq promoter and the degradation rates of the frq transcripts. The values d_f and d_w are the respective ratios of FRQ and WCC degradation rates over the degradation rates of the frq transcripts.
   
2. Amplitude of the oscillation. In the limit of small , the amplitude of the oscillations of free proteins scales as 1/ (Fig. 8).
   

View larger version (9K): [in this window] [in a new window]   FIGURE 9  Effect of strong parameter variations on oscillations. Each point represents a simulation where one parameter has been doubled or halved, keeping the other parameters constant. Amplitude and period are plotted relative to the amplitude and the period of the control. Measured amplitude is the amplitude of total FRQ oscillation. (solid triangle, RNA; solid circle, ; solid square, FRQ; open triangle, WCC; open square, FRQ; open circle, ß; open diamond, WCC; and x symbol, other parameters.) (A) One-loop model (parameters of Fig. 3). (B) One-loop model (with FRQ = 0.25 h–1, FRQ = 20 mol h–1, = 0.23 h–1, and = 4 mol–1 h–1, with other parameters the same as Fig. 3). (C) First two-loop model (parameters of Fig. 5). Five parameters destroy oscillations. (D) Second two-loop model (parameters of Fig. 6).

we find b = /exp(–(E_r–E)/RT), so that b remains constant if E_r = E = E. If we impose that the period does not vary >5% when temperature varies between T₀ = 301 K and T₁ = 311 K, we find E < 4000 J/mol, which seems a reasonable assumption and in the range of values proposed by Ruoff and Rensing (1996). Taking into account the degradation of the proteins (if they are not negligible) imposes a similar condition: the activation energy of the degradation rate must be of the same order as the activation energy of RNA degradation rate.

Another possibility is to consider that temperature has very little influence on dissociation of WCC from frq gene and over the degradation rates. This assumption is discussed below.

Finally, if production rates are modified while keeping the parameter a constant, the phase of the clock shifts, so that the oscillator can be entrained to a 24-h period (data not shown). So, importantly, the clock can be temperature-compensated and still entrained by changes in temperature.

Comparison between models: robustness to parameter variation
To gain insight into parameter dependence and to compare models, we doubled and halved constant rates for each reaction, one at a time, keeping the others constant. Results of this computation are given in Fig. 9 for two typical sets of parameters of a one-loop model and for both of the two-loop models.

For the first set of parameters for the one-loop model, the most sensitive parameters are transcription and translation rates (_WCC, _FRQ, ß), as well as the degradation rate of frq transcripts and of FRQ proteins and the detachment rate of WCC from the frq promoter, with period variations from 13% to 63%. For other parameters, the period never varied by >3% from control. This dependence is explained by the mathematical analysis summarized before. For the first one-loop model, amplitude variations are correlated to period variations: in most cases, period and amplitude of the oscillations vary with the same order of magnitude.

We also tested a set of parameters with a higher FRQ degradation constant. The main dynamical consequence of this second choice of parameters is that the binding of WCC to FRQ promoter is much weaker (i.e., no quasistatic in the activation phase) than in the first set of parameters. With this second set of parameters, WCC proteins do not saturate the FRQ promoter. As can be seen on Fig. 9 B, this enhances the global robustness to parameter variations, and the most sensitive parameters are _RNA, , and _FRQ.

The first two-loop model is not very robust to parameter variations. The oscillations disappear after modifications of _WCC, wc-1 transcription rate, and of four new parameters: the delay for the activation of translation, the wc-1 transcript degradation rate the enhanced translation rate ß₊ of these transcripts, and the rate of interaction between FRQ proteins and wc-1 transcripts. All these parameters are implicated in the delaying processes in the positive feedback-loop. The other sensitive parameters are the same as in the one-loop model.

The second two-loop model is far more robust to parameter variations. The three most sensitive parameters for the period are , the rate for the release of WCC from the frq promoter; and the degradation rates of frq transcripts and FRQ dimers, with variations of the period from 20% to 27%. Doubling the WCC proteins' degradation rate also gives a period that is 15% shorter. For all other parameter variations, the period does not vary by >9%, as can be seen on Fig. 9 D. Contrary to the one-loop model, it is possible to have large amplitude variations without modifying the period of the clock. The amplitude of the second two-loop model still depends on synthesis rates, but its period is much less sensitive (compare parts A and D of Fig. 9). For instance, a doubling of the WCC transcription rate modifies the amplitude of the oscillations without modifying the period, as can be seen in Fig. 6 D.

Possible role of the positive feedback loop
Experiments show that in Neurospora, it is possible to have large variations of the amplitude of the oscillations while keeping the period constant (Liu et al., 1998; Cheng et al., 2001b). According to the previous analysis in the one-loop model, if transcription or translation rates for both proteins are multiplied by the same factor f, a remains constant, and the period does not change; but from the expression of , the amplitude of the adimensioned variable is multiplied by and the real amplitude of the protein by f. Cheng et al. (2001b) proposed that one possible role of the positive feedback loop was to precisely adjust protein production rates to keep the period constant. The second two-loop model supports this suggestion: when synthesis rates of proteins are modified, the oscillator adjusts itself to keep the period constant, as can be seen in Fig. 6 D and Fig. 9 D. When the parameter _WCC is multiplied by 4, the diminution of the period is <2% of the reference period, whereas the amplitude is more than three times higher; when it is multiplied by eight, the diminution of the period is <7% and the amplitude is approximately five times higher (data not shown). For comparison, for the second set of parameters for the one-loop model, when the same parameter is multiplied by 4, the period is 25% lower than the reference period, and when it is multiplied by 8, the period is 50% lower than the reference period. The second two-loop model consequently shows that specific biochemical mechanisms in the positive feedback loop could help in keeping the period constant, despite changes in the amplitude.

Phase response curve
The precise biochemical processes mediating the response of Neurospora circadian clock to light pulses are still not completely known. It was shown first that the effect of light pulses is to switch off the negative feedback loop (Crosthwaite et al., 1995). A light pulse first greatly increases the production of the frq transcripts (fourfold to 25-fold as compared to the average level during one cycle). Then, these newly formed transcripts are quickly degraded compared to normal transcripts (with half-lives of the order of 1 h). It was more recently shown that there are two specific binding sites for light response in the frq promoter called light-response elements (Froehlich et al., 2002). And finally, another negative feedback loop, implicating at least one gene, called vivid, has been discovered. VIVID seems to negatively regulate (but not to fully control) the gating of light input, probably via hyperphosphorylation of WC-1 (Heintzen et al., 2001).

A precise model of this feedback loop is not possible yet because of the lack of more precise experimental data. It is possible, however, to test some light response properties of the Neurospora circadian clock. Dunlap and co-workers hypothesized that the phase of the clock was given by the concentration of the frq transcripts (Crosthwaite et al., 1995; Loros and Dunlap, 2001), and that the effect of light was to switch the concentration of the frq transcripts to its maximum. We tested this heuristic model by suddenly raising the concentration of the frq transcripts to its maximum at different times of the cycle. We also introduced a supplementary effect due to the other negative feedback loop: we supposed that one role of vivid was to trigger the degradation of WC-1 (as proposed by Heintzen et al., 2001) and set WCC total concentration to zero, also including a degradation of the WCC bound to the frq promoter. Light pulses at different times of the cycle delay (negative-phase shift) or advance (positive-phase shift) the oscillations by different amounts. Phase response curves (PRCs) show the phase shift that corresponds to light pulses at different times of the cycle.

The PRC was computed for the one-loop model and for the second two-loop model. These PRCs are very similar, and only the PRC for the two-loop model is shown in Fig. 10.

View larger version (12K): [in this window] [in a new window]   FIGURE 10  Phase response curve. (+ symbols) Response due to a sudden rise of RNA concentration to its maximum level with a simultaneous degradation of all WCCs at different times in the cycle. (Solid line) The variation of frq transcripts during a circadian cycle (without light pulses) is also shown for reference.

	DISCUSSION

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION ACKNOWLEDGEMENTS REFERENCES

 
The model explains oscillations, and can be improved to have robustness
The one-loop model (Eqs. 1–5) describes basic features of Neurospora circadian oscillations. Both qualitative and quantitative behaviors of frq transcripts and FRQ protein are reproduced by this model. Reactions at the level of the transcription were supposed to be slower than protein-protein reactions, which completely explain phase shifts and qualitative behaviors. The two phases described in Merrow et al. (1997) clearly appear: the repression phase corresponds to the phase where FRQ protein concentration is high, and the de-repression phase corresponds to the phase where free WCC protein concentration is high. Rhythmic binding is a natural consequence of this switch at the protein level.

The explicit distinction in the equations between transcription and translation of the RNA is necessary to explain experimental curves. Modeling protein production by a single effective step leads to the destruction of the oscillations, in the present one-loop model. Actually, the FRQ protein and the frq transcripts are seen experimentally to have different behaviors, since FRQ concentration variation does not simply reproduce messenger RNA variation after a time delay (Fig. 1 A). This experimental result clearly stresses the need for separate modeling of transcription and translation. The level of RNA has an unexpected consequence on the protein concentration: for instance, in the de-repression phase, the rate at which free WCC proteins are sequestered is controlled by the concentration of the frq transcripts.

Protein-protein interactions are at the core of the system, and should be taken into account explicitly. Such reactions are not equivalent to simple repression at the transcriptional level: heterodimerization is essential to couple RNA and proteins dynamics, and therefore plays a dynamical role different from a simple repression.

The present one-loop model also helps us to understand the clock temperature compensation as described in Temperature Dependence (above); when the production rates of proteins are varied the same way, in the limit of small , the period does not change despite a change of the amplitude of oscillation.

The fact that FRQ oscillations are observed for the one-loop model shows that the second positive feedback loop is not necessary for the occurrence of oscillations, confirming a previous study (Smolen et al., 2002). To take into account WC-1 oscillations, another feedback loop is required. We tested the hypothesis of a direct activation at the post-transcriptional level, without any delays, mediated, for instance, by some intermediate enzymes (models not shown). However, it was impossible to obtain realistic out-of-phase oscillations for FRQ and WC-1 with such models. To obtain out-of-phase relationships between FRQ and WCC, it was necessary to suppose that FRQ represses WCC production at high concentration, and activates it at low concentrations. This was modeled by taking into account the homodimerization of FRQ (Cheng et al., 2001a). This hypothesis gives behaviors in agreement with experimental observations. Besides, if the transcription rate of WCC is raised, the oscillation amplitudes are higher but the period changes only slightly (Fig. 6 D). In this model, one of the roles of the positive feedback loop could therefore be to adjust protein production rates to keep the period constant, despite amplitude changes, as proposed by Cheng et al. (2001b).

A light-PRC was also computed. It was shown that taking into account both the production of the frq transcripts and the hyperphosphorylation of WC-1 proteins explains the shifts of the clocks. This model of light influence confirms Crosthwaite et al. (1995), and the role of frq transcripts in the phase determination of the clock.

Comparison with other works
Even if major components of circadian clocks have been well-described experimentally, the dynamical origin of the oscillations remains quite unclear. Actually, most models of circadian clocks can be classified in two categories: models where delays are necessary to oscillations; and models where oscillations only depend on the specific assumptions made about the genetic interactions.

Examples of models with delays have been proposed by Smolen et al. (2001, 2002) and Goldbeter and co-workers (Gonze et al., 2000; Leloup and Goldbeter, 1998, 2003).

Smolen et al. (2001) hypothesize that the delays observed in circadian clocks are consequences of slow biological processes (due to transcription, translation, or cellular transport, for instance) and use delayed differential equations to model circadian oscillations with phenomenological delays accounting for these mechanisms. The main conclusion of their models is that a positive feedback loop is not necessary to have oscillations, but that long time-delays (7 h in the Neurospora case) are necessary to account for oscillations.

In the present one-loop model, we explicitly modeled transcription and translation. This one-loop model shows that explicit delays in the equations are not necessary to produce oscillations. Then we took into account the second feedback loop to better explain the biological data. Two models were formulated: a model with explicit delays and a model without delays but with supplementary biophysical interactions. These two models present the same qualitative behavior. However, their properties are different. We showed, for instance, that our second two-loop model is far more robust to parameter variations than our first two-loop model. This means that it may not be possible to reduce the second two-loop model to a simplified version such as the first two-loop model with delays, without destroying some important properties of the model.

Goldbeter and co-workers (Gonze et al., 2000; Leloup and Goldbeter, 1998, 2003) have also intended to explicitly model the delaying mechanisms. In the Neurospora case, WCC activity has not been considered, but in the mammalian case, the corresponding proteins dynamics (BMAL, CLOCK) has been modeled. Goldbeter and co-workers have hypothesized that nuclear transports and successive phosphorylations observed in most of circadian clocks are at the origin of delays and are necessary for the oscillation. Only the hyperphosphorylated form of the proteins has been supposed to form heterodimers to repress transcription. Actually, experimental studies showed that FRQ is quickly phosphorylated (Garceau et al., 1997) and that its phosphorylation rate determines its degradation rate (Liu et al., 2000). Also, hypophosphorylated FRQ is also known to be able to bind to WCC (Yang et al., 2002). Therefore, in the present model, phosphorylation has been supposed to fix the degradation rate of FRQ, which is one of the most important parameters for the determination of the period length.

Some other models do not introduce slow processes, and these models suppose that specific interactions in the genetic network help in destabilizing the fixed point.

The Ruoff-Rensing model (Ruoff and Rensing, 1996) is essentially based on the Goodwin model (Goodwin, 1965). Transcription and translation are explicitly modeled. As in the present model, a slow amplification process (transcription and translation) is coupled with a rapid switch, at the transcriptional level. Dynamics of repression is modeled by a Hill function accounting for fast kinetics. For the system to oscillate, a high Hill exponent (>9) is needed, which implies a very high cooperativity. This acts as a phenomenological switch, accounting for possible mechanisms of undescribed origin. The dynamics of the present model is close in spirit to the mechanism suggested by the Ruoff-Rensing model, with a slow accumulation of RNA and proteins coupled to a rapid repression. For instance, the influence of degradation rates predicted by the present model is very similar to what was proposed before for the Goodwin oscillator (Ruoff et al., 1999) and was confirmed experimentally by Liu et al. (2000). However, the present model bypasses the need of high cooperativity by taking into account the interaction between FRQ and WCC at the post-transcriptional level.

Another interesting model was proposed for the Drosophila circadian clock (Tyson et al., 1999) with a goal similar to that of the present model: to provide a minimal model, simple to analyze and to improve. As in the present model, dimerization played an essential role, but in Tyson et al. (1999), the crucial positive feedback loop was a consequence of stabilization of PER induced by this dimerization. However, the two models are quite different: Tyson and co-workers concluded that a positive feedback loop was required to explain oscillations. Such a positive feedback loop is not needed when one does not make any quasiequilibrium assumptions on the dynamics of the proteins as shown by the one-loop model of the present article. We propose that the role of the positive feedback loop is, rather, to improve robustness to variations in parameters.

Finally, in the previously described models of circadian clocks, regulation of transcription was modeled by Hill or Michaelis-Menten kinetics, modeling fast binding between DNA and protein, and in most of the models, quasiequilibrium is also assumed for RNA dynamics. Kinetics at the level of transcription is therefore supposed to be very quick. This focus on proteins requires us to make specific assumptions on the dynamics of the networks to have oscillations. The present models show that if one models transcription and translation and does not make any quasiequilibrium assumptions, both oscillations and biological delays can be explained without any further hypothesis.

The models raise experimental questions
Testing the model
First, the protein-protein reaction between FRQ and WCC plays a crucial role in the dynamics of the system. This reaction should be fast. Irreversibility is not necessarily needed, but multimerization should be greatly favored. For instance, for the parameters of Fig. 3, the dissociation rate of the complex must not be higher than 3.4 h^–1 to have sustained oscillations (data not shown). A possible experimental indication of this fact would be to measure the ratio between complexed and total proteins. For the protein with the lower concentration, this ratio should be close to 1.

Second, a consequence of this dimerization is the influence of frq transcripts on the dynamics of the system. The decay of messenger RNA plays a major role in the period determination. Raising or lowering the degradation time-constant of frq messenger RNA significantly changes the repression phase length and the qualitative behavior of the proteins. An alternative possibility is to modify the dynamics of binding of WCC to frq promoter. Mathematical analysis reveals that a lower detachment rate should lengthen the cycle if this rate is lower than the frq transcript's degradation rate, whereas a higher detachment rate should shorten the repression phase if this rate is higher than the frq transcript's degradation rate. The influence of seems difficult to test experimentally. However, the influence of transcript degradation rate could be easily tested, since it is, in principle, possible to alter the stability of frq transcripts by polyadenylation. One could thus test the correlation between the transcript's degradation rate and the period. One could also, for instance, imagine restoring the function of short period mutants (such as frq¹ or frq²; Feldman and Hoyle, 1973) by raising the transcript's stability. The present models also give an indication on the influence of FRQ degradation rate on the period of the clock: a lower degradation rate produces a longer period. This was already predicted and confirmed experimentally for the Goodwin oscillator applied to model the frq⁷ mutant (Ruoff et al., 1999; Liu et al., 2000). For the second two-loop model, dividing the degradation rate of FRQ dimers by 2 changes the period from 22 h to 29 h, as in frq⁷ strains, and frq transcript levels are also 32% higher than reference, qualitatively in agreement with experiments (Aronson et al., 1994) (data not shown).

Finally, to explain the phase shift between proteins, we proposed that only the monomer form of FRQ is active to enhance WC-1 translation. A consequence of this hypothesis is that in mutants where this homodimerization is switched off, constitutive levels of WC-1 should be very high, and even higher than the maximum level of WC-1 in normal cells. A possible test would be to vary FRQ levels and see that the enhancement of WC-1 translation does not vary monotonically with FRQ total concentration. For a low production rate, when FRQ mostly exists as a monomer, WC-1 translation rate should be high. For a high FRQ production rate, when FRQ mostly exists as a dimer, levels of WC-1 should be lower. This hypothesis could be tested experimentally. Strains have been artificially designed, where frq promoter is under the control of quinic acid (QA)(Aronson et al., 1994), and it is therefore possible to continuously vary FRQ protein production rate, while evaluating WC-1 concentration. One should observe a high WC-1 response only for a medium concentration of QA.

Improving the model
The measure of absolute concentrations would provide important data to refine the modeling. As some evolutions seem more or less linear, the present experimental data is not sufficient to fit the parameter values without this important information. This is also very important for understanding the mechanism of repression: if the repression mechanism is based on the heterodimerization which sequesters the WCC (Denault et al., 2001; Froehlich et al., 2003), the stoichiometry imposes constraints on the relative concentration of WCC and FRQ. In the present two-loop model, the WCC peak is approximately two-and-one-half times lower than FRQ peak. However, global extract (nucleus+cytosol) seems to show that FRQ and WCC peaks are approximately of the same order of magnitude (Denault et al., 2001).

To our knowledge, no nuclear extracts have been measured to evaluate this precise stoichiometry. Further hypotheses are therefore needed to explain this observed ratio. First, there could be specific different nuclear localization for the proteins, explaining a different ratio within the nucleus. Second, there could also be other negative feedback regulating the WCC level. Third, the stoichiometry