Constraint Identiﬁcation Using Modiﬁed Hoare Logic on Hybrid Models of Gene Networks

We present a new hybrid Hoare logic dedicated for a class of linear hybrid automata well suited to model gene regulatory networks. These automata rely on Thomas’ discrete framework in which qualitative parameters have been replaced by continuous parameters called celerities. The identiﬁcation of these parameters remains one of the keypoints of the modelling process, and is diﬃcult especially because the modelling framework is based on a continuous time. We introduce Hoare triples which handle biological traces and pre/post-conditions. Observed chronometrical biological traces play the role of an imperative program for classical Hoare logic and our hybrid Hoare logic, deﬁned by inference rules, is proved to be sound. Furthermore, we present a weakest precondition calculus (a la Dijkstra) which leads to constraints on dynamical parameters. Finally, we illustrate our “constraints generator” with a simpliﬁed circadian clock model describing the rhythmicity of cells in mammals on a 24-hour period.


Introduction
Formal methods from computer science have been largely applied to model and analyse biological systems [5,17].In particular, verification tools like model-checking have been used to verify dynamical properties of discrete models [3,7] in which the temporal aspects are only present through the succession of events: the delay between two successive events is not taken into account.Unfortunately, continuous time turns out to be important in most biological systems, in particular for gene regulatory networks.Gene regulatory networks are models representing influences between genes leading to the modification of the synthesis of associated proteins.Because proteins can regulate their target genes, positive or negative feedback loops emerge making possible a large variety of behaviours.These gene regulatory networks are designed to apprehend and predict effects of a component on the system but such models are also useful to confront hypotheses with the up-to-date collected knowledge on the gene interactions.
Several modelling framework have been devoted to gene networks.These frameworks differ by the aspects they highlight.Stochastic models emphasize non-determinism, differential models represent a system with a lot of details (transcription, traduction, transports...) [14] and give precise trajectories in terms of concentrations; qualitative models [16,15] focus on the major features that explain the observations (only main causalities are taken into account); and hybrid models link qualitative aspects with continuous variables such as time.
Whatever the modelling framework, the main difficulty of building such networks is the identification of parameters governing the dynamics of the system.The determination of these parameters is crucial to observe a behaviour consistent with biological knowledge.We already showed [3] that formal methods can help in this parameter identification step in the René Thomas' discrete modelling framework [16].Unfortunately, this framework abstracts temporal information, often necessary for tuning the models.This discrete approach is based on the splitting of homogeneous concentration areas which have the same effects on other components.In order to refine this complete discrete framework, we associated with each such area a celerity describing the evolution speed of each component.These parameters lead to a class of linear hybrid automata.We also showed how the known experimental traces can be used to establish constraints on dynamical parameters (celerities) and to restrain the set of admissible parameters [1].
Numerous works already focus on the study of hybrid automata [9].Especially, several tools already exist to tackle the model checking of linear hybrid automata, either with classical exhaustive approaches [10] or using abstract interpretation [8].Communicating sequential processes (CSP), a process algebra for describing patterns of interaction in concurrent systems, has also been extended to hybrid systems and hybrid Hoare logic has been proved to be useful in such context [18].These methods among others propose parameter synthesis in some ways.Nevertheless, these tools are tailored for a general purpose and will not take into account the specificities coming from biological contexts.
We divert Hoare logic (whose aim is to rigorously reason about the correctness of imperative programs) in order to determine constraints on celerities in such models.This approach was already developed, for the discrete framework [2,4] and we extend it in the present paper for hybrid automata.Hoare logic has already been extended to real time systems [11] in which continuous evolutions are not taken into account whereas they are important in our biological context.Hoare logic relies on triples of the form {Pre} p {Post} where Pre and Post are conditions on states of the system and p is a path of the system.A Hoare triple is considered true if, whenever the system is reset at a state satisfying the Pre condition, the path p is possible and leads the system into a state which satisfies the condition Post.Following E. Dijkstra [6], the aim of the game is then to determine, for each path p and postcondition Post, the weakest precondition Pre, thus covering the largest set of states, making the Hoare triple true.In our temporal approach, the time spent in each state becomes crucial to determine the constraints on parameters.
We illustrate our formalism with a biological process named circadian clock which synchronises all cells in mammals at a 24-hour rhythmicity.We design a hybrid automaton modelling this biological cycle and, from the observed trace of this process, we build constraints for each celerity of this hybrid automaton using the aforementioned hybrid Hoare logic.We finally show that simulations, run parameters values satisfying the constraints, exhibit curves which are similar to experimental data.
The paper is organised as follows.We first define in Section 2 the formalism of the hybrid modelling framework.Then Section 3 focuses on the syntax and the semantics of the modified Hoare logic, and the weakest precondition calculus, whereas Section 4 is devoted to the soundness of our hybrid Hoare logic.We illustrate in Section 5 the use of this Hoare logic for identifying the constraints on parameters of the simplified circadian clock model.Finally in Section 6, we discuss the limits of this approach and we expose some possible extensions.

Hybrid Modelling Framework of Gene Network
A gene network is visualised as a labelled directed graph (interaction graph) in which vertices are either variables (within circles) or multiplexes (within rectangles), see Fig. 1.Variables abstract genes and their products, and multiplexes contain formulas that encode situations in which a variable or a group of variables (inputs of multiplexes, dashed arrows) influences the evolution of some other variables (output of multiplexes, plain arrows).A multiplex can encode the formation of molecular complexes, phosphorylation by a protein, competition of entities for the activation of a promoter, etc. Definition 1 gives the formal details of a gene network.
Definition 1 (Hybrid Gene Regulatory Network).A hybrid gene regulatory network (GRN for short) is a tuple R = (V, M, E, C) where: V is a finite set whose elements are called variables of the network.With each variable M is a finite set whose elements are called multiplexes.With each multiplex m ∈ M is associated a formula ϕ m in the multiplex language formed of the atoms "(v n)", where n ∈ 1, b v 1 , and the usual logical connectives ¬, ∧, ∨ and ⇒.E is a set of edges of the form (m is called the celerity of v for ω at the level n and these celerities have to satisfy the following constraints: Let us remark that the dashed arrows of Fig. 1 are not present in the previous definition.When representing a gene network, it is convenient to visualise the variables contributing to a particular multiplex, but from a formal point of view, this information is redundant with the formula of the considered multiplex. Celerities (noted C v,ω,n ) give the evolution of each variable v when it is under the active regulation of the set ω of its predecessors and when it is in the qualitative state n.They code for the dynamics of the system and we aim at the identification of these celerities.The constraints on celerities given in the previous definition link the signs of celerities to the underlying dynamics and may need some explanations.The first one deals with the case where a celerity value is null for a given set of active predecessors ω of a variable v.This models an equilibrium state, thus the related constraint states that celerities around this

5:4
Constraint Identification Using Modified Hoare Logic equilibrium, for the same set of active predecessors ω, must force v to reach this equilibrium.As a consequence, there is a single null celerity at most for a given couple of v and ω.If, on the other hand, no celerity is null for these v and ω, a consequence of the second constraint is that they are all of the same sign.This models that v either always decreases until full degradation or always increases up to saturation.
In the remainder of this section, we focus on the dynamics of a gene network.Definition 2 introduces the hybrid states whereas Def. 3 explains the crucial notion of resources of a variable in a particular state.
π is a function from V to the interval [0, 1] of real numbers.η is called the discrete state or qualitative state of h whereas π is called its fractional part.For simplicity, we note in the sequel η v = η(v) and π v = π(v).We denote S the set of hybrid states of R. When there is no ambiguity, we often use η and π without explicit mention of h.Definition 3 (Resources).Let R = (V, M, E, C) be a GRN and let v ∈ V .The satisfaction relation h ϕ, where h = (η, π) is a hybrid state of R and ϕ is a formula of the multiplex language, is inductively defined as follows: If The usual meaning of the logical connectives is used.The set of resources of a variable v at a state h is defined by: ρ that is, the multiplexes that are predecessors of v and whose formula is satisfied.We note that the set ρ(h, v) only depends on the discrete state of h: all hybrid states having the same discrete part thus have the same resources.Indeed, inside a discrete state η, the dynamics of v is controlled in the same manner, thus the celerity is the same for all states h = (η, π), that is: C v,ρ(h,v),ηv .By abuse of language, we also use the notation ρ(η, v).From this celerity, and given a particular hybrid state, one can compute the touch delay (Def.4) of each variable, which measures the time necessary for the variable to reach a border of the discrete state.Definition 4 (Touch Delay).Let R = (V, M, E, C) be a GRN, v be a variable of V and h = (η, π) be a hybrid state.We note δ h (v) the touch delay of v in h for reaching the border of the discrete state.More precisely, δ h is the function from V to R + ∪ {+∞} defined by: Nevertheless, reaching the border of a discrete state is not sufficient to go beyond the frontier: there may be no other qualitative level beyond the border (we call such a border an external wall) or the celerity in the neighbour state may be of the opposite sign (internal wall), as given in Def. 5.

Definition 5 (External and Internal Walls)
. Let R = (V, M, E, C) be a GRN, let v ∈ V be a variable and h = (η, π) a hybrid state.1. v is said to face an external wall at state h if: Figure 2 Continuous transitions.Inside each state, a continuous transition (h0 → h 0 ) goes from the initial point h0 to the unique point h 0 from which a discrete transition takes place (h 0 → h1).Left: The celerity vector allows, without sliding mode, the trajectory to directly reach a border which is crossed.Center: From h 0 two possible discrete transitions can occur: h 0 → h1 or h 0 → h2.Moreover (πg, πpc) corresponds to the fractionnal coordinates of a hybrid state h along the path.Right: The grey zone depicts an external or internal wall.The only possible discrete transition is h 0 → h1.
where sgn is the classical sign function.We note sv(h) the set of sliding variables, that is, variables that face an internal or external wall in the qualitative state of h.
We note that a sliding variable v ∈ sv(h) may not be actually sliding if it has not reached its border yet.However, if in addition δ h (v) = 0, then v is located on an internal wall or external wall.In this case, its fractional part cannot evolve anymore in the current qualitative state (see Fig. 2-Right where variable g reaches an external wall).We introduce the notion of first changing variables in Def.6 which are the first variables able to change their discrete levels.
Definition 6 (First Changing Variables).Let R = (V, M, E, C) be a GRN and h = (η, π) be a hybrid state.The set of first changing variables is defined by: Moreover, δ first h denotes the time spent in the qualitative state of h when starting from h: The set first(h) represents the set of variables whose qualitative coordinates can change first.If a variable is on an external or internal wall, it cannot evolve as long as other variables do not change, thus: first(h) ∩ sv(h) = ∅.Similarly, if the celerity of v in the current state is null, its qualitative value cannot change because of an infinite touch delay.
Figure 2 illustrates several evolutions of a gene network.From a particular hybrid state h 0 , the dynamics alternates continuous transitions (within the discrete state) and discrete transitions (when changing the discrete state).When the trajectory reaches an external or internal wall (see Fig. 2-Right), the variable slides along the wall only if the celerity of some other variable can drive the trajectory in such a direction.This description leads to the following definition.
Definition 7 (Hybrid State Space).Let R = (V, M, E, C) be a GRN.We note R = (S, cT, dT ) the hybrid state space of R where S is the set of hybrid states, and cT (resp.dT ) is the set of continuous (resp.discrete) transitions:

5:6 Constraint Identification Using Modified Hoare Logic
1.There exists a continuous transition in cT from state h = (η, π) to state h = (η , π ) iff: a.Either first(h) = ∅ and there exists a variable v ∈ first(h) such that: i. δ h (v) = 0, where δ h (v) is called the duration of the (continuous) transition, Or first(h) = ∅ (meaning that each variable v either reaches an equilibrium state: C v,ρ(h,v),ηv = 0; or faces a wall: v ∈ sv(h)) and: There exists a discrete transition in dT from state h = (η , π ) to state h = (η , π ) iff there exists a variable v ∈ first(h ) such that: The states from which there do not exist any transitions (discrete or continuous) are called steady states.
The continuous transitions lead to the last hybrid state inside the current discrete state, at which point a qualitative change can happen.The instantaneous discrete transitions make the system evolve, as soon as the system can (that is, when δ h (v) = 0), into the next qualitative state by going through a border.These two different kinds of transitions can be observed on Fig. 3 where the discrete transitions are in dotted lines and the continuous transitions are in plain lines.Let us remark that there is a unique continuous transition starting at a given hybrid state.Indeed, assuming that there exist two continuous transitions h → h 1 and h → h 2 from the same hybrid state h, the item 1 of the previous definition leads to the equality h 1 = h 2 regarding the ends of the continuous transitions (whatever the value of the set first(h)).
Let us notice that the defined linear hybrid automata leads to an undeterministic behaviour: when the celerity vector allows the trajectory to reach more than one border at the same time, several discrete transitions can be considered (see Fig. 2-Center).Some of these discrete transitions can be forbidden in case of internal wall.

Hybrid Hoare Logic
This section is dedicated to the presentation of the Hoare logic adapted to our hybrid formalism.Hoare logic is based on Hoare triples noted {Pre} p {Post} meaning that if a program p is executed from a state satisfying a precondition Pre, then after execution, the postcondition Post is true.In our case, the program p is replaced by a biological trace characterising biological knowledge on the chronometrical qualitative behaviour of the system.For this, we define the property language used for pre-and postconditions in Subsec.3.1 and the path language used to describe observed traces in Subsec.3.2.Then, Hoare logic is defined using these languages and we give in Subsec.3.3 an adaptation of the weakest precondition calculus, that is, the computation of the weakest (the most general) precondition that makes the trace possible and such that the postcondition Post is satisfied afterwards.
In the rest of this section, we denote by any of the usual comparison symbols on integers or real numbers: <, ≤, >, ≥, =, =.

Property Language
We first define the property language describing pre-and postconditions.Definition 8 (Property Language L C ).The terms of the property language L C are inductively defined as follows: A discrete term is a variable η u with u ∈ V , or a constant of N; and n ∈ 0, b u , or a constant of R; The connectives +, −, × and / create new terms by composition, the latter being only valid for continuous terms.We use their usual semantics.
Discrete atoms are of the form n n where n and n are discrete terms and continuous atoms are of the form f f where f and f are continuous terms.The discrete conditions are defined by: The hybrid conditions are defined by: where a d and a c are respectively a discrete atom and a continuous one.
A property is a couple (D, H) formed by a discrete and a hybrid condition.All such couples (D, H) form the property language L C .
A hybrid state h satisfies a property ϕ = (D, H) ∈ L C iff both D and H hold in h, by using the usual meaning of the connectives; in this case, we note h ϕ.

Path Language
The path language given in Def. 12 takes the role of an imperative program in a Hoare triple by describing a biological behaviour.Such a path consists in explicit discrete transitions as given in Def. 9, but also in continuous transitions described by duration and some information, see Def. 10.The characterisation of continuous transitions is based on two kinds of atoms: C v c constrains the value of the current celerity of v, and slide(v) constrains v to slide.
Definition 9 (Discrete Path Atom).The (discrete) path atoms are defined by: For any states h = (η, π) and h = (η , π ), the transition In the following, if v ∈ V is a variable, v± refers indistinctly to v+ or v−.Definition 10 (Assertion Language L A ).The assertion language L A is defined by the following grammar: A couple (∆t, a) ∈ R + × L A of a non-negative real number and an element of the assertion language is called an assertion couple.
The following definition gives the semantics of such assertion couples.From an informal point of view, for any states h = (η, π) and h = (η, π ) in the same qualitative state η, the continuous transition h → h satisfies the assertion couple (∆t, a) if the continuous transition exists and if it lasts ∆t units of time and it respects the assertion a: is always true; C v c T I M E 2 0 1 7

5:8
Constraint Identification Using Modified Hoare Logic ,ηv is the celerity of v in the current qualitative state; slide + (v) (resp.slide − (v)) is satisfied iff v faces and reaches a wall at the top of the domain (resp.at the bottom); slide(v) is a shorthand for slide + (v) ∨ slide − (v); and logical connectives have their usual meanings.We note indifferently (h, h ) (∆t, a) or h (∆t,a) −→ h .Regarding Def. 10, the special case where ∆t equals 0 characterises a situation where the system enters a qualitative state in a "corner" and no continuous transition is required between two successive discrete transitions.
Definition 11 (Semantics of the Assertion Couple (∆t, a)).Let us consider a hybrid state h = (η, π) and the unique continuous transition starting from h and ending in h = (η, π ).The satisfaction relation between the continuous transition h −→ h and an assertion couple (∆t, a) ∈ R + × L A is noted (h, h ) (∆t, a), by overloading of notation, and is defined as follows: If The path language allows the modeller to express experimental biological traces as sequences of elementary paths.The next section shows how such information can be formally taken into account in order to help the identification of celerities compatible with such paths.

Hoare Triples and Weakest Precondition
We are now able to give the definition (Def.13) of a Hoare triple in the scope of our hybrid formalism which is a natural extension of the classical definition.Figure 3 gives an example of a valid Hoare triple.
Definition 13 (Hybrid Hoare Triples).A Hoare triple for a given GRN is an expression of the form {Pre} p {Post} where Pre and Post, called precondition and postcondition respectively, are properties of L C , and p is a path from L P .A Hoare triple {Pre} p {Post} is satisfied iff for all state h 1 Pre, there exists another state h 2 so that h 1 p −→ h 2 and h 2 Post.Now that the semantics of this new Hoare logic is defined, we aim at adapting the weakest precondition calculus as proposed by Dijkstra [6] to our hybrid framework (Def.14).Edsger Dijkstra introduced a predicate transformer semantics: the semantics of an imperative programming language is defined by assigning to each instruction in this language a corresponding predicate transformer.For each elementary instruction EI of the imperative programming language, the weakest precondition of EI is a function mapping any postcondition Post to a precondition Pre.Actually, this function returns the weakest precondition on the initial state ensuring that the execution of EI terminates in a final state satisfying Post.For each sequential imperative program "P ; EI " whose last instruction is EI , and for each postcondition Post, the predicate transformer of EI allows us to first determine the weakest precondition just before the last instruction and by iterating the same process, it becomes possible to determine the weakest precondition of the whole imperative program (loops are treated in a particular way with the help of invariants).
In our setting, the same approach leads to build the minimal constraints on the celerities insuring that starting from a state satisfying the precondition Pre, the model exhibits the known path p (corresponding to a biological trace) leading to a state satisfying the postcondition Post.Each constraint depends on each elementary path which is defined by the time ∆t spent in the current qualitative state, the assertion a and the discrete path atom v±.Each elementary path takes the role of an elementary instruction.Definition 14 (Weakest Precondition).Let p be a path program and Post = (D, H f ) be a post-condition parameterized by a final state index f .The weakest precondition attributed to p and Post is a property: ), parameterized by a fresh initial state index i and the same final state f , and whose value is recursively defined by: If p = ε is the empty sequence program, then D ≡ D and which is parameterized by a fresh intermediate state index m; T I M E 2 0 1 7

5:10
Constraint Identification Using Modified Hoare Logic ∆t), A(∆t, a) and J v are sub-properties detailed in Appendix A.
We note that in the cases corresponding to atoms p = (∆t, a, v±), the formula H i,f contains H f with substitutions, in conjunction with Φ ± v (∆t), ¬W ± v , F(∆t), A(∆t, a) and J v , which makes the weakest precondition of a sequence of instructions very big and difficult to compute or analyse by hand.Nevertheless, each of the previous subformulas corresponds to a condition which has to be met to allow the execution of an atomic instruction (∆t, a, v±): The sign of the celerity of v in the current state is given by v±; Traversing the qualitative state lasts ∆t units of time (Φ ± v (∆t)); There is no internal or external wall preventing v to increase or decrease qualitative state (¬W ± v ); All components other than first changing variables must either reach their border after v, or face an internal or external wall (F(∆t)).The assertion a is verified along the continuous transition (A(∆t, a)); The continuous transition inside a discrete state links the fractional parts of v, its celerity and time spent in the current discrete state.Similarly the discrete transition indicates that the fractional parts of states before and after a discrete transition are the same except for the variable v changing its discrete level (J v ).
Finally, the computation of the weakest precondition for a given Hoare triple {Pre} p {Post} is automated using the classical backward proof strategy: If p is of the form (∆t, a, v±) or ε, then we compute the precondition.If p = p 1 ; p 2 with p 2 = (∆t, a, v±), we compute the precondition before p 2 and we iterate for path p 1 (we never consider p 2 as ε).These two items are mutually exclusive which means that the proof strategy generates a unique proof tree.
An implementation of this weakest precondition calculus has been realised2 .Section 5 details its result on a model of the circadian clock and before that, next section gives the theorem of its soundness.

4
Soundness of the Hybrid Hoare Logic

Inference Rules and Axioms
The considered Hoare logic for hybrid gene regulatory networks is defined by the following inference rules: where v is a variable, η v its expression level, D (resp.H) the discrete (resp.hybrid) condition,

5:11
both detailed in Appendix A, ∆t the time spent inside the current discrete state and a an assertion.The last inference rule is the sequential composition rule: Sequential composition rule: where Q 1 , Q 2 , Q 3 are properties of the form (D, H) having the role of pre-and postconditions, and p 1 and p 2 are particular paths deduced from biological experiments.
The two following axioms, based on the semantics of the hybrid model, complement the inference rules: the expression level has to be in its definition domain), C v,ω,ηv × C v,ω,ηv+1 ≥ 0 (for two neighbour qualitative states, if the variable v is controlled by the same resources, then the celerities of v cannot be of opposite signs).

Soundness of the Hoare logic
The following lemmas are useful for the proof of soundness.Lemma 15 states that the time spent in the current discrete state is equal to the time mandatory, for the variable which changes first, to reach its border.Lemma 16 expresses the fact that the truth value of a formula remains the same after a continuous transition.
Lemma 15 (Time Spent in a Discrete State).Let h be a hybrid state.If h (D 1 , H 1 ) and From the definition of the sub-property Φ + v (∆t), see Appendix A.2, variable v reaches its upper border (π i v = 1) and its celerity is positive (C v,ω,n > 0).Let h = (η, π ) be the hybrid state where v first touches this border.The time spent in the current qualitative state η corresponds to the time necessary to reach the border where v changes its qualitative level.Indeed, from Φ + v (∆t) we deduce: Cv,ω,n .From Def. 4, we have ∆t = δ h (v).Finally, the sub-property F(∆t), see Appendix A.4, expresses the fact that if a variable u different from v reaches its border before v, u faces an internal or external wall (see Appendix A.3). Thus, since δ h (v) = ∆t and according to Def. 6, we deduce

Proof.
Since h and h belong to the same discrete state, the expression levels of all variables are the same.The evaluation of D in h is then the evaluation of D in h.The atoms of H i,f (see subformulas in Appendix A) concern either celerities or discrete or continuous coordinates of different points of the trajectory (entrance and arrival points in different discrete states).These points are either outside the current discrete state, or are the points h or h .Moreover, the celerities are constants, the points of the trajectory which do not belong to the current discrete state as well as the points h and h do not change.The interpretation of (D , H i,f ) is therefore the same in h and h .

Constraint Identification Using Modified Hoare Logic
The soundness of the modified Hoare logic, adapted for the hybrid modelling framework, means that if a Hoare triple is built in agreement with the inference rules (Def.14) then this Hoare triple is satisfied according to the semantics of Hoare triples (Def.13).
Theorem 17.The hybrid Hoare logic is sound.
The proof is detailed in Appendix B.

Example: Simplified Circadian Cycle
The circadian rhythm is a biological process regulating cells of an organism with a 24-hour period and controlling the electrical and metabolic processes.

Presentation of the Circadian Cycle
In mammals, the main circadian cycle is located in the suprachiasmatic nucleus and regulates the peripheral clocks.It is affected by light, acting like a synchronizer called Zeitgeber, which means "giver of time".The circadian rhythm is mainly controlled by two protein complexes which are PER/CRY and BMAL1/CLOCK.When light appears, the BMAL1/CLOCK complex activates the per and cry genes by binding the E-box response element in the promoter upstream these genes [12].The PER and CRY proteins are synthesized and dimerised in the cytoplasm.During night, this complex is found inside the nucleus and inhibits BMAL1/CLOCK implying a negative feedback of PER/CRY on its genes.Finally, PER/CRY is degraded by proteasome.The circadian cycle completes and a new one begins with the transcription of genes bmal1 and clock.
We decided to use an interaction graph which focuses on the per and cry components, as represented in Fig. 1.The Light/Day cycle (whose duration is 12h/12h) is represented by the node named L. (Let us notice that the node labelled X is a modelling artefact to get an oscillating feature for light.)This node enhances the per and cry genes (modelled with node g) when the light is activated, that is, when the qualitative value of L is at level 1.These genes synthesize their proteins which complex and spread inside the nucleus.When the complex is activated (which is modelled by an expression level of 1 for pc), those genes are inhibited, blocking the synthesis during night.Because genes are disabled, the protein complex will be degraded by proteasome and a new cycle begins with the reactivation of the genes.All four nodes of this model have two qualitative levels of expression named 0 (not active) and 1 (active).

Hoare Triple and Results
The steps of the circadian clock explained in Subsec.5.1 are represented in the Hoare Triple below.The time spent in each qualitative state comes from biological information obtained during a light/day cycle of 12h/12h.The assertions slide + (g) and slide + (pc) (resp.slide − (L)) characterize a saturation (resp.complete degradation) of g and pc (resp. of L, corresponding to the beginning of the night).This information is summed up in the following Hoare triple:

Constraints on celerities of g and pc
Constraints on celerities of L and X Using this Hoare triple, we compute via the backward strategy the weakest precondition iteratively by crossing the intermediate states of the path.In addition, because the behaviour is cyclic, we know that the starting hybrid state of this path is equal to the finishing one.The provided implementation allows us to identify and simplify the constraints of the celerities obtained through our weakest precondition calculus.After some automatic and manual simplifications, we finally obtain the constraints summed up in Table 1 which make the known cyclic behaviour possible.
In order to illustrate the validity of this process, we used the constraint solver Ibex 3 to extract values satisfying the previous constraints.This constraint solver takes as input only conjunctions.Thus the obtained constraints are transformed in a disjunctive normal form (we obtained 3 terms in disjunction) and all terms of this disjunction are successively given to Ibex to extract possible values for variables.Amongst the values returned by Ibex, we arbitrarily chose one set of possible values to be injected in the model for simulation.Interestingly, the obtained constraints fully characterize one of the hybrid states along the limit cycle, which gives us an initial state for the simulation.
The simulated traces have then to be compared to biological experimental data.Because such data for the PER/CRY protein complex inside the nucleus are not published, the simulation (Fig. 4) is compared to experimental data of genes per1, per2, cry1 and cry2 as well as their respective proteins taken separately [13].We noticed that the maximal activity of per genes of our simulation and experimental data occurs at the end of a day, and the curves of proteins are maximal during night at the same time slot for both curves.Thus our simplified circadian cycle model is consistent with biological experimental data although we arbitrarily used parameter values satisfying constraints.This simulation reinforces reliability of our formalism for determining constraints.

Conclusion
In this paper, we have developed a suitable approach to determine constraints on the parameters of a linear hybrid automaton.Our hybrid Hoare logic combined with experimental biological traces including precise chronometrical information leads to constraints on celerities which have to be satisfied to allow the model to represent the observed behaviour (HCSP [18] is a formalism similar to ours but does not include chronometrical information along the path).The obtained constraints via the weakest precondition calculus are analysed using the solver Ibex which extracts all admissible intervals of celerities.Choosing celerity values satisfying these constraints leads to a model which exhibits simulation traces similar to the aforementioned experimental data, this approach has been tested on the simplified circadian clock model.
The soundness of our hybrid Hoare logic is proved which means that simulations obtained with parameters satisfying the computed precondition leads to simulated traces which are in concordance with the path representing the experimental data.
This work opens many outlooks.Generally, it is useful to prove the completeness of the weakest precondition calculus.Because of the continuous terms (real numbers), our hybrid Hoare logic cannot be complete regarding all possible formulas.Nevertheless, we think that our hybrid Hoare logic should be complete regarding closed propositional formulas constructed from polynomial (in)equations and logical connectors.The decidability of the theory of real closed fields implies that the precondition constraints can be analysed; it does not mean that for each semantically correct Hoare triple, there exists a proof tree built on our inference rules.Completness of our framework would mean that if a Hoare triple is semantically correct and if pre-and postconditions are expressions in the first order language of real closed fields, then there exists a proof tree for this Hoare triple.Finally, because of the combinatorial explosion of the size of the weakest precondition formula and despite some on-the-fly simplifications, it could be interesting to investigate other ways to simplify the result in some particular cases.

A Appendix: Sub-properties of the Weakest Precondition Calculus
In this appendix, we detail each subformula of the weakest precondition in Def.14.

A.1 Weakest Precondition
In order to fully compute the weakest precondition, it is required to label the fractional parts of the states mentioned in the properties.For this, we use labels called below f (final), i (initial) and m (intermediate).Moreover, by convention, we use π (resp.π) to specify the fractional part of the exit from the current discrete state (resp.entrance into the current discrete state).
Let us notice that all the following properties depend on the indices i and f used in Def.14, although for readability issues we did not mention them on the names of each sub-property.Furthermore, for a given index i, we call by convention π i u (resp.π i u ) the fractional part of the entering (resp.exiting) state inside the discrete state i.
Finally, for all variable u ∈ V and all ω ⊂ R − (v) subset of predecessors of u, we define: In other words, Φ ω v is true in a state h if and only if the resources of u are exactly ω, that is, ρ(h, v) = ω.

A.2 Discrete Transition to the Next Discrete State
For all component v ∈ V , Φ + v (∆t) (resp.Φ − v (∆t)) describes the conditions in which v increases (resp.decreases) its discrete expression level after ∆t units of time: its celerity in the current state must be positive (resp.negative) and its fractional part only depends on ∆t in the way given at the very end of Section 2.

A.3 Internal and External Walls
For all component u ∈ V , W + u (resp.W − u ) states that there is a wall preventing u to increase (resp.decrease) its qualitative state.This wall can either be an external wall EW + u (resp.EW − u ) or an internal wall IW + u (resp.IW − u ).Furthermore, Φ ω u+ (resp.Φ ω u− ), which is required in these subformulas, is true if and only if the set of resources of u is exactly ω in the state where u is increased (resp.decreased) by 1. where:

A.4 First Changing Variables
F(∆t) states that all components that are not first changing variables must either reach their border after the first changing variables, or face an internal or external wall.

A.5 Hybrid Assertions
The sub-property A(∆t, a) allows one to translate all assertion symbols given about the continuous transition related to the instruction (celerities and slides) into a property: where a is the assert part of the instruction P = (∆t, a, v±), and, for all variable u ∈ V : These sub-properties indicate that the exit position of the corresponding variable u is located on a threshold.In addition, the constraints • ∆t mean that the duration before reaching the border is lower that the one spent inside the current state (∆t).The sign of the celerity of the sliding variable u is constrained by the sub-property F and the constraint π T I M E 2 0 1 7

A.6.1 Continuous Junctions Inside Discrete States
For all component v ∈ V , and for a continuous transition between two hybrid states h = (η, π) and h = (η, π ), CJ v establishes a relationship between the fractional parts and the celerity of the variable v.If the exit fractional part of v is 0 or 1, the sign of the celerity can be deduced and the time mandatory to v to reach the border is lower than the time spent in the current discrete state.If v does not reach its border, the exact position of the entrance fractional part of v can be deduced from the exit position, the time spent in the current discrete state and the celerity.

A.6.2 Discrete Junctions Between Discrete States
For all component v ∈ V , and for a discrete transition happening on component v between an initial and a final state corresponding to the indices i and f , DJ v establishes a junction between the fractional parts of these states.This formula states that the fractional part of v switches from 1 to 0 for an increase, or from 0 to 1 for a decrease, whereas all other fractional parts are unchanged: Finally, we define: These relationships can be easily observed on Fig. 3 on the discrete transition in the centre: all fractional parts are left the same, except for the variable performing the transition.

B Soundness Proof
The soundness proof is made for each inference rule which depends on its corresponding assertion (Def.11).Each of them is treated according to the assertion type.We focus here on the proof of the soundness of the incrementation rule since the proof of the soundness of the decrementation rule is similar, and that the one for the sequential composition rule is classical.In this subsection, we consider the Hoare triple associated with the incrementation rule, described in Subsection 4.1, and a hybrid state h = (η, π) satisfying the precondition.

B.1 First Case: a =
α) Let us first prove the existence of the continuous transition.According to the sub-property Φ + v (∆t), π v = 1 (arrival at the top border of v), C v,ω,n > 0 and π v = π v + C v,ω,n × ∆t (the time spent in the current state is ∆t) if Φ ω v and η v = n are satisfied.Let us also consider the unique hybrid state h = (η, π ) such that the continuous transition h → h exists.Thus, according to Lemma 15 and Definition 11, (h, h ) (∆t, ).

5:19 β) Let us now prove the existence of the discrete transition.
Let us simplify the subformula ¬W + v ≡ ¬EW + v ∧ ¬IW + v satisfied at h .We have: which is evaluated to true because v increases its level (v+ is the discrete path atom) and thus is not already at its maximal discrete value.Thus, ¬W + v ≡ ¬IW + v : Amongst all premises of the remaining disjunctions, only one is satisfied because the current qualitative state and the next state have a unique qualitative level (η v = n and m = n + 1) and a unique set of resources (Φ ω v and Φ ω v+ ).Replacing ω and ω by the right resources of the corresponding states ρ(η, v) and ρ(η , v) and naming η the next state, we deduce: However, since Φ + v (∆t) is true at h , we have C v,ρ(η,v),ηv > 0. Thus ¬W + v is equivalent to C v,ρ(η ,v),η v ≥ 0 and the previous inequation is true since ¬W + v is satisfied at h .Thus the variable v reaches its threshold in ∆t time (Φ + v (∆t)) and crosses it (¬W + v ) allowing a discrete transition h → h which increases v because the signs of the celerities of v in h and in h are the sames.
γ) Let us finally prove that the postcondition is satisfied after the elementary path.We previously proved that there exists a unique continuous transition h → h and a discrete one h → h .Since h (D[η v \η v + 1], H i,f ), we deduce with Lemma 16: h (D[η v \η v + 1], H i,f ).The discrete transition increases the variable v (η v = η v + 1), we deduce that: The hybrid condition So, the discrete and hybrid conditions D and H f are satisfied at h and the postcondition is verified.α) Similarly to the first case, we consider the unique hybrid state h = (η, π ) such that the continuous transition h → h exists.The time spent in the current qualitative state is also ∆t (sub-property Φ + v (∆t)).Since a = , the sub-property A plays a crucial rule:

T I M E
Amongst all premises of these conjunctions, only one is satisfied because the current qualitative state has a unique qualitative level for each variable v l (η v l = n l ) and a unique set of resources for each v l (Φ ω l v l ).We can then replace slide + (u) by the sub-property S + : where ω is the resources of u and n its current qualitative level.This formula means that the exit position of the current qualitative state is on the top border (π i u = 1).We then deduce: According to Lemma 15, we have δ h (v) = δ first h = ∆t and so δ h (v) > δ h (u).In other words, u reaches its top border before v reaches its one.Thus, the continuous transition h → h is such that (h, h ) (∆t, a), see Definition 11. β and γ) The proof of the discrete transition existence from h is similar to the first case.This transition leads to h which satisfies the postcondition h (D, H) (see the stages β and γ of the first case).

B.3 Third Case
The sub-property A(∆t, a) allows one to replace the celerity C u in the assertion a by the celerity indexed by the relevant set of resources of the current qualitative state.So, we deduce C u,ρ(η,u),ηu c.From Definition 11, the unique continuous transition h → h where h = (η, π ) is such that (h, h ) (∆t, a).
β and γ) The proof of the discrete transition existence from h is similar to the first case.This transition leads to h which satisfies the postcondition h (D, H) (see the stages β and γ of the first case).
β and γ) The proof of the discrete transition existence from h is similar to the first case.This transition leads to h which satisfies the postcondition h (D, H) (see the stages β and γ of the first case).
This proof is generalisable for all logical connectives and recursively to all formulas.

Figure 1
Figure 1 Simplified model of the circadian clock in mammals.L is a zeitgeber (see Section 5).

Figure 2 -
Figure 2-Centre illustrates an example of hybrid state.The tuple of all fractional parts represents coordinates inside the current qualitative state.

3Figure 4
Figure 4 Model simulation based on arbitrarily chosen celerities satisfying the deduced constraints.The plain (resp.dashed) line represents the PER/CRY complex (resp.gene) activity.

2
Second Case: a = slide + (u) Hoare triple example: {PreC } (∆ 1 t , , g+) {PostC }.Starting from the hybrid state h1 PreC , and considering the path in bold line, it is possible to chain a continuous transition (h1 → h 1 ) of duration ∆ 1 t and a discrete transition (h 1 → h2) leading to a h2 PostC : this Hoare triple is therefore satisfied.

Table 1
Constraints obtained by computation of the weakest precondition.Left: Constraints on celerities of g and pc.Right: Constraints on celerities of L and X.