Pauline Traynard
email: pauline.traynard@gmail.com

Céline Feillet
email: celine.feillet@unice.fr

Sylvain Soliman
email: sylvain.soliman@inria.fr

Franck Delaunay
email: franck.delaunay@unice.fr

François Fages $
email: francois.fages@inria.fr

Model-based Investigation of the Circadian Clock and Cell Cycle Coupling in Mouse Embryonic Fibroblasts: Prediction of RevErb-α Up-Regulation during Mitosis

Keywords: Quantitative biology, Circadian clock, Cell cycle, Model coupling, Data fitting, Oscillations, Formal methods, Model Checking

Experimental observations have put in evidence autonomous self-sustained circadian oscillators in most mammalian cells, and proved the existence of molecular links between the circadian clock and the cell cycle. Some mathematical models have also been built to assess conditions of control of the cell cycle by the circadian clock. However, recent studies in individual NIH3T3 fibroblasts have shown an unexpected acceleration of the circadian clock together with the cell cycle when the culture medium is enriched with growth factors, and the absence of such acceleration in confluent cells. In order to explain these observations, we study a possible entrainment of the circadian clock by the cell cycle through a regulation of clock genes around the mitosis phase. We develop a computational model and a formal specification of the observed behavior to investigate the conditions of entrainment in period and phase. We show that either the selective activation of RevErb-α or the selective inhibition of Bmal1 transcription during the mitosis phase, allow us to fit the experimental data on both period and phase, while a uniform inhibition of transcription during

Introduction

In most organisms, from bacteria to plants and animals, spontaneous gene expression oscillations with a period close to 24 hours have been observed. A biochemical circadian clock present in each cell is responsible for maintaining these oscillations at this period, generally in the form of a self-sustained genetic oscillator entrained by the day/night cycle through various input pathways.

This circadian clock has many effects on cell signaling and metabolism [START_REF] Partch | Molecular architecture of the mammalian circadian clock[END_REF].

Experimental results have also shown a regulation of the cell division cycle by the circadian clock [START_REF] Matsuo | Control mechanism of the circadian clock for timing of cell division in vivo[END_REF][START_REF] Barnes | Requirement of mammalian timeless for circadian rhythmicity[END_REF][START_REF] Ünsal-Kaçmaz | Coupling of human circadian and cell cycles by the timeless protein[END_REF], in particular in mammalian cells with possible applications to cancer chronotherapies [START_REF] Ballesta | A combined experimental and mathematical approach for molecular-based optimization of irinotecan circadian delivery[END_REF][START_REF] De Maria | Design, optimization, and predictions of a coupled model of the cell cycle, circadian clock, dna repair system, irinotecan metabolism and exposure control under temporal logic constraints[END_REF]. Molecular links between these two cycles have been exhibited to explain this regulation. In particular the regulation of Wee1, an inhibitor of the G2/M transition, by the clock genes has been proposed to explain the circadian gating of mitosis during the liver regeneration process [START_REF] Matsuo | Control mechanism of the circadian clock for timing of cell division in vivo[END_REF] with 48 hours period doubling phenomena for the cell cycle [START_REF] Nagoshi | Circadian gene expression in individual fibroblasts: cell-autonomous and selfsustained oscillators pass time to daughter cells[END_REF]. Other similar molecular links going in the same direction, through p21 [START_REF] Gréchez-Cassiau | The circadian clock component bmal1 is a critical regulator of p21WAF1/CIP1 expression and hepatocyte proliferation[END_REF] and Chk1 and Chk2 [START_REF] Ünsal-Kaçmaz | Coupling of human circadian and cell cycles by the timeless protein[END_REF][START_REF] Gery | The circadian gene per1 plays an important role in cell growth and dna damage control in human cancer cells[END_REF], have been shown in different cells in the literature. A few models have also been developed to further investigate those hypotheses, by coupling a model of the cell cycle with a model of the circadian clock through those direct molecular links, and analyzing the conditions of entrainment in period [START_REF] Gérard | Entrainment of the mammalian cell cycle by the circadian clock: Modeling two coupled cellular rhythms[END_REF][START_REF] Calzone | Coupling the cell cycle and the circadian cycle[END_REF][START_REF] Glass | Synchronization and rhythmic processes in physiology[END_REF]. However, in mouse embryonic fibroblasts NIH3T3, several studies using large-scale time-lapse microscopy to monitor circadian gene expression and cell division events in real time and in individual cells during several days have unveiled unexpected behaviors, hinting that the relationship might be more complex. Nagoshi et al. [START_REF] Nagoshi | Circadian gene expression in individual fibroblasts: cell-autonomous and selfsustained oscillators pass time to daughter cells[END_REF], have first shown that circadian gene expression in fibroblasts continues during mitosis, but with a consistent pattern in circadian period variation relatively to the circadian phase at division, leading them to hypothesize that mitosis elicits phase shifts in circadian cycles. A more recent study of Bieler et al. [START_REF] Bieler | Robust synchronization of coupled circadian and cell cycle oscillators in single mammalian cells[END_REF] relating the same experiments on dividing fibroblasts found the two oscillators synchronized in 1:1 mode-locking leading the authors to hypothesize a predominant influence of the cell cycle on the circadian clock in NIH3T3 cells. This is in agreement with another detailed experimental study of Feillet et al. [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF] which found several synchronization states in NIH3T3 fibroblasts in different conditions of culture. In particular, it was observed in [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF] that enriching the culture medium with growth factors by increasing the concentration of Foetal Bovine Serum (FBS) not only accelerates the cell division cycle but also the circadian clock. For cells cultured in 10% FBS, both distributions of the cell cycle length and the circadian clock are centered around 22h. For cells cultured in 15 % FBS, both the cell cycle and the circadian clock accelerate, with period distributions centered around 19 hours. However, when cells reach confluence and stop dividing, the circadian clock slows down and the period distribution is then centered around 24 hours. None of the currently available models coupling the cell cycle and the circadian clock can explain these observations since they are based on an unidirectional influence of the circadian clock on the cell cycle [START_REF] Gérard | Entrainment of the mammalian cell cycle by the circadian clock: Modeling two coupled cellular rhythms[END_REF][START_REF] Calzone | Coupling the cell cycle and the circadian cycle[END_REF] and not on the other direction.

In this paper, in order to explain these observations, we investigate the reverse influence of the cell cycle on the circadian clock, using computational modeling tools. We develop a mathematical model of the influence of the cell cycle on the circadian clock through the differential regulation of clock genes around the mitosis phase, and study the conditions in which the cycles are entrained in period and phase as observed in [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF]. For this, we use the circadian clock model of Relogio et al. [START_REF] Relógio | Tuning the mammalian circadian clock: robust synergy of two loops[END_REF] which has been carefully fitted to phase data on suprachiasmatic cells, and a simple model of the cell cycle by Qu et al. [START_REF] Qu | Dynamics of the cell cycle: checkpoints, sizers, and timers[END_REF] which focuses on the mitosis phase. In [START_REF] Traynard | Model-based investigation of the effect of the cell cycle on the circadian clock through transcription inhibition during mitosis[END_REF], we have already shown that the uniform inhibition of transcription during mitosis, as observed in eukaryotes [START_REF] Weisenberger | A possible mechanism for the inhibition of ribosomal RNA gene transcription during mitosis[END_REF], could explain the acceleration of the circadian clock in non-confluent cells when the concentration of FBS increases [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF]. In particular, our model could reproduce the same periods for the cell cycle and the circadian clock for different levels of FBS, modeled by different values for the synthesis parameters of the cell cycle model, but with an incorrect time delay between the cell division and the peak of Reverb-α, which seemed impossible to fix under the hypothesis of a uniform inhibition of transcription during mitosis.

Here, we show that these difficulties can be resolved, using a different hypothesis of selective regulation of one clock gene during the M phase, either the activation of Reverb-α or the inhibition of Bmal1. Our coupled model under one of these hypotheses is able to reproduce the experimental measures on periods and phases made by Feillet et al. [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF] in individual unperturbed fibroblasts. Furthermore we argue that the complex behaviors observed with high variability after a treatment by dexamethasone to synchronize cellular clocks, modeled by the induction of a high level of Per and the inhibition of the other clock core genes, can be explained by the perturbation of the clock after this treatment. Indeed, our model shows that the stabilization time after that pulse appears to be greater than the time horizon of 72 hours used in those experiments. This computational model has been built using the Biocham modeling software [START_REF] Calzone | BIOCHAM: An environment for modeling biological systems and formalizing experimental knowledge[END_REF] for 1. importing and exporting models in SBML, and modeling the molecular interactions of the coupling of the models, 2. specifying the observed behavior in quantitative temporal logic using pattern formulae for periods and phases [START_REF] Fages | Temporal logic modeling of dynamical behaviors: First-order patterns and solvers[END_REF][START_REF] Traynard | Trace simplifications preserving temporal logic formulae with case study in a coupled model of the cell cycle and the circadian clock (best student paper award)[END_REF],

3. searching parameter values [START_REF] Rizk | Continuous valuations of temporal logic specifications with applications to parameter optimization and robustness measures[END_REF] and measuring robustness and parameter sensitivity indices [START_REF] Rizk | A general computational method for robustness analysis with applications to synthetic gene networks[END_REF] with respect to the temporal logic specification of the dynamical behavior1 .

Experimental Data and their Formal Specification in Temporal Logic

Experimental Observations and Measurements

In this section we explain the single cell experiments and analyses performed in [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF] and the conclusions drawn by the authors. The reported experiments have been done using time lapse videomicroscopy and cell tracking using different fluorescent reporters for the cell cycle and the circadian clock observed during 72 hours in proliferating NIH3T3 embryonic mouse fibroblasts. This cell line was modified to include three fluorescent markers of the circadian clock and the cell cycle: the RevErb-α::Venus clock gene reporter [START_REF] Nagoshi | Circadian gene expression in individual fibroblasts: cell-autonomous and selfsustained oscillators pass time to daughter cells[END_REF] for measuring the expression of the circadian protein2 RevErb-α, and the Fluorescence Ubiquitination Cell Cycle Indicators (FUCCI), Cdt1 and Geminin, two cell cycle proteins which accumulate during the G1 and S/G2/M phases respectively, for measuring the cell cycle phases [START_REF] Sakaue-Sawano | Visualizing spatiotemporal dynamics of multicellular cell-cycle progression[END_REF].

The cells were left to proliferate in regular medium supplemented with different concentrations of FBS (10% and 15%). Long-term recording was performed in constant conditions with one image taken every 15 minutes during 72 hours.

The lengths of the cell cycles were measured as the time interval between two consecutive cell divisions.

The expression traces of RevErb-α proteins were detrended and smoothed.

Spectrum resampling was used to estimate the clock period. Cells with less than two RevErb-α peaks within their lifetime, a period length outside the interval between 5 hours and 50 hours or a relative absolute error (RAE) bigger than 0.25 (showing a confidence interval wider than twice the estimated period) were classified as non-rhythmic and discarded, assuming that they do not have a functioning clock. Finally, the delay between cell division and the next clock marker peak was measured. It revealed that RevErb-α-Venus peaked about 7 h after cell division in all conditions, quite consistently with the delay of 5 hours for Reverb-α without Venus observed in [START_REF] Feillet | Coupling between the circadian clock and cell cycle oscillators: Implication for healthy cells and malignant growth[END_REF] and [START_REF] Bieler | Robust synchronization of coupled circadian and cell cycle oscillators in single mammalian cells[END_REF].

The quantitative data on the periods of the cell cycle and the circadian clock and the phase between them are summarised in 1: Estimated periods of the circadian molecular clock and the cell division cycle measured in [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF] in fibroblast cells without treatment by dexamethasone, for two concentrations of FBS. The time delay is between the cell division time and the next peak of RevErb-α protein.

Experimental Observations after Treatment by Dexamethasone

Furthermore, a series of experiments were done with a pulse of dexamethasone (abbreviated in the rest of this paper as Dex) before recording. This glucocorticoid agonist is known to exert a resetting/synchronizing effect on the circadian molecular clocks in cultured cells through the induction of Per genes.

In that case, the cells were incubated for 2 hours in the same medium supplemented with Dex, just before returning to a Dex-free medium for the recording.

The resulting dynamics in Dex pulsed cells are more complex. Table 2 summarizes the observed period and phase values reported in [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF]. The cells in 10% FBS show an increased clock period and a low cell cycle period, with an overall ratio of 5:4. In 20% FBS the cell lineages are dominated by two groups. The first group shows close periods, i.e. a 1:1 mode-locking similarly to the experiments without dexamethasone. The second group shows a high clock period and a fast cell cycle, with an overall ratio close to 3:2 between the clock and cell cycle, explaining the three-peaks distribution of the circadian phase at division, as already observed by Nagoshi et al. [START_REF] Nagoshi | Circadian gene expression in individual fibroblasts: cell-autonomous and selfsustained oscillators pass time to daughter cells[END_REF] ten years before. It has to be noted that the 20% FBS dexamethasone-synchronized experiment was repeated with similar results available in the Supplementary Information of [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF], although the distribution of the period ratios for the second group is wider in the interval ranging from 1.2 to 2. [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF] in fibroblast cells after treatment by dexamethasone, for two concentrations of FBS. The time delay is between the cell division and the next peak of RevErb-α protein. The experiment done with 20% FBS have been clustered by the authors of [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF] in two groups with different periods.

Medium

In [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF], the authors suggest that these observations might be interpreted by the existence of distinct oscillatory stable states coexisting in the cell populations, in particular with 5:4 and 1:1 phase-locking modes for the condition 10% FBS, and 3:2 and 1:1 phase-locking modes for the condition 20% FBS, and that the dexamethasone could knock the state out of the 1:1 mode toward other attractors.

Formal Specification of Oscillation Properties in Quantitative Temporal Logic

For the analysis of the dynamical behavior of this complex system, we shall make use of a temporal logic language which allows us to express the relevant system's oscillatory properties to fit, instead of over-specifying them by providing a precise curve to fit. This allows us to combine qualitative properties of oscillations and quantitative properties on the shapes of the traces such as distances between peaks or peak amplitudes. This is useful to capture the periods on either experimental and simulated traces, even when the traces are irregular and noisy. We use formal constraints on the amplitudes and regularity of the peaks to filter out ambiguous traces, keeping only sustained oscillations with small irregularities, as depicted for instance in Fig. 10.

More precisely, we use the temporal logic formula patterns described in [START_REF] Fages | Temporal logic modeling of dynamical behaviors: First-order patterns and solvers[END_REF][START_REF] Traynard | Trace simplifications preserving temporal logic formulae with case study in a coupled model of the cell cycle and the circadian clock (best student paper award)[END_REF] and implemented in our modeling software Biocham [START_REF] Calzone | BIOCHAM: An environment for modeling biological systems and formalizing experimental knowledge[END_REF] to specify the constraints about the successive peaks of concentrations between either the same molecular species (period constraints) or different molecular species (phase constraints). Biocham then provides commands for automatically

• extracting periods and phases from either simulation or experimental numerical data time series [START_REF] Fages | On temporal logic constraint solving for the analysis of numerical data time series[END_REF],

• searching the space of the unknown parameters of the model for satisfying period and phase constraints [START_REF] Rizk | Continuous valuations of temporal logic specifications with applications to parameter optimization and robustness measures[END_REF],

• measuring parameter sensivity indices and robustness with respect to period and phase constrainst [START_REF] Rizk | A general computational method for robustness analysis with applications to synthetic gene networks[END_REF].

For instance, the following command (used in Section 3) computes the validity domain of the variables of a formula pattern used to extract the period of the mitosis promoting factor (MPF) for the cell cycle, and of RevErb-α for the circadian clock, in a trace, together with their relative phase:

validity_domain( Exists([e1,e2,e3], periodErrors([RevErb_nucl],[periodselect,e1,e2,e3],100) & e1<2 & e2<2 & e3<2)) & period([MPF],[periodMPF]) & phase([MPF,RevErb_nucl],[phase])).
The result for the simulation trace displayed on the bottom of Fig. 6 is periodselect = 20.023, periodMPF = 20.046, phase = 6.765

Since the trace of MPF shows sustained and regular oscillations in all simulations, the simple predicate period is used here to extract the mean of the last two peak-to-peak intervals. The predicate phase similarly captures the mean of the last two time intervals between MPF and RevErb-α peaks. On the other hand, the period constraint on the oscillations of RevErb-α is expressed by the predicate periodErrors which deals with irregular traces. The validity domains of its variables provide the mean of the last two RevErb-α peak-to-peak intervals in the variable periodselect, along with several variables characterizing irregularity features of the trace, namely e1 for the irregularities in distances between peaks (it denotes the maximum difference between two intervals), e2 for the irregularities in the amplitudes of the peaks (it quantifies the differences between the amplitudes of the peaks), and e3 being a non-null error if the concentration amplitude is too small (below 0.1). Setting thresholds on these variables ensures that irregular traces are filtered out. Furthermore, a transient time of 100h is used to avoid the irregularities caused by the initial state. The logical quantifier

Exists projects the validity domain of the variables on the single dimension for

periodselect for the answer, and eliminates the other dimensions for e1, e2 and e3.

It is worth noting that for the purpose of parameter search, irregular traces after transient time should not be filtered out in order to orient the search algorithm in a promising direction, when oscillations begin to appear for instance.

The formula below specifies a constraint on period and phase for parameter search, without constraint on the period errors:

add_search_condition( periodErrors([RevErb_nucl],[period,e1,e2,e3],100) & Exists([phase], phase([MPF,RevErb_nucl],[phase]) & phase>minphase & phase<maxphase), [period,e1,e2,e3,minphase,maxphase],[21.3,0,0,0,3,5.5], 300,[(kdie,kdie21)]).
This formula is used by the parameter search procedure for testing the parameter values in different conditions. The condition expressed here, kdie=kdie21, fixes a particular value for the parameter kdie which determines the period of the cell cycle model (detailed in the next Section) in accordance with one of the experimental values for FBS. A time horizon of 300h is specified for the simulation.

The validity domains for the variables [period,e1,e2,e3,minphase,maxphase] are compared to the objective values [21.3,0,0,0,3,5.5] in order to determine a score for the constraint (as the Euclidean distance to the objective values). The addi-tional variables minphase and maxphase are used to score the distance from the phase to the objective interval [3,5.5], enabling some flexibility in the searched value for the phase. For instance, the best set of parameters found as a solution after the calibration procedure detailed in 4.2 has the score 0.54 for this constraint.. Furthermore, the difference between the periods of the cell cycle and the circadian clock is used to score the entrainment in period of the circadian clock by the cell cycle. This is achieved by the following formula with the variable diff :

satisfaction_degree( Exists([p1,p2], Exists([e1,e2,e3], periodErrors([RevErb_nucl],[p1,e1,e2,e3],100) & e1<3 & e2<3 & e3<3) & period([MPF],[p2]) & diff=p2-p1), [diff],[0],300).
This formal specification of entrainment in period is used to compute the satisfaction degrees displayed in Fig. 7. Each value scores the difference between the two periods with a continuous satisfaction degree between 0 (no entrainment) to 1 (null difference, perfect entrainment).

Cell Cycle and Circadian Clock Models

Model of the Cell Cycle

The cell cycle of mammalian cells is composed of five phases: the quiescent phase G0 where cells can stay without dividing, the growth phase G1 for entering the cell cycle, the DNA replication phase S, the gap phase G2, and the chromosome segregation and mitosis phase M phase. Each phase is characterized by a particular protein of the cyclin family, which forms a complex with a cyclin-dependent-kinase (CDK) and determines the activity of the phase. The mouse embryonic fibroblast cells considered in this paper are quickly dividing cells. However, those cells with growth factors also reach confluence and the G0 quiescent phase when they have no more space to divide, i.e. G0 by contact inhibition.

For our purpose of investigating the hypothesis of a regulation of some clock genes around mitosis, it is sufficient to use a cell cycle model focusing on the mitosis phase. We use a model proposed by Qu et al. [START_REF] Qu | Dynamics of the cell cycle: checkpoints, sizers, and timers[END_REF] in which the cell cycle is divided in two different phases: the G1/S/G2 phase and the M phase. Of course, more detailed models distinguishing the four phases of the cell cycle exist, e.g. [START_REF] Gérard | Entrainment of the mammalian cell cycle by the circadian clock: Modeling two coupled cellular rhythms[END_REF], making possible to represent various regulations of the cell cycle by the circadian clock genes, for instance through p21 and c-Myc on G1, and Wee1 on the G2/M transition. However, since the consequences of those regulations have not been observed in the experimental data considered in this paper, those extra details are irrelevant for our focus on the reverse effect of the cell cycle on the circadian clock by the regulation of clock genes around mitosis. The simpler two phase model of Qu et al. [START_REF] Qu | Dynamics of the cell cycle: checkpoints, sizers, and timers[END_REF] is thus sufficient to investigate this hypothesis.

In this model, the M phase is triggered by the complex Cdk1/Cyclin-B. This complex appears in two forms, an active form called MPF (M-phase Promoting Factor) and a phosphorylated, inactive form called preMPF. MPF is phosphorylated and inactivated by the kinase Wee1, and dephosphorylated and activated by the phosphatase Cdc25. Both the kinase and phosphatase activities are themselves regulated by MPF, respectively inactivated and activated by the complex, as depicted in Fig. 1. In this model, we assessed the effect of the different reaction rate constants on the period of the cell cycle using sensitivity analysis. We found that two parameters are able to change widely the range of the cell cycle period without changing significantly the strength of the coupling: kdie, the degradation rate constant of the intermediary enzyme involved in the negative feedback loop between MPF and the proteasome APC, particularly important in G1/S, and kampf, the activation rate constant of MPF by Cdc25p, which plays a role in G2/M. In the supplementary material of [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF], both the phases G1 and S/G2/M seem to be shortened in enriched FBS. Therefore there is no reason to prefer kdie (active in G1) or kampf (active in G2/M) to modulate the cell cycle period. We choose kdie as varying parameter affected by FBS because its activity is during the cell growing phase G1. The same results could be obtained with kampf.

MPF preMPF Cdc25P Cdc25 Wee1 Wee1P IE APCi APC
A simple parameter search gives the following values for kdie: 0.147 for a cell division period of 21.3 hours corresponding to 10% FBS, and 0.23 for a period of 18.6 hours corresponding to 15% FBS).

Model of the Circadian Clock

It has been shown in mice that in absence of synchronisation by a central clock, autonomous circadian oscillators are maintained in peripheral tissues, although they are progressively desynchronized [START_REF] Yoo | LUCIFERASE real-time reporting of circadian dynamics reveals persistent circadian oscillations in mouse peripheral tissues[END_REF] because each one of them has its own period, that varies slightly from 24h. This has been confirmed in cultured NIH3T3 cells first in [START_REF] Nagoshi | Circadian gene expression in individual fibroblasts: cell-autonomous and selfsustained oscillators pass time to daughter cells[END_REF] and then in [START_REF] Bieler | Robust synchronization of coupled circadian and cell cycle oscillators in single mammalian cells[END_REF] and [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF]. In each of those studies, confluent fibroblasts have a circadian clock period close to 24 hours regardless of the medium concentration.

In this paper we use the circadian clock model of Relogio et al. It is worth noting that this way of coupling the models enforces the fact that for quiescent cells, whatever the FBS concentration is, the transcription rate is unaffected by mitosis. As a consequence, the clock necessarily returns to a period of 24 hours in confluent cells, as observed in the experiments.

Coupling Parameters

In our coupled model, the influence of the cell cycle on the circadian clock is modeled by a set of parameters that express the regulation coefficient (activation or inhibition) of clock gene synthesis during mitosis. We use the parameter search procedure of Biocham based on stochastic optimization [START_REF] Rizk | Continuous valuations of temporal logic specifications with applications to parameter optimization and robustness measures[END_REF] to find the sets of values for the coupling parameters that reproduce the entrainment in period and phase observed in the data. Using the behavior specification detailed in Section 2.3, we define a multi-condition objective for the four conditions kdie=0.1, kdie=0.147, kdie=0.18 and kdie=0.23: the period of the circadian clock must be equal to the period of the cell cycle, respectively 24 hours, 21.3 hours, 20h and 18.6 hours. Furthermore, in each condition, the delay between MPF and RevErb-α peaks must be between 6.5 hours and 8.6 hours. The simulation of the model with any of these two sets of parameters shows a delay between the starting time of the mitosis effect and the circadian clock consistent with the experimental data (close to 7 hours for the first set, and between 7 and 8.5 hours for the second one), as well as consistent period values.

Coupling parameters

The stochastic optimization procedure for parameter search returns numerical values for all coupling parameters, however one must check which parameter values are necessary for the correct entrainment, and which parameters have no impact on the satisfaction of the specification of the behavior. To this end, we perform a sensitivity analysis by computing the response curves for the period of the circadian clock and the delay between mitosis and the next RevErb-α peak by varying each coupling parameter, in the condition where the cell cycle has a period of 21 hours.

The results for the first set of parameters are shown in Fig. 4. These simulation results reveal that the entrainment in period and phase of the circadian clock depends only on the effect of mitosis on RevErb-α and Bmal1, and on the duration of this effect. In the first parameter set, the activation of RevErb-α is crucial for the entrainment in period. On the other hand, varying the synthesis coefficients during mitosis for Per, Cry or Ror has no significant effect on the entrainment in period and phase. More specifically, the clock is entrained to the cell cycle period of 21 hours if the coefficient multiplied to the synthesis rate of RevErb-α is at least 1.7, and the inhibition lasts at least 2 hours. In this condition the time delay between divisions and RevErb-α peaks is consistent with the data, with RevErb-α peaks occuring 7 hours after division. An inhibition of RevErb-α would also preserve the entrainment in period, but the circadian phase at division would become inconsistent with the data, with RevErb-α peaks occurring just after mitosis. One can notice that an activation on Bmal1 would also increase the time delay of RevErb-α peaks after division.

The results for the second set of parameters are shown in Fig. 5. These simulation results reveal similarly that the entrainment in period and phase of the circadian clock depends only on the effect of mitosis on Bmal1 and RevErbα, and not on Per, Cry or Ror. However in this case the entrainment in period depends on the effect on Bmal1 which has to be inhibited with a coefficient smaller than 0.2 in order for the circadian clock to be entrained at 21 hours with a correct delay after division.

Simulation Results for the Selective Activation of RevErb-α during Mitosis

The first set of coupling parameter values points toward the hypothesis of a selective activation of the circadian clock gene RevErb-α triggered by mitosis.

Thanks to the previous sensitivity analysis, the first coupling parameter set can be simplified to consider a single activation of RevErb-α (corresponding to a coefficient equal to 2), and no effect on the other clock genes.

It is worth noting that an up-regulation of RevErb-α translation during the M phase has already been observed in a genome-wide experimental study interrogating the translational lanscape during cell cycle progression using ribosome profiling [START_REF] Stumpf | The translational landscape of the mammalian cell cycle[END_REF]. This provides one first argument for considering this hypothesis in the first place. Table 4 shows the periods of the circadian clock and the cell division cycle and the delay between the starting time of the inhibition of RevErb-α, when the peak of MPF overtakes the threshold 0.5, and the following peak of RevErb-α in our model with different values of kdie corresponding to the different culture conditions. As shown in the right panel of Fig. 2, it is possible to simulate the experimental medium enrichment with 10 or 15% FBS by varying the parameter kdie of the cell cycle model to obtain the same values for the period of the cell division cycle. In all cases, the cell division manages to entrain the circadian clock (that has a free period around 24 hours) to its period, simply through this mechanism of selective transcription activation, as depicted in Fig. 6. These simulation results reproduce quite well the data of Table 1 when there is no treatment by Dex. Interestingly, our model can also have a cell division time higher than 24 hours, for instance with kdie=0.077 which might correspond to a concentration of FBS around 5%. In that case the model predicts that the cell cycle will still entrain the circadian clock, lowering its period. Moreover, RevErb-α peaks are predicted to stay just after the mitosis. The landscape in Fig. 7 is computed to assess the role of the inhibition or activation duration. It shows the variation of the difference between the periods of RevErb-α for the circadian clock and MPF for the cell cycle when the two parameters kdie and duration vary. The value of each period is captured with a temporal logic specification as seen in the subsection 2.3. The result for the activation of RevErb-α. is shown in Fig. 7. Three domains can be distinguished in this parameter space. In the domain in yellow, the circadian clock is entrained to the same period as the cell cycle. This domain of entrainment is wider for a long duration of activation. For a short duration, the circadian clock can only be entrained by the cell cycle if the entraining period is close to 24 hours (corresponding to kdie= 0.1). In the purple domain at the bottom (for a low value of duration), the difference between the two periods is high because the clock is not entrained, hence it keeps its period constant and close to 24 hours.

Finally, these two domains are separated by a black domain where the clock oscillations are partially entrained and become irregular.

One can notice that the longer the activation of RevErb-α, the wider the range of values of kdie over which the circadian clock can be entrained. In particular, the clock can be entrained by the cell cycle when kdie = 0.23, corresponding to the smallest period (18 hours) reported in the data, if the duration of the activation is at least 3 hours.

The entrainment both in period and phase with an activating effect during 3 hours is visualized in Fig. 8. It shows the response curve for the periods of the cell cycle and the circadian clock, and the time delay between the peaks of MPF and RevErb-α in the simulations when the parameter kdie varies.

Predictions on the phases in the clock

The experimental data on the phases between clock components in proliferating cells are sparse. The model allows us to investigate whether the coupling from the cell cycle affects the phases between clock mRNAs and proteins. The following in silico experiment has been performed to this end: in the coupled model with a fast cell cycle (21 hours), the strength of the activation of RevErbα is changed in a set of simulations and the phases between clock components are captured in each simulation, normalised by the clock period. The results are shown in Fig. 9. The simulations reveal that in the entrained condition (when the synthesis rate of RevErb-α is activated with a coefficient higher than 1.75), the phases between clock components are impacted by the periodic activation of RevErb-α resulting from the coupling with the cell cycle, compared to their values in the free clock (when the synthesis rate of RevErb-α is activated with a coefficient close to 1). In particular, the phase between Bmal1 and RevErb-α mRNAs shows a small advance, that impacts similarly the phases between Bmal1 and RevErb-α in the cytoplasm and in the nucleus. Notably, the other clock genes and proteins targeted by Bmal1 exhibit a phase delay when the synthesis of RevErb-α is activated during mitosis.

• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • 0.
• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •
The prediction is thus that in dividing cells, the phases between the clock proteins slightly but significantly differ from the phases in quiescent cells.

Comparison to Experimental Data after Treatment by Dexamethasone

In order to take into account the experiments with dexamethasone, the model can be extended with an event, lasting for two hours, and inducing Per mRNA while inhibiting the other clock genes.

Fig. 10 shows that in our models, regardless of the growth factors in the medium (i.e. of the value of kdie), the Dex pulse results in a perturbation of the clock and then returns to the observed entrainment. The pulse alters the clock before returning to the previously observed entrainment regime. In the left panel the pulse is from time 162 to 164 while on the right it is from 170 to 172. The left panel's peak-to-peak distance is in the [18.8; 22.7] interval, while the right one remains in the [20.9; 21.7] interval. This might correspond to the two groups observed in [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF]. The time to recover normal entrainment varies but is often larger than 72 hours.

These simulations point us to the possibility that the noisy data reported in Table 1 after the Dex pulse might simply be due to the various cellular states in which the pulse happened and to the time necessary for the cells to recover their clock entrainment, rather than to two different oscillatory attractors of the system. A pulse at time 170h disrupted only slightly our clock, leading to mostly remaining in mode-locking 1 : 1, whereas advancing that same pulse by 8 hours (corresponding to giving the pulse to a cell in a different state) leads to a bigger disturbance, some peak-to-peak distances close to 23 hours, others to 18 hours, and even if this is transitory, this might correspond to the type of data observed in the Group 2 of Table 1.

One can remark that an alternative hypothesis could be that Dex causes a stress response linked to the cell cycle and the circadian clock. This would still look like noise, but there could be a mechanistic explanation for this effect which could be investigated with a modeling approach. However, to our knowledge, a short treatment by Dex is not known to cause a stress response.

Simulation Results for the Selective Inhibition of Bmal1 during Mitosis

The second solution found by the calibration procedure in section 4.2 is discussed here. Like the activation of RevErb-α transcription, the inhibition of Bmal1 during 3 hours, triggered at mitosis, is able to entrain the circadian clock in a wide range of periods (18-28 hours). The results summarized in Table 5 show a similarly good fit to the experimental data on periods and phases. The resulting traces are displayed in Fig. 11 for different cell cycle length conditions.

The entrainment in phase, ie. the duration between mitosis and the following RevErb-α peak, differs between the two couplings, as shown in Fig. 12. In both cases, when the cell cycle has a period smaller than 24 hours, the circadian clock marker RevErb-α peaks 6 hours to 8 hours after mitosis. When the cell cycle period is greater than 24 hours, a notable difference can be seen for the predicted phase: with an activation of RevErb-α, RevErb-α still peaks just after the mitosis. But with an inhibition of Bmal1, RevErb-α peaks 18 to 24 fours after the mitosis. No experimental observations exist in a slowed down cell cycle condition, but [START_REF] Nagoshi | Circadian gene expression in individual fibroblasts: cell-autonomous and selfsustained oscillators pass time to daughter cells[END_REF] and [START_REF] Bieler | Robust synchronization of coupled circadian and cell cycle oscillators in single mammalian cells[END_REF] report some cells dividing not long after the circadian peak where the circadian clock was found to be slowed down. On the other hand, the large translational landscape of gene expression during cell-cycle progression reported in [START_REF] Stumpf | The translational landscape of the mammalian cell cycle[END_REF] does not consider Bmal1 nor any clock gene other than RevErb-α.

Uniform Inhibition of Transcription during Mitosis

It has been shown that in eukaryotes, gene transcription can be significantly inhibited during mitosis [START_REF] Weisenberger | A possible mechanism for the inhibition of ribosomal RNA gene transcription during mitosis[END_REF]. The impact of a global transcription inhibition of clock genes during mitosis on the circadian oscillator has been studied by modeling in [START_REF] Kang | Modeling the effects of cell cycle m-phase transcriptional inhibition on circadian oscillation[END_REF]. In this study, the authors found that a periodic inhibition of transcription during one hour was able to entrain a model of the mammalian circadian clock, but only when the inhibition period was close to one half, twice or equal to the intrinsic circadian model period. In these cases, a phase locking between the circadian clock and the periodic inhibition was observed, albeit with one or two preferential circadian phases for the inhibition and values that varied greatly with the inhibition period. The discrepancies with the recent data could come from the short inhibition duration considered or from the arbitrarily parameterised model used for the circadian clock, taken from [START_REF] Leloup | Toward a detailed computational model for the mammalian circadian clock[END_REF].

In [START_REF] Traynard | Model-based investigation of the effect of the cell cycle on the circadian clock through transcription inhibition during mitosis[END_REF], we also investigated the uniform inhibition of all clock genes and found that it was sufficient to reproduce the entrainment of the circadian clock by the cell cycle in period, but not in phase. The delay between MPF and RevErb-α remained inconsistent with the data, as depicted in Fig. 13,i.e. the mitosis triggered by MPF occurs just after the peaks of RevErb-α, while the experimental studies consistently report the opposite: peaks of RevErb-α 5-7 hours after divisions. Furthermore, it seemed impossible to find parameter values to reproduce the observed delay under that hypothesis of a uniform inhibition of transcription during mitosis, which thus cannot explain the experimental data in mouse embryonic fibroblasts NIH3T3 measured in [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF].

Conclusion

By hypothesizing a selective activation of RevErb-α or a selective inhibition of Bmal1 triggered at the beginning of mitosis, we have been able to build a mechanistic dynamical model which reproduces the somewhat surprising numerical data reported in [START_REF] Bieler | Robust synchronization of coupled circadian and cell cycle oscillators in single mammalian cells[END_REF][START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF] about the acceleration of the circadian clock observed in dividing fibroblasts with high FBS concentrations. These observations suggest that the primary coupling between the cell division cycle and the circadian clock results from an influence of the cell cycle on the circadian clock in those cells. While considering a uniform inhibition of the transcription during mitosis [START_REF] Weisenberger | A possible mechanism for the inhibition of ribosomal RNA gene transcription during mitosis[END_REF] was shown to be sufficient to fit the period data in [START_REF] Traynard | Model-based investigation of the effect of the cell cycle on the circadian clock through transcription inhibition during mitosis[END_REF],

the phase data reported in [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF] seemed to be impossible to reproduce under that uniform inhibition hypothesis. The use of Biocham search algorithms for computing transcription inhibition parameters satisfying the period and phase observations formalized in quantitative temporal logic, led us to the hypothesis that around mitosis, either the transcription of RevErb-α has to be strongly activated, or Bmal1 strongly inhibited with no inhibition of RevErb-α. These two hypotheses differ by their predictions on slow cell cycle cells, possibly obtained with low levels of FBS, but for which no quantitative data are currently available.

Since RevErb-α is a repressor of Bmal1 and Bmal1 an activator of RevErbα, these two hypotheses correspond to two alternative mechanisms for a similar effect. Interestingly, an experimental study investigating the translational landscape of dividing cells reported an up-regulation of RevErb-α translation during the M phase [START_REF] Stumpf | The translational landscape of the mammalian cell cycle[END_REF].

This is one reason to favour this hypothesis. Furthermore, the activation of RevErb-α might be caused by various mechanisms [START_REF] Loyer | Role of CD-K/cyclin complexes in transcription and RNA splicing[END_REF]. For instance, the transcription factor c-Myc displays bursts of transcriptional activity during G1 phase just after mitosis and the S to G2/M transition of the cell cycle [START_REF] Seth | Cell cycle regulation of the c-Myc transcriptional activation domain[END_REF].

The c-Myc protein regulates its target genes through the same E-box DNA response element as the Clock/Bmal1 heterodimer. It is therefore conceivable that during the G2/M phase of the cell cycle, RevErb-α is positively regulated by c-Myc leading to the transcriptional repression of Bmal1 [START_REF] Altman | MYC disrupts the circadian clock and metabolism in cancer cells[END_REF]. In this scenario the E-box regulated Per and Cry genes are expected to be also up-regulated by the higher transcriptional activity of c-Myc. This is compatible with our simulations which show that the phase and period are resilient to variation of the coefficient synthesis for Per and Cry. An additional potential link between RevErb-α and the cell cycle machinery is provided by the recent observation that the Cyclin B-CDK complex can phosphorylate RevErb-α [START_REF] Zhao | Circadian Amplitude Regulation via FBXW7-Targeted REV-ERBα Degradation[END_REF].

Our model also postulates a different interpretation of some of the results presented in [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF] when cells are treated by a 2 hours pulse of dexamethasone.

Indeed, instead of different autonomous cycling regimes, the model predicts temporary perturbations leading to shorter or longer peak-to-peak distances, but returning to the previous entrainment regime after some time, longer than the horizon used in the experiments.

Furthermore, in our coupled model, the phases between some of the clock gene products are shifted when entrained by a fast cell cycle. We are able to quantify these phase shifts and show that they concern mainly RevErb-α whose mRNA peaks are advanced by the periodic activation during mitosis. Other clock mRNAs and proteins are slightly delayed compared to their activator Bmal1. A prediction of the model is therefore that in quickly dividing cells, these protein peaks are shifted with respect to quiescent cells where such a phenomenon should not be observed.

Figure 1 :

 1 Figure 1: Schema of the cell cycle model of Qu et al. [16]. Solid arrows represent biochemical reactions while dashed arrows denote enzyme catalysis. Red arrows denote an inhibition through activation of the degradation.

Figure 2 :

 2 Figure 2: Left: Simulation of the cell division cycle model of Qu et al. Right: Period of the cell division cycle (measured as the distance between successive peaks of MPF) as a function of the parameter kdie for MPF activation by Cdc25p in the model of Qu et al.

Figure 3 :

 3 Figure 3: Simulation trace of the concentrations of the clock gene products in the nucleus over a time horizon of 100h in the model of Relogio et al. [15].

Figure 4 :

 4 Figure 4: Period of the circadian clock (red) and time delay between MPF and RevErb-α proteins peaks (blue) in the coupled model calibrated with the first set of parameters, and when the cell cycle has a period of 21 hours. The peaks on the blue curves on the right figures characterize irregular oscillating traces (due to a partial entrainment in period).

Figure 5 :

 5 Figure 5: Period of the circadian clock (red) and phase between the division and RevErb-α protein (blue) in the coupled model calibrated with the second set of parameters, and when the cell cycle has a period of 21 hours.

Figure 6 :

 6 Figure 6: Simulation of the model with the activation of the transcription of RevErb-α triggered by mitosis, during 72 hours. Top: the cell cycle has a period of 21.3 hours. Middle: the cell cycle has a period of 20.1 hours. Bottom: the cell cycle has a period of 26 hours.

Figure 7 :

 7 Figure7: Absolute difference between the periods of the circadian clock and the cell cycle, as a function of kdie for varying the cell cycle period, duration, the duration of the activation of RevErb-α transcription triggered by mitosis. The landscape is computed as the satisfaction degree of the third formula detailed in 2.3. The color translates the distance from the value found for the period difference diff to the objective 0. Full satisfaction in yellow indicates equal periods for MPF and RevErb-α proteins, while the other colours indicate the absolute difference. Black indicates an absence of result for the specification, meaning that the regularity constraints set on the trace of RevErb-α with the function periodErrors were not met.

Figure 8 :

 8 Figure 8: Entrainment in period and phase of the circadian clock when the period of the cell cycle varies with the parameter kdie, with the activation of RevErb-α triggered by mitosis. The blue curve depicts the time delay between peaks of MPF and RevErb-α in the simulations.

Figure 9 :

 9 Figure 9: Circadian clock phases in the coupled model when the cell cycle has a period of 21 hours and the activation of RevErb-α varies. Note that the deep variations between synthesis rates 1.25 and 1.75 can be explained by irregular oscillations caused by partial entrainment, on which measured phases are not reliable.

Figure 10 :

 10 Figure10: Effect of a dexamethasone pulse on the entrainment resulting from the periodic activation of RevErb-α synthesis by the cell cycle. The pulse alters the clock before returning to the previously observed entrainment regime. In the left panel the pulse is from time 162 to 164 while on the right it is from 170 to 172. The left panel's peak-to-peak distance is in the[18.8; 22.7] interval, while the right one remains in the[20.9; 21.7] interval. This might correspond to the two groups observed in[START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF]. The time to recover normal entrainment varies but is often larger than 72 hours.

Figure 11 :

 11 Figure 11: Simulation of the model with the the periodic transcription of Bmal1 activated by mitosis, during 72 hours. Top: the cell cycle has a period of 21.3 hours. Middle: the cell cycle has a period of 20.1 hours. Bottom: the cell cycle has a period of 26 hours.

  0

Figure 12 :

 12 Figure 12: Entrainment in period and phase of the circadian clock when the period of the cell cycle varies with the parameter kdie, with the inhibition of Bmal1 (left) or the activation of RevErb-α (right) triggered by mitosis. In the right panel, the circadian clock period is missing for low kdie values because the oscillations are irregular. The blue curves depict the time delay between peaks of MPF and RevErb-α in the simulations.

Figure 13 :

 13 Figure 13: Entrainment to a period around 21.3 hours with kampf = 3.75 corresponding to a medium enriched with 10% FBS, found in [17].

Table 1 [

 1 [START_REF] Feillet | Phase locking and multiple oscillating attractors for the coupled mammalian clock and cell cycle[END_REF]. Surprisingly, increasing FBS from 10% to 15%, not only decreases the mean period of the cell cycle from 21.3 hours to 18.6 hours, but also the clock period from 21.9 hours to 19.4 hours, i.e. to essentially the same period. This shows that both oscillators remain unexpectedly in 1:1 mode locking. While the speedup of the cell cycle can be directly attributed to the growth factors in increasing concentration of FBS, it can not account for the speedup of the clock the same way, since Table

	confluent cells keep a 24-hours period for the circadian clock independently of
	the FBS concentration [14].
	Medium	Clock period Division period Mean delay
	FBS 10% 21.9 h ± 1.1 h 21.3 h ± 1.3 h	8.6 h
	FBS 15% 19.4 h ± 0.5 h 18.6 h ± 0.6 h	7.1 h

Table 2 :

 2 Estimated periods of the circadian molecular clock and the cell division cycle measured in

Table 3 :

 3 Two sets of synthesis regulation coefficients during mitosis, found by Biocham's calibration procedure for satisfying the temporal logic specification of period and phase data. The first set was found with initial values equal to 1 for all synthesis coefficients, while the second was found with null initial values.Per transcription. Cry and Ror transcriptions are activated while the RevErb-α transcription is mostly unaffected.

	Essentially, two types of solutions, shown in Table 3, are found. If the initial
	values are 1 for all synthesis coefficients, corresponding to no inhibition during
	mitosis, the best result found after 55 iterations on a population of 95 sets of

Table 4 :

 4 Periods and time delays measured in the coupled model with different values of kdie for modeling the different culture conditions (the correspondance with 5% FBS is speculative since no significant experiment was done in this condition). The delays are the time observed by simulation between the peaks of concentration of MPF and RevErb-α.

Table 5 :

 5 Periods and time delays reproduced by the coupled model with different values of kdie for modeling the different culture conditions (the correspondance with 5% FBS is speculative since no experiment was done in this condition). The delays are the time observed by simulation between the peaks of concentration of MPF and RevErb-α.

The models and the formal specifications used in this paper are available on http:// lifeware.inria.fr/wiki/software/biosystems16.

In this paper, the genes are distinguished from the proteins by writing the names of the genes in italics and the proteins in normal text.

Acknowledgements.

This work has been partially funded by the ANR HyClock project (ANR-14-CE09-0011), ANR "Investments for the future" LABEX SIGNALIFE program (ANR-11-LABX-0028-01) and the former ERASysBio+ C5Sys project (ANR grant 2009-SYSB-002-02). This work was also granted access to the HPC resources of the CINES under the allocation c2015036437 made by GENCI.


F Ladieu

J Le Cochec

P Pari

P Trouslard

P Ailloud

Dielectric Constant of Glasses: Evidence for Dipole-Dipole Interactions

Keywords: numbers: 61.43.Fs, 72.20.Ht, 72.20.My, 77.22.Gm

The 1 kHz real part 0 of the dielectric constant of a structural glass (a-SiO 2x C 1y H z ) was measured at low temperature T. Reducing the sample thickness h below 100 nm weakens the slope j@ 0 =@Tj for T & 0:1 K, for all measuring fields E. This contrasts with the predictions of the two-level system (TLS) model but is in agreement with the recently proposed delocalization of excitations derived from a field-induced TLS-TLS interaction mechanism. For small h this interaction is screened, which explains the h effects on 0 . Hence, interactions must play a key role in standard thick samples as soon as T & 0:1 K.

For more than 30 years, amorphous solids have been widely investigated at low temperature T [1]. Their properties turn out to be both ''universal'' (i.e., quasiindependent from their chemical composition) and strongly different from their crystalline counterparts. This is explained within the tunneling two-level system (TLS) model [START_REF] Anderson | [END_REF] where some atomic systems fluctuate between two neighboring energy minima separated by a potential barrier, which, at low T, is crossed by tunneling.

Even if the low density of glasses (in comparison with their crystalline counterpart) may justify [3] the existence of ''voids'' and thus the TLS framework, the ability of the standard noninteracting TLS model to account for experiments is very intriguing. Indeed, by extracting from measurements the coupling strength of a TLS to phonons, one finds [4] that TLS's are strongly coupled via virtual phonons. A tunnel transition on a given TLS deforms elastically the neighboring matrix, thus inducing an energy change U / 1=r 3 for a TLS located at a distance r. For two neighboring TLS's [5] one gets Ur 1 nm ' 10 K. Similarly, since many TLS's are charged, they interact through dipolar interaction, mediated by virtual photons, giving an interaction energy [5] in the same range of 10 K between neighboring TLS's. This large energy scale contrasts with the (very few) experimental results that suggest weak interactions between TLS's. In the 100 mK range, only the small nonequilibrium effects of the kHz dielectric susceptibility [6], as well as its small unexpected H dependence [7], have been related to interactions. In the mK range, somewhat larger effects, such as the ultralow T plateau [8] of the dielectric constant, and the internal friction behavior [9], were ascribed to TLS interactions.

This work provides strong evidence of the key role played by TLS's interactions in the 100 mK range. We show that reducing the thickness h of our glassy samples below 100 nm strongly affects the real part 0 of the 1 kHz dielectric constant. The 0 T minimum, which occurs at a given T rev for h > 100 nm, progressively shifts to lower T as h is decreased (and finally disappears for h ' 12 nm). We show that this behavior contradicts the noninteracting TLS model and is most likely related to TLS interactions. Many theories dealing with interactions in the TLS model have been proposed. In Ref. [10], it is suggested that interactions are renormalized to zero by frustration, whereas in Ref. [11] interactions are assumed to supersede disorder. More recently Burin et al. [12] proposed a mechanism that accounts for the reported data, once properly factored in the numerical determination of 0 [13].

The samples were all produced similarly: on a vitreous a-SiO 2 0.1 mm thick substrate, a layered Cu=glass= Cu=glass structure was deposited, by using a multichamber system excluding exposure to air during the whole process [14]. The 15 nm thick Cu electrodes were evaporated (at 0:1 nm=s). The glass layers of a-SiO 2x C 1y H z (with x ' 2:9, y ' 0:9, z ' 1:2) were made from TetraEthylOrthoSilane with a 13 MHz vapor plasma, and x-ray analysis confirmed their amorphous nature [START_REF]X-ray experiments[END_REF]. The autopolarization was set to ÿ100 V, the incoming flux to 2 sccm, and the pressure to 0.80 Pa. The resulting glass deposition rate was 0:1 nm=s, allowing one to select the thickness h of the internal glass layer by choosing the deposition time. Since the value of h is crucial, it was independently measured by three other methods that all gave compatible results: (i) in situ laser interferometry was realized onto a Si substrate placed close to the sample during the glass deposition; (ii) the glass layer grown onto the Si substrate was irradiated by a deuteron 0.91 MeV beam allowing, using the nuclear reaction on 16 O, an h estimate; (iii) the value of the capacitance C / 0 was checked to scale with the expected h. The 15 nm thick top glass layer hinders any spurious atmospheric effect during the cryogenic experiment. The sample was glued inside a copper box connected to the mixing chamber of the 3 He= 4 He dilution refrigerator. Semirigid coaxial wiring was ensured from the cold copper box up to the 2500-Andeen capacitance bridge. A capacitance c F , twice as large as that of the sample, was set in parallel of each cable, so as to filter high frequency parasitic fields. The data remained unchanged when c F was halved, proving the filtering efficiency. The 0 T; E behavior of six samples, with h ranging from 12 to 800 nm, is reported in Figs. 123, where CT CT ÿ C0:5 K and C max C0:5 K ÿ C min , while C min is the minimum value of C in our T experimental range (14 mK; 0.5 K). Figure 1 shows that, for a given measuring field E 40 kV=m, the temperature T rev , where 0 T is minimum, gradually shifts towards lower values as h is decreased. For the thinnest sample, the minimum is suppressed for fields lower than 250 kV=m. At E 40 kV=m, the simplest way to describe the CT behavior of the thinnest sample is CT / T , where is equal to 1:35 0:15 above 70 mK and to 2:75 0:4 below. This strongly contrasts with the usual logarithmic behavior CT / lnT, valid for thicker samples as well as for all previously studied glasses [START_REF] Rogge | [END_REF]. Figure 2 shows that the lowering of T rev when decreasing h occurs for all E (except for E 1 MV=m where some heating arises; see below). However, this decrease occurs only for h & 100 nm since the 200 and 800 nm samples behave similarly. Finally, Fig. 3 shows, for the smallest h, that a minimum of 0 T is recovered when E * 1 MV=m, but the value of the slope j@ 0 =@Tj below T rev is smaller than that of thick samples.

P H Y S I C A L R E V I E W L E T T E R S
The above properties of 0 T; E; h are intrinsic to the glassy state. Indeed, as we show below, they come neither from a heating effect nor from a variation of the glass composition with h:

(i) The E-induced dissipated power P might heat the sample to a temperature above the measured T. This explanation is ruled out since the data of Fig. 1 are taken at constant E. Note that one has P / E 2 h; i.e., P increases with h at a given E. Thus, if the measurements were affected by heating, the observed trend would be the opposite of that of Figs. 1 and2. The presence of heating can be suspected only for the two thickest samples and for E * 1 MV=m, where T rev is slightly lower at h 200 nm than at h 75 nm. To check this point, the thermal resistance R between the sample and the cryostat was estimated as follows: first for E > 10 kV=m, the law T rev / E 1=2 was assumed to hold (see below) for a perfectly well thermalized sample; second, the (small) difference between this law and the measured T rev in the range 10-300 kV=m was attributed to a (small) heating effect. From the two thickest samples one gets R ' 50 MK=W at 100 mK: this is the good order of magnitude to account for the above-mentioned behavior of T rev when E * 1 MV=m and h 200 nm. Finally, one finds R / T ÿ3 , as expected for boundary Kapitza resistances.

(ii) Besides, one may imagine that the effect of Figs. 1-3 comes from the fact that the glass composition depends on h, even for a given set of plasma parameters. The components of the plasma may be mixed with the Cu of the first electrode and yield a ''boundary'' layer, of thickness b 0 , with a chemical composition that strongly differs from the rest of the ''normal'' glass (of thickness h ÿ b 0 ). Since the value of T rev slightly depends on the glass composition [START_REF] Rogge | [END_REF], one may wonder if the data could be explained by combining the susceptibilities of these two dielectrics in series. This is not the case. Indeed, taking for the boundary layer the 12 nm sample, and for the normal glass the 800 nm sample, one finds, for the 22 nm thick sample, that CT should, by far, be dominated by the ''thin'' T behavior, for any E. This is clearly contradicted by Fig. 1 where the standard lnT trend holds for the 22 nm thick sample. Even if the 12 nm thick sample is discarded because of its unusual T behavior, the scenario of two dielectrics in series can be further tested by using the 22 nm sample as the new ''boundary layer.'' For any given E, the resulting T rev for the 75 nm thick sample should be only 5 mK below that of the thickest sample. This is clearly contradicted by the data of the 100 kV=m range of Fig. 2. At this step, it is clear that decreasing h below 100 nm qualitatively changes the physics of the glass sample. This contradicts the standard (i.e., noninteracting) TLS picture. This model accounts for the decrease of 0 above T rev by the progressive freezing of the diagonal (or relaxational) part 0 z of the susceptibility, whereas the increase of 0 below T rev comes from the off-diagonal (or resonant) part 0

x (see Ref. [START_REF]The susceptibility 0 is given by 0 0 z 0 x ; see Ref[END_REF]). Because of its purely quantum nature, 0

x grows when T decreases, as do all quantum effects. In other words, according to the standard model, 0 is dominated by the dynamics of individual TLS's (see Ref. [START_REF]More precisely, in the standard TLS model, TLS-TLS interactions are assumed to be[END_REF]). Since the TLS's size is below 1 nm, no h effect can be explained in the range of tens of nm, contrary to the data reported above.

This failure of the standard TLS model strongly contrasts with its former success, in particular, its ability to account for the 0 / lnT behavior as well as for the order of magnitude of T rev [1,[START_REF] Rogge | [END_REF]. In fact, it has been very recently discovered that the agreement holds only in the ''linear regime'' (i.e., small E's) where 0 , and thus T rev , do not depend on E. Indeed, as shown in Ref. [13], the standard TLS model does not account for the nonlinear behavior of 0 T; E (such as that reported in Ref. [START_REF] Rogge | [END_REF] and that of our h > 100 nm samples). Indeed, due to the quantum nature of 0 below T rev , 0 is strongly depressed by a strong measuring electric field E at a given T. This is due to the quantum saturation phenomenon coming from the fact that increasing E decreases the population difference between the two energy levels. As the Rabi oscillations produced by E on the upper level are in phase opposition with respect to those produced on the ground level, the quantum response, once averaged on many independent TLS's, tends to zero when E is increased. This was checked by solving numerically the Bloch equations of TLS's [13] with a nonperturbative method (see the dotted lines in Fig. 4).

Both the nonlinear measurements and the h effect on 0 can be accounted for by using the field-induced TLS-TLS interaction proposed by Burin et al. [12]. The idea is to consider a thermally excited TLS of gap and to seek for the probability that such an excitation delocalizes by hopping to another TLS of gap 0 located at a distance from the first one. The probability is large only if j ÿ 0 j U, where U is the interaction coupling the two TLS's. With the usual density P P(independent on ) of gaps in glasses, one finds that is extremely small at E 0, but that a finite E strongly favors the delocalization of excitations. This is due to the fact that E modulates any by an amount e dip pE= r , where p ' 1 D is the TLS dipole and r ' 5 is the glass dielectric constant [6,[START_REF] Rogge | [END_REF]. Thus, all the TLS's whose gaps differ from each FIG. 4. 0 T; E was computed in Ref. [13] by solving Bloch equations. The noninteracting TLS model is shown in dotted lines: due to the quantum saturation effect, 0 / C is strongly depressed when raising E at low T. This suppresses the minimum of 0 T for experimental fields E > 40 kV=m and could model the 12 nm sample where TLS's interactions are strongly reduced. Including TLS's interactions through the E-induced mechanism of Ref. [12] restores the 0 increase with E (at any given T), as well as that of T rev . This models the usual 0 T; E of thick samples. Finally, since strong E's decrease the distance between interacting TLS's, one recovers the usual behavior of 0 T when E < h. This explains, for h 12 nm, the appearance of a 0 T minimum for strong enough E, as in Figs. 2 and3. other by less than e dip form an ''infinite cluster'' within which the delocalization of excitations occurs.

As far as the nonlinear behavior of 0 is concerned, this mechanism leads to the following predictions: (i) Energy diffusion enhances the diagonal part 0 z of the susceptibility, by decreasing the associated relaxation time 1 . This is of great qualitative interest since, as stated above, 0

x is mainly depressed by E and thus cannot account for the increase of 0 with E at low T. (ii) The rate of delocalization (i.e., 1= 1 ) increases with E, and the associated increase of 0 z with E overcomes the E-induced depression of 0

x . (iii) The rate 1= 1 grows when T decreases, since quantum coherence is needed for the delocalization of excitations. Figure 4, drawn from Ref. [13], shows that the mechanism of Burin et al. allows one to recover a nonlinear 0 T; E behavior in qualitative agreement with that of our two thickest samples (or that of Ref. [START_REF] Rogge | [END_REF]). Along the same lines [13], one predicts T rev / E p , as seen in Fig. 2 in the nonlinear regime (E 10 kV=m).

The delocalization of excitations also accounts for the thickness effects reported here. Indeed, dipolar interactions are strongly screened by the electrodes [19], i.e., by the fact that their numerous electrons intercept and cancel the electric field responsible for the interaction between TLS's. This weakens U for distances * h. Since the distance P Pe dip ÿ1=3 between two neighboring TLS's of the infinite cluster is typically 60 nm for E 40 kV=m, the mechanism of Burin et al. is expected to be progressively switched off as h decreases in this range. Since this mechanism is responsible for the behavior of 0 T T rev , one understands that T rev should decrease with h, as seen on Fig. 2. As an extreme consequence, just as in Fig. 1, T rev could even disappear for the thinnest sample due to the complete suppression of the delocalization of excitations (and to the E-induced depression of 0 x valid at any h).

Besides, due to the E-induced modulation of , in the scenario of Burin et al., decreases for strong E, opening the possibility to recover < h for large enough E. This could account for Fig. 3 where, even for h 12 nm, one recovers a minimum for 0 T.

In conclusion, decreasing the thickness of glass samples below 100 nm changes the physics of the real part of the kHz dielectric susceptibility. This experimental result strongly supports the key role of dipole-dipole interactions up to 100 mK. The main features of the reported data can be accounted for by assuming that TLS's interactions occur via an E-induced delocalization of excitations between quasisimilar TLS's.

2 is 1 K, which is much smaller than any relevant energy scale. That is why, according to the standard model,'' 0 depends only on what happens on the scale of individual TLS's.' ' [19] The distance between the two TLS's turns from / e dip ÿ1=3 , as derived in the text for < h, to / e dip ÿ1=2 when > h. This ''simple geometric effect'' somewhat enlarges and thus reduces both the elastic and the dipolar parts of U. It is not easy to evaluate the importance of the ''geometric effect'' with respect to that of screening by the electrodes. Whatever it may be, this does not change our conclusion that the h effects reported here come from a reduction of TLS-TLS interactions.

FIG. 1 .

 1 FIG.1. Relative variation with T of the capacitance, for a given measuring field E 40 kV=m. As the sample thickness h decreases the minimum shifts towards lower T and disappears for the smallest h. This contradicts the standard noninteracting TLS model. Inset: same plot (T in linear scale); the 24 nm curve is vertically moved ( 0:2) for clarity. Note that for h 12 nm, the standard lnT behavior is not obeyed.

FIG. 3 .

 3 FIG.3. Relative variation of the capacitance with T, for various E and for h 12 nm. No minimum occurs in 0 T down to 14 mK for the four smallest E's, but a minimum reappears for larger E's. At E ' 0:4 MV=m (gray symbols) one may have T rev 19 mK. The vertical scale is reduced compared to that of Fig.1, since the variation of 0 T < T rev is weaker than in thick samples.

We thank M. Ocio and E. Vincent for lending us their cryostat and L. Le Pape for experimental support. We especially thank J.-Y. Prieur and J. Joffrin (CNRS, Orsay) for very encouraging scientific discussions, D. L'Ho ˆte and A. Forget for experimental help, and J. P. Bouchaud for helping us to write this Letter.