Christophe Chorro
email: christophe.chorro@univ-paris1.fr

Dominique Guegan
email: dguegan@univ-paris1.fr

Florian Ielpo
email: florianielpo@ensae.org

Martingalized Historical approach for Option Pricing

Keywords: Generalized Hyperbolic Distribution, Option pricing, Incomplete market, CAC40, Stochastic Discount Factor, Martingale Correction. JEL classification: G12, G17, C5, C8

published or not. The documents may come

Introduction

Once the time series properties of financial asset's returns have been assessed, the biggest challenge to option pricing is to thoroughly select the right pricing measure. This article questions the importance of this step. We propose to answer this question by comparing two different modeling approaches: first, the affine Stochastic Discount Factor (SDF), as presented in [START_REF] Cochrane | Asset Pricing[END_REF] and then the Empirical Martingalization Correction (EMC) presented in [START_REF] Duan | Empirical Martingale Simulation for Asset Prices[END_REF]. The rationale behind this choice is the following: for the SDF framework, when moving from the historical distribution to the pricing one, each of the first four moments are changed. This is consistent with the existence of volatility trading and with the idea of skewness premia. In the second approach, the moment modification only involves the expectation, as in the [START_REF] Black | The Pricing of Option and Corporate Liabilities[END_REF]'s model. By comparing the results on a comparable ground with each methodology, we aim at answering the question of the real impact on option prices of the modifications of higher order moments when moving from the historical to pricing probability measure. This work is quite different in its philosophy of recent works focusing on a closed form SDF for pricing under the risk neutral measure such as [START_REF] Heston | A Closed-Form GARCH Option Valuation[END_REF], [START_REF] Mercuri | Option pricing in a garch model with tempered stable innovations[END_REF] or Chorro et al. (2008a).

To achieve such a purpose, we use a time series model that properly fits the key known features of equity returns. The [START_REF] Nelson | Conditional Heteroskedasticity in Asset Returns[END_REF]'s EGARCH model is used to handle the asymmetric reaction of the volatility to negative returns. To match the returns' leptokurticity and asymmetry, the conditional distribution of the innovations is assumed to be the Generalized Hyperbolic distribution. Using a CAC 40 dataset, we show that the option mispricing errors for both the SDF and EMC methodology are fairly closed: what seems to be the key to option pricing is more the precise modeling of the conditional distribution of the returns' innovations than the change in probability measure.

The Modeling under P

Under the historical measure P, the discrete time economy we consider is characterized by the time series dynamics of the stock price process (S t ) t∈{0,1,...,T } given by

Y t = log S t S t-1 = r + m t + h t z t , S 0 = s ∈ R + , (1)
where r is the risk free rate expressed on a daily basis and supposed to be constant and where z t is a sequence of independent and identically distributed centered random variables with variance 1.

In equation (1), we consider a general conditional time varying excess of return m t that depends on the constant risk premium λ 0 . In practice, it will be fixed for the empirical study and we retain the following classical form:

m t = λ 0 h t - 1 2 h t . (2)
In order to capture asymmetry phenomena (such as leverage effects) observed in equity returns datasets, we restrict to the Nelson (1991)'s EGARCH(1, 1) model:

log(h t ) = a 0 + a 1 (| z t-1 | +γ z t-1 ) + b 1 log(h t-1 ). (3)
So as to deal with the leptokurticity and the asymmetry of the CAC 40 returns, the conditional distribution of the EGARCH model is chosen to be the Generalized Hyperbolic distribution. This distribution has been introduced by [START_REF] Barndorff-Nielsen | Exponentially Decreasing Distributions for the Logarithm of Particle Size[END_REF]. It is driven by five parameters: λ being a shape parameter, α controlling the kurtosis, β controlling the skewness, δ and µ being respectively linked to the volatility and the expectation of the process.

For (λ, α, β, δ, µ) ∈ R 5 with δ > 0 and α >| β |> 0, the one dimensional GH(λ, α, β, δ, µ) distribution is defined by the following density function

d GH (x, λ, α, β, δ, µ) = ( √ α 2 -β 2 /δ) λ √ 2πK λ (δ √ α 2 -β 2 ) ×e β(x-µ) K λ-1/2 α √ δ 2 +(x-µ) 2 √ δ 2 +(x-µ) 2 /α 1/2-λ , (4)
where K λ is the modified Bessel function of the third kind. For λ ∈ 1 2 Z, the basic properties of the Bessel function allow to find simpler forms for the density. In particular, for λ = 1, we get the Hyperbolic distribution whose log-density is a hyperbola. For λ = -1 2 , we obtain the Normal Inverse Gaussian distribution which is closed under convolution. This family has already been suggested as a model for financial price processes because its exponentially decreasing tails seem to fit some fundamental statistical behaviors of asset returns remarkably, as presented in [START_REF] Barndorff-Nielsen | Normal Inverse Gaussian Processes and the Modeling of Stock Returns[END_REF] or [START_REF] Eberlein | The Generalized Hyperbolic Model: Financial Derivatives and Risk Measures[END_REF]. For instance, in [START_REF] Chorro | Option Pricing Under GARCH Models with Generalized Hyperbolic Innovations (II): Empirical Results[END_REF] with a similar dataset, the joint use of an EGARCH model with this GH distribution offers a very realistic assessment of the returns stochastic behavior.

Option pricing modeling

Using the previous time series framework, we now discuss how to deal with the change in probability measure toward the risk neutral distribution. To serve our purpose, we consider two ways out of this problem. The first (the SDF approach) leads to a change in every moment of the distributions, whereas the second one (The EMC approach) states that the only difference between the risk neutral distribution and the historical one is the expectation.

The Stochastic Discount Factor approach

The approach that is briefly overviewed here received different presentations, the first one being reviewed in [START_REF] Cochrane | Asset Pricing[END_REF]. The second one can be found in [START_REF] Gerber | Option Pricing by Esscher Transforms[END_REF]. The methodology unfolds as follows. The stochastic discount factor is the ratio of the risk neutral distribution to the historical one. We assume for the stochastic discount factor a particular parametric form:

∀t ∈ {0, ..., T -1} M t,t+1 = e θt+1Yt+1+ξt+1 (5)
where

Y t+1 = log St+1
St and where θ t+1 and ξ t+1 are

F t = σ(z u ; 0 ≤ u ≤ t) measurable random
variables. This hypothesis is implicit in many asset pricing models such as the Black and Scholes (1973)'s model.

We compute explicitly (θ t+1 , ξ t+1 ). Considering the bond and the risky asset, the pricing relation

P t = E P [Φ T M t,T | F t ] (6)
with Φ T the financial asset future random payoff, gives for T = t + 1 the following restrictions for the stochastic discount factor

     E P [e r M t,t+1 | F t ] = 1 E P [e Yt+1 M t,t+1 | F t ] = 1. (7)
Using the previous restrictions with Generalized Hyperbolic innovations, the distribution of Y t given

F t-1 under Q is again a Generalized Hyperbolic one GH λ, α √ h t , β √ h t + θ q t , δ h t , M t + µ h t (8)
with M t = r + m t and where (θ q t+1 , ξ q t+1 ) is the unique solution of (7). This essential property is proved in Chorro et al. (2008a) where the corresponding Monte Carlo pricing methodology is also explained.

Empirical Martingale Correction approach

This second approach to "risk neutralization" is inspired from the variance reduction tool presented in [START_REF] Duan | Empirical Martingale Simulation for Asset Prices[END_REF]. A similar idea used together with an optimization procedure is followed by [START_REF] Barone-Adesi | GARCH Option Pricing Model in Incomplete Markets[END_REF] to overcome the classical option pricing framework. Here, we make no assumption on the shape of the pricing kernel and we compute prices for options with time to maturity e Y k where Y k is given by (1) .

(T -t)
The Empirical Martingale Correction works by replacing the previous sampled prices by:

ST,i = S T,i 1/N N i=1 S T,i S t e -r(T -t) . (9)
The sampled average of ST,i is exactly equal to S t e r(T -t) , that is the risk neutral conditional expectation. With this approach, we only shift the historical distribution in a way that prevents arbitrage opportunities by implicitly changing the drift of this distribution.

Results

The dataset that we use contains the following time series. We consider daily log returns of the French CAC 40 whose price at time t is denoted S t . The sample starts on January 2, 1988 and ends on October 26, 2007.

Concerning the estimation we are naturally led to a two steps methodology. We first estimate the EGARCH parameters driving the historical distribution with a Quasi Maximum Likelihood Estimation using a rolling window of 4'000 observations and ending at time t. At the second stage, since we exactly know the form of the density function of a GH distribution (7) we adopt a classical maximum likelihood approach to estimate the unknown remaining parameters (λ, α, β, δ, µ) using the residuals obtained at the previous stage. The key point in our approach is to maintain as much outliers in the dataset as possible: these extremal events are essential to fit the GH parameters and to control the tail behaviors.

The average values obtained for each parameters driving the EGARCH model are presented in table 1 and those driving the GH distribution are presented in 2. Then, we price options using either the SDF or the EMC methodologies.

We get the approximated option price C(.) as the sample average of the simulated final prices:

C(t, T, K) = 1 N N i=1 (S T,i -K) + e -r(T -t) . (10)
In practice, the number of sampled paths N is equal to 10 000.

Finally, to compare the two option pricing methods, we use the average absolute relative pricing errors criterion for the working days t between January 2, 2006 (τ 1 ) and October 26, 2007 (τ 2 ). Let C(t, T j , K i ) be the estimated call option price with a time to maturity equal to T j -t and a strike price worth K i and C(t, T j , K i ) be the corresponding quoted market option price. Then the criterion is

AARP E = τ2 t=τ1 Jt j=1 Gt,j i=1 C(t, T j , K i ) -C(t, T j , K i ) C(t, T j , K i ) , (11)
where J t is the number of call option maturities T j available at time t and G t,j the number of strikes K i available at time t for this particular maturity T j .

This criterion is robust to the well known fact that option pricing errors are proportional to the moneyness: out of the money call option prices are very low, and so are the usual errors found. The converse is true for deep in the money option prices. Table 3 presents the value obtained for this criterion with the two different approaches to produce this change in distribution.

First, the scale of errors obtained are very closed to the one obtained in [START_REF] Barone-Adesi | GARCH Option Pricing Model in Incomplete Markets[END_REF], whereas they directly optimize the GARCH parameters to match the option prices. Second, for time to maturities beyond 1 year, the pricing kernel based option prices are clearly dominated by the simpler method. Third, for time to maturities between .5 and 1 year, the results get trickier to analyze: depending on the moneyness, one method may dominate the other. Globally, the errors obtained for each approach are very close. The main conclusion of these estimation results can be summarized as follows: when taking into account the returns' time series properties with much care, the change in distribution from the historical to the risk neutral probabilities is not that essential anymore. The simple martingalization produced by the EMC approach is enough to ensure small pricing errors.

Conclusion

We compare two different approaches to price options in a discrete time setting. The first one lies on classical theoretical results derived from the [START_REF] Black | The Pricing of Option and Corporate Liabilities[END_REF]'s hedging methodology to change from the historical to the risk neutral distribution. The second one is only based on no arbitrage requirements avoiding any technical arguments. By doing this, we end up with the intuition that a simple change in the drift clearly challenges a more sophisticated methodology that is widely used in the financial economics literature. The use of a time series model with innovations that encompasses most of the known features of financial assets' returns seems to be more essential to option pricing.

A forthcoming paper is ongoing to deeply explore this feature with other datasets and also to compare the (EMC) to the pricing methodology proposed by [START_REF] Elliott | A Discrete Time Equivalent Martingale Measure[END_REF].

  by simulating sampled paths of the previous stochastic model under the historical distribution

P.

To rule out arbitrage opportunities, we directly impose risk neutral constraints. The i th sampled historical final price for the underlying is denoted by S T,i = S t T k=t+1

Table 1 :

 1 Average -5.51 0.47 0.01 0.91 0.08 Std. Dev. 0.33 0.04 0.02 0.02 0.03 Descriptive statistics on the average estimated parameters for the EGARCH model.

		a 0	b 1	a 1	γ	λ 0
		α	β	δ	µ	λ
	Average	308.42 -2.77 335.64 3.04 -78.46
	Std. Dev. 139.90 5.08 152.02 5.51 76.68

Table 2 :

 2 Descriptive statistics on the average estimated parameters for the GH distribution.

				.5 < Maturity < 1
			[.8-.9]	[.9-1]	[1-1.1]	[1.1-1.2]	> 1.2
	SDF correction	Nelson GH	0.072508156	0.034220274	0.08973236	0.211844412	0.473990647
	Mart. Correction	Nelson GH	0.070375688	0.032051441	0.112397312	0.225184538	0.40081137
				Maturity > 1
			[.8-.9]	[.9-1]	[1-1.1]	[1.1-1.2]	> 1.2
	SDF correction	Nelson GH	0.120757162	0.088165393	0.069513071	0.118837264	0.261635523
	Mart. correction	Nelson GH	0.119581221	0.087505694	0.067820071	0.119243876	0.266680359

Table 3 :

 3 Absolute average mispricing errors for the CAC 40 french index disaggregated by maturities and moneyness, using the stochastic discount factor restrictions (SDF correction) or the empirical martingale correction (Mart correction).


Keywords: Réseaux génétiques, Commande en feedback, Linéaire par morceaux, Solutions périodiques

Commande qualitative de solutions périodiques de systèmes affines par morceaux ; application aux réseaux génétiques Résumé : Les systèmes hybrides, en particulier affines par morceaux (APM), sont souvent employés comme modèles de réseaux génétiques. Dans ce rapport nous approfondissons des travaux antérieurs sur la commande de tels systèmes, utilisant également des résultats récents garantissant l'existence et l'unicité de cycles limites, sur la seule base d'une abstraction discrète du système et de sa structure d'interaction. L'objectif est de contrôler le graphe de transitions d'états du système APM pour obtenir un comportement périodique, ce qui est une propriété très importante de nombreux systèmes biologiques. Nous montrons comment commander l'apparition ou la suppression d'un cycle limite unique, par une action qualitative sur les taux de dégradation d'un système APM, aussi bien par commande statique que dynamique, c'est-à-dire par le couplage adéquat d'un sous-réseau contrôleur. Ceci est illustré sur deux réseaux de gènes classiques, présentant une structure de boucles de retro-action imbriquées.

Introduction

Gene regulatory networks often display both robustness and steep, almost switchlike, response to transcriptional control. This motivates the use of an approximation of these response laws by piecewise affine differential (PWA) equations, to build hybrid models of genetic networks. PWA systems are affine in each rectangular domain (or box) of the state space. They have been introduced in the 1970's by Leon Glass [START_REF] Glass | Combinatorial and topological methods in nonlinear chemical kinetics[END_REF] to model genetic networks. It has led to a long series of works by different authors, dealing with various aspects of these equations, e.g. [START_REF] Casey | Piecewise-linear models of genetic regulatory networks: Equilibria and their stability[END_REF][START_REF] Edwards | Symbolic dynamics and computation in model gene networks[END_REF][START_REF] Farcot | Geometric properties of piecewise affine biological network models[END_REF][START_REF] Glass | Combinatorial and topological methods in nonlinear chemical kinetics[END_REF][START_REF] Gouzé | A class of piecewise linear differential equations arising in biological models[END_REF]. They have been used also as models of concrete biological systems [START_REF] Jong | Qualitative simulation of the initiation of sporulation in bacillus subtilis[END_REF].

From an hybrid system point of view, the behaviour of PWA systems can be described by a transition graph, which is an abstraction (in the hybrid system sense) of the continuous system. This transition graph describes the possible transition between the boxes. It is also possible to check properties of the transition graph by model checking techniques [START_REF] Batt | Model checking genetic regulatory networks with parameter uncertainty[END_REF].

Nowadays, the extraordinary development of biomolecular experimental techniques makes it possible to design and implement control laws in the cell system. The authors have recently developed a mathematical framework for controlling gene networks with hybrid controls [START_REF] Farcot | A mathematical framework for the control of piecewise-affine models of gene networks[END_REF]; these controls are defined on each box. It is easy to see that this amounts to change the transition graph to obtain the desired one.

From another point of view, more oriented towards dynamical systems, it is also possible to obtain results concerning the limit cycles in PWA systems (see [START_REF] Glass | Stable oscillations in mathematical models of biological control systems[END_REF] and the recent generalisation in [START_REF] Farcot | Periodic solutions of piecewise affine gene network models: the case of a negative feedback[END_REF]). For example, one can show that a simple negative loop in dimension greater that two produces a unique stable limit cycle. It is clear that biological oscillations play a fundamental role in the cell ( [START_REF] Elowitz | A synthetic oscillatory network of transcriptional regulators[END_REF]).

Our aim in this paper is to control PWA systems to make a single stable limit cycle appear or disappear. To fulfil that goal, after some recalls concerning the PWA systems, we use some results enabling to deduce the existence of a single stable limit cycle in the state space from a periodic behaviour in a box sequence (section 3), then the results on the control of the transition graph in the space of boxes (section 4), to obtain our main results illustrated by 2 examples (section 5).

Related works on control aspects concern the affine or multi-affine hybrid systems ( [START_REF] Habets | A control problem for affine dynamical systems on a full-dimensional polytope[END_REF][START_REF] Belta | Controlling a class of nonlinear systems on rectangles[END_REF]). The authors derive sufficient conditions for driving all the solutions out of some box. Other related works study the existence of limit cycles in the state space [START_REF] Glass | Stable oscillations in mathematical models of biological control systems[END_REF][START_REF] Lu | Structural principles for periodic orbits in glass networks[END_REF]. We are not aware of works linking control theory and limit cycle for this class of hybrid systems.

Piecewise affine models 2.1 General formulation

This section contains basic definitions and notations for piecewise affine models [START_REF] Glass | Combinatorial and topological methods in nonlinear chemical kinetics[END_REF][START_REF] Edwards | Analysis of continuous-time switching networks[END_REF][START_REF] Farcot | Geometric properties of piecewise affine biological network models[END_REF][START_REF] Jong | Qualitative simulation of genetic regulatory networks using piecewise-linear models[END_REF]. The general form of these models can be written as:

dx dt = κ(x) -Γ(x)x (1)
The variables (x 1 . . . x n ) represent levels of expression of n interacting genes, meaning in general concentrations of the mRNA or protein they code for. We will simply call genes the n network elements in the following. Since gene transcriptional regulation is often considered to follow a steep sigmoid law, an approximation by a step function has been proposed to model the response of a gene (i.e. its rate of transcription) to the activity of its regulators [START_REF] Glass | Combinatorial and topological methods in nonlinear chemical kinetics[END_REF]. We use the notation:

s + (x, θ) = 0 if x < θ, s + (x, θ) = 1 if x > θ,
This describes an effect of activation, whereas s -(x, θ) = 1s + (x, θ) represents inhibition. Unless further precision are given, we leave this function undefined at its threshold value θ.

The maps κ : R n + → R n + and Γ : R n + → R n×n + in (1) are usually multivariate polynomials (in general multi-affine), applied to step functions of the form s ± (x i , θ i ), where for each i ∈ {1, • • • , n} the threshold values belong to a finite set

Θ i = {θ 0 i , . . . , θ qi i }. (2)
We suppose that the thresholds are ordered (i.e. θ j i < θ j+1 i ), and the extreme values θ 0 i = 0 and θ qi i represent the range of values taken by x i rather than thresholds. Γ is a diagonal matrix whose diagonal entries Γ ii = γ i , are degradation rates of variables in the system. Obviously, Γ and the production rate κ are piecewiseconstant, taking fixed values in the rectangular domains obtained as Cartesian products of intervals bounded by values in the threshold sets [START_REF] Batt | Model checking genetic regulatory networks with parameter uncertainty[END_REF]. These rectangles, or boxes, or regular domains [START_REF] Plahte | A methodological basis for description and analysis of systems with complex switch-like interactions[END_REF][START_REF] Jong | Qualitative simulation of genetic regulatory networks using piecewise-linear models[END_REF], are well characterised by integer vectors: we will often refer to a box D a = i (θ ai-1 i , θ ai i ) by its lower-corner index a = (a 1 -1 . . . a n -1). The set of boxes is then isomorphic to

A = n i=1 {0, • • • , q i -1}, (3)
Also, the following pairs of functions will be convenient notations:

θ ± i : A → Θ i , θ - i (a) = θ ai-1 i and θ + i (a) = θ ai i .
Let us call singular domains the intersections of closure of boxes with threshold hyperplanes, where some x i ∈ Θ i \ {θ 0 i , θ qi i }. On these domains, the right-hand side of (1) is undefined in general. Although the notion of Filippov solution provides a generic solution to this problem [START_REF] Gouzé | A class of piecewise linear differential equations arising in biological models[END_REF], in the case where the normal of the vector field has the same sign on both side of these singular hyperplanes, it is more simply possible to extend the flow by continuity. In the remaining of this paper, we will only consider trajectories which do not meet any singular domain, a fact holding necessarily in absence of auto-regulation, i.e. when no κ i depends on x i . This leads to the following hypothesis: ∀i ∈ {1, • • • , n}, κ i and γ i do not depend on x i .

(H1)

On any regular domain of index a ∈ A, the rates κ = κ(a) and Γ = Γ(a) are constant, and thus equation ( 1) is affine. Its solution is explicitly known, for each coordinate i :

ϕ i (x, t) = x i (t) = φ i (a) + e -γit (x i -φ i (a)) , (4)

INRIA

where t ∈ R + is such that x(t) ∈ D a , and

φ(a) = (φ 1 (a) • • • φ n (a)) = κ 1 (a) γ 1 (a) • • • κ n (a) γ n (a) .
It is clearly an attractive equilibrium of the flow (4). It will be called focal point in the following for reasons we explain now. Let us first make the generic assumption that no focal point lies on a singular domain:

∀a ∈ A, φ(a) ∈ a ′ ∈A D a ′ . (H2)
Then, if φ(a) ∈ D a , it is an asymptotically stable steady state of system (1). Otherwise, the flow will reach the (topological) boundary ∂D a in finite time. At this point, the value of κ (and thus, of φ) changes, and the flow changes its direction, evolving towards a new focal point. The same process carries on repeatedly. It follows that the continuous trajectories are entirely characterised by their successive intersections with the boundaries of regular domains (extending them by continuity, as mentioned previously). This sequence depends essentially on the position of focal points with respect to thresholds. Actually, {x |

x i = θ - i (a)} (resp. {x | x i = θ + i (a)}) can be crossed if and only if φ i (a) < θ - i (a) (resp. φ i (a) > θ + i (a)). Then, let us denote I + out (a) = {i ∈ {1, • • • , n}| φ i > θ + i (a)}, and similarly I - out (a) = {i ∈ {1, • • • , n}| φ i < θ - i (a)}. Then, I out (a) = I + out (a) ∪ I - out (a)
is the set of escaping directions of D a . Also, we call walls the intersections of threshold hyperplanes with the boundary of a regular domain. When it is unambiguous, we will omit the dependence on a in the sequel. Now, in each direction i ∈ I out the time at which x(t) encounters the corresponding hyperplane, for x ∈ D a , is readily calculated:

τ i (x) = -1 γ i ln φ i -θ ± i φ i -x i , i ∈ I ± out . (5)
Then, τ (x) = min i∈Iout τ i (x), is the exit time of D a for the trajectory with initial condition x. Then, we define a transition map T a : ∂D a → ∂D a :

T a x = ϕ (x, τ (x)) = φ + α(x)(x -φ). (6)
where α(x) = exp(-τ (x)Γ). The map above is defined locally, at a domain D a . However, under our assumption (H1), any wall can be considered as escaping in one of the two regular domains it bounds, and incoming in the other. Hence, on any point of the interior of a wall, there is no ambiguity on which a to choose in expression [START_REF] Jong | Qualitative simulation of genetic regulatory networks using piecewise-linear models[END_REF], and there is a well defined global transition map on the union of walls, denoted T . On the boundaries of walls, at intersections between several threshold hyperplanes, the concept of Filippov solution would be required in general [START_REF] Gouzé | A class of piecewise linear differential equations arising in biological models[END_REF]. This problem will either be solved on a case by case basis, or we implicitly restrict our attention to the (full Lebesgue measure) set of trajectories which never intersect more than one threshold hyperplane.

To conclude this section let us define the state transition graph TG associated to a system of the form (1) as the pair (A, E) of nodes and oriented edges, where A is defined in (3) and (a, b) ∈ E ⊂ A 2 if and only if ∂D a ∩ ∂D b = ∅, and there exists a positive Lebesgue measure set of trajectories going from D a to D b . It is not difficult to see that this is equivalent to b being of the form a ± e i , with i ∈ I ± out (a) and e i a standard basis vector. From now on, it will always be assumed that (H1) and (H2) hold, at least in some region of state space (or transition graph) on which we focus.

Illustrative example

Let us now illustrate the previous notions on a well-known example with two variables, in order to help the reader's intuition. Consider two genes repressing each other's transcription. In the context of piecewise-affine models, this would be described by the system below:

1 2 ẋ1 = κ 0 1 + κ 1 1 s -(x 2 , θ 1 2 ) -γ 1 x 1 ẋ2 = κ 0 2 + κ 1 2 s -(x 1 , θ 1 1 ) -γ 2 x 2
where inhibition is modeled by s -(x, θ), as already mentioned. A usual notation for the interaction graph uses to denote inhibition, and to denote activation. The two constants κ 0 i represent the lowest level of production rates of the two species in interaction. It will be zero in general, but may also be a very low positive constant, in some cases where a gene needs to be expressed permanently. In the given equation, arbitrary parameters may lead to spurious behaviour, in particular an inhibition which would not drop its target variable below its threshold. To avoid this, it suffices to assume the following conditions on focal points' coordinates:

κ 0 i γ i < θ 1 i and κ 0 i + κ 1 i γ i > θ 1 i , for i = 1, 2.
This might be called structural constraints on parameters. The phase space of this system is schematised on Figure 1.

Then, the transition graph of the system takes the form:

01 11 00 10
where circled states are those with no successor. It appears in this case that TG constitutes a reliable abstraction of the system's behaviour. In general, things are not as convenient, and some paths in the transition graph may be spurious.

In particular, cyclic paths may correspond to damped oscillations of the original system, but even this cannot be always ascertained without a precise knowledge of the parameter, see section 3 for related results. However, one general goal of

INRIA θ 1 1 θ 2 1 θ 1 2 θ 2 2 κ 0 1 γ1 κ 0 1 +κ 1 1 γ1 κ 0 2 γ2 κ 0 2 +κ 1 2 γ2
Figure 1: The dashed lines represent threshold hyperplanes, and define a rectangular partition of state space, and dotted lines indicate focal points' coordinates.

Arrows represent schematic flow lines, pointed toward these limit points. Note that pieces of trajectories are depicted as straight lines, which is the case when all degradation rates γ i coincide, a fact we never assume in the present study.

the present study will be to search for feedback control laws ensuring that given systems are indeed well characterised by their abstraction. Such a property can be deduced from the shape of the abstraction TG itself, whence the term 'qualitative control'.

Stability and limit cycles

Periodic solutions have soon been a prominent topic of study for systems of the form (1) [START_REF] Glass | Stable oscillations in mathematical models of biological control systems[END_REF][START_REF] Snoussi | Qualitative dynamics of piecewise-linear differential equations: a discrete mapping approach[END_REF][START_REF] Mestl | Periodic solutions of piecewiselinear differential equations[END_REF][START_REF] Edwards | Analysis of continuous-time switching networks[END_REF][START_REF] Lu | Structural principles for periodic orbits in glass networks[END_REF]. With the notable exception of [START_REF] Snoussi | Qualitative dynamics of piecewise-linear differential equations: a discrete mapping approach[END_REF], all these studies focused on the special case where Γ is a scalar matrix, which greatly simplifies the analysis, since trajectories in each box are then straight lines towards the focal point, as in Figure 1. In a recent work [START_REF] Farcot | Periodic solutions of piecewise affine gene network models: the case of a negative feedback[END_REF][START_REF] Farcot | Periodic solutions of piecewise affine gene network models: the case of a negative feedback loop[END_REF][START_REF] Farcot | Limit cycles in piecewise-affine gene network models with multiple interaction loops[END_REF], we have shown that the local monotonicity properties of transition maps can be used to prove existence and uniqueness of limit cycles in systems like [START_REF] Andrianantoandro | Synthetic biology: new engineering rules for an emerging discipline[END_REF]. In this section we recall without proof some of these results.

In the rest of this section we consider a piecewise-affine system such that there exists a sequence C = {a 0 . . . a ℓ-1 } of regular domains which is a cycle in the transition graph, and study periodic solutions in this sequence. We abbreviate the focal points of these boxes as φ i = φ(a i ). Let us now define a property of these focal points: we say that the points

φ i are aligned if ∀i ∈ {0, • • • , ℓ -1}, ∃!j ∈ {1, • • • , n}, φ i+1 j -φ i j = 0, (7)
where φ ℓ and φ 0 are identified. Since C is supposed to be a cycle in TG, for each pair (a i , a i+1 ) of successive boxes there must be at least one coordinate at which their focal points differ, namely the only s i ∈ I out (a i ) such that a i+1 = a i ± e si . We keep on denoting s i this switching coordinate in the following. Hence condition [START_REF] Deveaux | The ethanol switch : a tool for tissue specific gene induction during plant development[END_REF] means that s i is the only coordinate in which φ i and φ i+1 differ. This implies in particular that a i+1 is the only successor of a i , i.e. there is no edge in TG from C to A \ C. It might seem intuitive in this case that all orbits in C converge either to a unique limit cycle, or to a point at the intersection of all crossed thresholds. However, this fact has only been proved for uniform decay rates (i.e. Γ scalar), [START_REF] Glass | Stable oscillations in mathematical models of biological control systems[END_REF], and its validity with distinct decay remains an open question. The condition [START_REF] Deveaux | The ethanol switch : a tool for tissue specific gene induction during plant development[END_REF] is of purely geometric nature. However, it can be shown that it holds necessarily when the interaction graph has degree one or less, see [START_REF] Farcot | Limit cycles in piecewise-affine gene network models with multiple interaction loops[END_REF] for more details.

If

{s i } 0 i<ℓ = {1, • • • , n}, i.
e. all variables are switching along C, then the intersection of all walls between boxes in C is either a single point, which we denote θ C , or it is empty. The latter holds when two distinct thresholds are crossed in at least one direction. When defined, θ C is a fixed point for any continuous extension of the flow in C, see [START_REF] Farcot | Limit cycles in piecewise-affine gene network models with multiple interaction loops[END_REF].

Let us now rephrase the main result from [START_REF] Farcot | Limit cycles in piecewise-affine gene network models with multiple interaction loops[END_REF].

Theorem 1. Let C = {a 0 , a 1 • • • a ℓ-1
} denote a sequence of regular domains which is periodically visited by the flow, and whose focal points satisfy condition [START_REF] Deveaux | The ethanol switch : a tool for tissue specific gene induction during plant development[END_REF]. Suppose also that all variables are switching at least once. Let W denote the wall ∂D a 0 ∩ ∂D a 1 , and consider the first return map T : W → W defined as the composite of local transition maps along C. A) If a single threshold is crossed in each direction, let λ = ρ(DT(θ C )), the spectral radius of the differential DT(θ C ). Then, the following alternative holds:

i) if λ 1, then ∀x ∈ W , T n x → θ C when n → ∞.
ii) if λ > 1 then there exists a unique fixed point different from θ C , say q = Tq. Moreover, for every x ∈ W \ {θ C }, T n x → q as n → ∞. B) If there are two distinct crossed thresholds in at least one direction, then the conclusion of ii) holds.

In [START_REF] Farcot | Periodic solutions of piecewise affine gene network models: the case of a negative feedback[END_REF][START_REF] Farcot | Periodic solutions of piecewise affine gene network models: the case of a negative feedback loop[END_REF] we have resolved the alternative above for a particular class of systems:

Theorem 2. Consider a negative feedback loop system of the form

ẋi = κ 0 i + s εi (x i-1 , θ i-1 ) -γ i x i , ε i ∈ {-, +} i ∈ {1, • • • , n},
with subscripts understood modulo n, and an odd number of negative ε i . It can be shown that there exists a cycle C in TG whose focal points satisfy [START_REF] Deveaux | The ethanol switch : a tool for tissue specific gene induction during plant development[END_REF]. Then, in Theorem 1, A.i) holds in dimension n = 2, and A.ii) holds for all n 3.

Piecewise Control

Feedback regulation is naturally present in many biological systems, as the widespread appelation 'regulatory network' suggests. Hence, it seems appropriate to take advantage of the important body of work developed in feedback control theory for decades, in order to study gene regulatory networks and related systems [START_REF]Control Theory and Systems Biology[END_REF][START_REF] Sontag | Molecular systems biology and control[END_REF].

In particular, the recent advent of so called synthetic biology [START_REF] Andrianantoandro | Synthetic biology: new engineering rules for an emerging discipline[END_REF][START_REF] Kobayashi | Programmable cells: interfacing natural and engineered gene networks[END_REF] , has led to a situation where gene regulatory processes are not only studied, but designed to perform certain functions. Hence, autonomous systems of the form (1) should to be extended, so as to include possible input variables. In [START_REF] Farcot | A mathematical framework for the control of piecewise-affine models of gene networks[END_REF],

INRIA
we have presented such an extension, where both production and decay terms have some additional argument u ∈ R p , of which they were affine functions. In this context, we defined a class of qualitative control problems, and showed that were equivalent to some linear programming problems. As in our previous work, we consider qualitative feedback laws, in the sense that they depend only on the box containing the state vector, rather than its exact value. This choice is motivated by robustness purposes. More pragmatically, it is also due to the fact that recent techniques allow for the observation of qualitative characteristics of biological systems, for instance by live imaging, using confocal microscopy, of GFP marker lines, where the measured state is closer to an ON/OFF signal than to a real number.

Recent experimental techniques allow furthermore for the reversible induction of specific genes at a chosen instant, for instance using promoters inducible by ethanol [START_REF] Deveaux | The ethanol switch : a tool for tissue specific gene induction during plant development[END_REF], or light [START_REF] Shimizu-Sato | A light-switchable gene promoter system[END_REF], to name only two. Also, degradation rates may be modified, either directly by introducing a drug [START_REF] Wyke | Induction of protein degradation in skeletal muscle by a phorbol ester involves upregulation of the ubiquitin-proteasome proteolytic pathway[END_REF], or via a designed genetic circuit [START_REF] Grilly | A synthetic gene network for tuning protein degradation in saccharomyces cerevisiae[END_REF], which might be induced using previously mentioned techniques.

To simplify the presentation, we focus in this paper on the particular case where decay rates can be linearly controlled by a scalar and bounded input u. For each i ∈ {1, • • • , n}, let us denote this as:

dx i dt = κ i (x) -(γ 1 i (x)u + γ 0 i (x))x, u ∈ [0, U ] ⊂ R + , (8)
where γ 0 i and γ 1 i are piecewise constant functions assumed to satisfy γ 0 i > 0 and γ 1 i > -

γ 0 i
U , in any box. This ensures that decay rates are positive, but yet can be decreasing functions of u (for γ 1 i < 0). Now, a feedback law depending only on the qualitative state of the system is simply a expressed as the composite of a map a D a → A indicating the box of the current state, with a function u : A → [0, U ] which represents the control law itself. In other words, in each box a constant input value is chosen. For a fixed law of this form, it is clear that the dynamics of ( 8) is entirely determined, and in particular we denote its transition graph by TG(u).

Let us now recall our definition of control problem.

Global Control Problem: Let TG ⋆ = (A, E ⋆ ) be a transition graph. Find a feedback law u : A → [0, U ] such that TG(u) = TG ⋆ .

Clearly, E ⋆ cannot be arbitrary in A 2 , and must in particular contain only arrows of the form (a, a ± e i ). Now in the present, restricted, context the equivalent linear programming problem described in [START_REF] Farcot | A mathematical framework for the control of piecewise-affine models of gene networks[END_REF] is very simple. For each a ∈ A, the control problem above requires that the focal point φ(a, u(a)) belongs to a certain union of boxes, i.e. its coordinates must satisfy inequalities of the form θ

j -(a) i < κ i (a)/(γ 1 i (a)u(a) + γ 0 i (a)) < θ j + (a) i
, or equivalently

κ i (a) -γ 0 i (a)θ j + (a) i γ 1 i (a)θ j + (a) i < u(a) < κ i (a) -γ 0 i (a)θ j -(a) i γ 1 i (a)θ j -(a) i (9)
if γ 1 i (a) > 0, and in reverse order otherwise. Hence, the solution set of the control problem is just the Cartesian product of all intervals of the form (9), when a varies in A. It is thus identical to a rectangle in R #A (where # denotes cardinality), which is of full dimension if and only if the problem admits a solution. Thanks to the explicit description [START_REF] Edwards | Symbolic dynamics and computation in model gene networks[END_REF], this set can be computed with a complexity which is linear in #A. The latter grows exponentially with the dimension of the system, but in practice, one will face problems where E and E ⋆ only differ on a subset of initial vertices, say A ⋆ , and then the actual complexity will be of order #A ⋆ .

In addition to this type of control, we introduce in this note some first hints toward dynamic feedback control, where instead of a direct feedback u, one uses some additional variable (here a single one), evolving in time according to a system of the form (1), and coupled to the initial system. This is suggested by the rectangular form of admissible inputs found in (9): instead of fixing an arbitrary value in a rectangle of an external input space, one increases the state space dimension, which has the effect of adding new boxes to the system. The dynamics of the supplementary variables is then defined by analogy with the direct feedback case: when the initial variables are in a box a, this makes additional variables tend to a box of the form ( 9). This raises a number of questions, in large part due to the fact that instead of applying an input u(a) instantaneously when entering box D a , the feedback now tends toward some value, which takes some time. Instead of fully developing a general theory, we thus have to chosen to illustrate it on a simple example, in section 5.1.

Examples

We now illustrate with examples how it is possible to combine results of the two previous sections, and compute qualitative feedback laws ensuring (or precluding) the existence and uniqueness of oscillatory behaviour of a system of the form (8).

Example 1: disappearance of a limit cycle

Consider the following two dimensional system:

ẋ1 (t) = K 1 s -(x 2 ) -(γ 1 1 u + γ 0 1 )x 1 ẋ2 (t) = K 2 [s + (x 1 , θ 1 1 )s + (x 2 ) + s + (x 1 , θ 2 1 )s -(x 2 )] -γ 0 2 x 2 (10)
where x 2 has a unique threshold, and

s ± (x 2 ) = s ± (x 2 , θ 1
2 ). We assume moreover that the following inequalities stand:

γ 1 1 > 0, K 1 > γ 0 1 θ 2 1 , K 2 > γ 0 2 θ 1 2 , (11)
so that the first decay rate increases with u. Also, the interactions are functional: an activation of a variable leads to the corresponding focal point coordinate being above a variable's threshold (chosen as the highest one for x 1 , since otherwise θ 2 1 cannot be crossed from below). Remark that in this system, x 2 violates (H1). However, it will appear soon that this autoregulation is only effective at a single wall, which is unstable, and thus can be ignored safely. This system corresponds to a negative feedback loop, where x 2 is moreover able INRIA to modulate its activation by x 1 : when x 2 is above its threshold, the interaction is more efficient, since it is active at a lower threshold θ 1 1 < θ 2 1 . Biologically, this may happen if the proteins coded by x 1 and x 2 form a dimer, which activates x 2 more efficiently than x 1 protein alone. This is reminiscent of the mixed feedback loop, a very widespread module able to display various behaviours [START_REF] François | Core genetic module: The mixed feedback loop[END_REF]. It might be depicted by this graph 1 2

As seen in the equations, the scalar input is assumed to affect the first decay rate, but not the second (i.e. γ 1 2 = 0). Now, one readily computes the focal points of all boxes: 00 01 10 11 20 21

φ 1 0 φ 1 0 φ 1 0 0 0 0 φ 2 φ 2 φ 2 (12)
where φ 1 is an abbreviation for K 1 /(γ 1 1 u + γ 0 1 ), and φ 2 for K 2 /γ 0 2 . Under the constraints [START_REF] Farcot | Geometric properties of piecewise affine biological network models[END_REF], this readily leads to the transition graph in absence of input (i.e. u = 0 in all boxes):

TG(0) = 01 11 21 00 10 20
The dotted line represents an unstable wall, for which Filippov theory would be required for full rigour. However, this wall is not reachable, and we ignore it afterwards. Now, since this graph has a cycle, the two thresholds θ 1,2 1 are crossed, and ( 12) is easily seen to imply condition [START_REF] Deveaux | The ethanol switch : a tool for tissue specific gene induction during plant development[END_REF], conclusion B) of Theorem 1 applies : there is a unique stable limit cycle. Now, in accordance with the section's title, let us look for a u leading to:

TG ⋆ = 01 11 21 00 10 20
Clearly from TG ⋆ , the box D 10 attracts trajectories from all other boxes, and contains its own focal point, which is thus a globally asymptotically stable equilibrium. The only states whose successors differ in TG(0), and TG ⋆ are 10 and 20, hence we assume u(a) = 0 for all other a ∈ A, or A ⋆ = {10, 20} to recall the notations of previous section. Then, Eq. ( 9) with θ

j -(a) 1 = θ 1 1 and θ j + (a) 1 = θ 2 1 gives: K 1 -γ 0 1 θ 2 1 γ 1 1 θ 2 1 < u(a) < K 1 -γ 0 1 θ 1 1 γ 1 1 θ 1 1 (13)
for both a ∈ A ⋆ . This defines a nonempty interval by θ 1 1 < θ 2 1 , hence the Control Problem of previous section can be solve under constraints [START_REF] Farcot | Geometric properties of piecewise affine biological network models[END_REF]. An illustration on a numerical example is shown Figure 2. Now, let us focus on the question of realising an extended network which solves the same problem, by adding a variable to system [START_REF] Elowitz | A synthetic oscillatory network of transcriptional regulators[END_REF]. In other words, one now seeks to impose the dynamics described by TG ⋆ using dynamic feedback. Biologically, this amounts to designing a genetic construct whose promoter depends transcriptionaly on x 1 and x 2 , and increases the degradation rate of x 1 . Let us denote by y the expression level of this additional gene. The most obvious version of such an extended system arises by increasing y production rate exactly at boxes in A ⋆ :

       ẋ1 (t) = K 1 s -(x 2 ) -(γ 1 1 υ s + (y) + γ 0 1 )x 1 ẋ2 (t) = K 2 [s + (x 1 , θ 1 1 )s + (x 2 ) + s + (x 1 , θ 2 1 )s -(x 2 )] -γ 0 2 x 2 ẏ(t) = s + (x 1 , θ 1 1 )s -(x 2 ) -γ y y (14)
υ a constant in the interval [START_REF] Farcot | A mathematical framework for the control of piecewise-affine models of gene networks[END_REF], so that forcing s + (y) = 1 would lead us back to a static feedback solution. This use of a single constant υ is possible in this particular example because constraints (13) are identical for the two boxes in A ⋆ , but it should be noted that in general several constants might be required. We consider without loss of generality that y ∈ [0, 1/γ y ], since higher values of y tend to 1/γ y or 0. Also, s + (y) is defined with respect to a threshold θ y ∈ (0, 1). We also assume θ y γ y < 1, ensuring that y may cross its threshold when activated. Now, ( 14) defines an autonomous systems of the form (1), whose transition graph has indeed a fixed point 101:

INRIA

This fixed point corresponds the fixed point 10 of TG ⋆ : in fact, the upper part of the graph above, where s + (y) = 1 is exactly TG ⋆ . However, it is not invariant, and some trajectories can escape to s + (y) = 0, where we see TG(0), and thus the possibility of periodic solutions. Besides, there are other cycles in this graph. Unlike static feedback control -and more realistically -the effect of y on γ 1 takes some positive time, explaining why the situation is not a direct translation of previous case. We will now show that under additional constraints of the parameters governing y's dynamics, it is possible to guarantee that D 101 contains a globally asymptotically stable equilibrium. To achieve this, let us rephrase a lemma, proved as Lemma 1 in [START_REF] Farcot | Geometric properties of piecewise affine biological network models[END_REF]:

Lemma 1. For any box, there is at most one pair of parallel walls successively crossed by solution trajectories of a system of the form [START_REF] Andrianantoandro | Synthetic biology: new engineering rules for an emerging discipline[END_REF].

In other words, there is at most one direction i such that opposite walls, of the form x i = θ - i and x i = θ + i , are crossed. Moreover, such an i is characterised, see [START_REF] Farcot | Geometric properties of piecewise affine biological network models[END_REF], by the condition

∀j = i, τ i (θ - i ) < τ j (θ - j ), (16)
under the assumption I out = I + out (which simplifies the description without loss of generality), i.e. all exiting walls occur at higher threshold values, of the form θ + i , which is thus the threshold involved in the definition of τ i , Eq. ( 5). This allows us to prove the following result: Proposition 1. Suppose that the parameters of ( 14 (1 -γ y θ y )

1 γy > φ 1 -θ 2 1 φ 1 -θ 1 1 1 γ 0 1
Then there the steady state in box D 101 attracts the whole state space of system [START_REF] Farcot | Limit cycles in piecewise-affine gene network models with multiple interaction loops[END_REF].

Proof. Since each box is either containing an asymptotic steady state, or has all its trajectories escaping it toward a focal point, all limit set must be contained in a strongly connected component of the transition graph, i.e. a collection of cyclic paths sharing some vertices. A visual inspection of the transition graph displayed in [START_REF] Farcot | Periodic solutions of piecewise affine gene network models: the case of a negative feedback loop[END_REF] shows that any cyclic path in the transition graph TG must visit the state 100. This state has only two successors: 200 and 101, the fixed state.

Hence it has two exit walls, which we denote by:

W + 1 = {θ 2 1 } × (θ 0 2 , θ 1 2 ) × (0, θ y ) and W + y = (θ 1 1 , θ 2 1 ) × (θ 0 2 , θ 1
2 ) × {θ y }. All other walls are incoming. Denoting them by obvious analogy with the two exiting walls, let us consider each of them. First, both walls W ± 2 are repelling: this has already been said for W + 2 when discussing auto-regulatory terms in [START_REF] Elowitz | A synthetic oscillatory network of transcriptional regulators[END_REF]. For W - 2 , this follows from the fact that D 100 contains only trajectories escaping in finite time, and can be extended to this wall by continuity. Moreover, it follows from TG that any trajectory escaping W ± 2 either reaches D 101 (where the fixed point lies), or enters D 100 again via the wall W - 1 . Thus, any trajectory which does not enter D 101 must cross the pair W ± 1 in succession. Now, from Lemma 1, among the two pairs of walls W ± 1 , W ± y , only one can be crossed in succession by trajectories. Moreover, the inequality in the statement is the exact translation of the condition [START_REF] François | Core genetic module: The mixed feedback loop[END_REF], in the case where W ± y is the crossed pair of walls, precluding any attractor but the known fixed point, φ(101). Some elementary calculus shows that the left-hand side in the inequality of proposition 1 is a increasing function of γ y when γ y ∈ (0, 1/θ y ), as assumed previously. Thus, this inequality is equivalent to requiring a lower bound to γ y , eventhough this bound does not have a simple explicit form. This fact can be given an intuitive explanation: γ y is inversely proportional to the characteristic time of the variable y, in each box. Hence, proposition 1 means that the dynamics of y must be fast enough in order to retrieve the behaviour of the static feedback control, which corresponds to the limit of an instantaneous feedback. See Figure 3 for a numerical example.

The results of this section can be summarised as Proposition 2. A system of the form [START_REF] Elowitz | A synthetic oscillatory network of transcriptional regulators[END_REF], with structural constraints [START_REF] Farcot | Geometric properties of piecewise affine biological network models[END_REF], has a unique, stable and globally attractive limit cycle in absence of input, i.e. u = 0. There exists a control law ensuring a unique, stable and globally attractive equilibrium point. This control can be achieved in two ways:

•) Using a scalar piecewise constant feedback u, such that u(a) satisfies ( 13) for a ∈ {10, 20}.

•) Using dynamic feedback with a single additional variable y, as in ( 14), whose decay rate satisfies the condition in proposition 1, and with υ a solution of (13).

INRIA

Example 2 : birth of a limit cycle

Let us now consider the following system

   ẋ1 (t) = K 1 s + (x 2 ) -(γ 1 1 u + γ 0 1 ) x 1 ẋ2 (t) = K 3 2 s -(x 3 ) + K 1 2 s -(x 1 , θ 2 1 ) -(γ 1 2 u + γ 0 2 ) x 2 ẋ3 (t) = K 3 s + (x 1 , θ 1 1 ) -(γ 1 3 u + γ 0 3 ) x 3 (17)
where s + (x i ) abbreviates s + (x i , θ 1 i ) for i = 2, 3. We assume the following constraints to be satisfied

γ 1 i > 0, i = 1, 2, 3, K 1 > θ 2 1 γ 0 1 > θ 1 1 γ 0 1 K 3 > θ 3 γ 0 3 , K i 2 > θ 2 γ 0 2 , i = 1, 3. ( 18
)
This system is a particular case of two combined negative feedback loops, of the form:

3 2 1
Since the behaviour of a single loop is well characterised by theorem 2, it can be considered as one of the simplest systems whose behaviour might be worth investigating.

Computing the focal points of all boxes, with the abbreviations

φ i = K i /(γ 1 i u + γ 0 i ) (
0 0 φ 1 φ 1 φ 1 0 0 0 φ 1 φ 1 φ 1 φ + 2 φ + 2 φ 3 2 φ + 2 φ + 2 φ 3 2 φ 1 2 φ 1 2 0 φ 1 2 φ 1 2 0 0 φ 3 φ 3 0 φ 3 φ 3 0 φ 3 φ 3 0 φ 3 φ 3
Under the indicated parameter constraints, the following transition graph is easily deduced:

TG(0) =
The region with bold arrows -i.e. the whole graph in this case -is invariant, and we see that depending on the parameter values, the actual solutions of ( 17) may have various behaviours: to each periodic path in TG(0), a stable periodic orbit may possibly correspond, and there is an infinity of such paths. Some examples have already been provided of such situations, where periodic paths of arbitrary length can be realised as stable limit cycles, by suitable choice of parameters [START_REF] Gedeon | Attractors in continuous time switching networks[END_REF]. Although it does not present a fixed box, it may also have a stable equilibrium, limit of damped oscillations, as will be illustrated soon.

In order to guarantee that the system oscillates, we fix the following objective: [START_REF] Gouzé | A class of piecewise linear differential equations arising in biological models[END_REF] RR n 7130 [START_REF] Grilly | A synthetic gene network for tuning protein degradation in saccharomyces cerevisiae[END_REF] for a ∈ A ⋆ . The problem is thus reduced to the satisfiability of the inequality between the two extreme terms above. This fact holds for some parameter values satisfying constraints [START_REF] Glass | Combinatorial and topological methods in nonlinear chemical kinetics[END_REF], and the results of this section can be summarised as Proposition 3. A system of the form [START_REF] Gedeon | Attractors in continuous time switching networks[END_REF], with structural constraints [START_REF] Glass | Combinatorial and topological methods in nonlinear chemical kinetics[END_REF], may present a large variety of asymptotic behaviours without input, i.e. when u = 0. This includes steady states, as shown Figures 4 and5 , as well as limit cycles (not shown). If u is a scalar piecewise constant feedback, such that u(a) satisfies ( 21) for a ∈ A ⋆ , and u(a) = 0 elsewhere, then there exists a unique, stable and globally attractive limit cycle.

TG ⋆ =

Conclusion

We have given, and illustrated by two examples, a control methodology to make unique stable limit cycles appear or disappear in hybrid PWA systems. The obtained feedback laws are termed qualitative control because they depend only on a qualitative abstraction of the original system; its transition graph. Future work suggested by this study are mostly related to the question of dynamic feedback. Actually, the first example shows the effective possibility of using an additional variable to control a system, i.e. to design a controller system to be coupled to the original one. Moreover, the design of this dynamic feedback relied in a simple way on the static feedback problem. This technique should be formalised in more general terms, and applied to other examples in the future. 4, vs time. The blue curve stops before the red one for the following reason. The numerical solutions are computed using the transition map, from wall to wall, and the uncontrolled trajectory crosses successive thresholds within time intervals tending to zero.

Figure 2 :

 2 Figure 2: Dashed lines: with feedback control. Plain lines: without. Two initial conditions, (0.95, 0.95) in box 21 (blue curves) and (0.85, 0.15) in 20 (red curves). The controlled and autonomous trajectories only diverge in box 10, 20 where the feedback is active. See parameters in Appendix A.1

1 :

 1 ) satisfy, denoting φ 1 = K1 γ 0

Figure 3 :

 3 Figure 3: Dashed line: inequality of proposition 1 satisfied. Plain line: inequality violated. Two common initial conditions, (x 1 , x 2 , y) = (0.95; 0.95, 0.1) (blue) and (0.95, 0.95, 0.95) (red). The value of y has been divided by 10 to keep all variables in [0, 1]. In both cases a limit cycle is controlled into an equilibrium point. See parameters in Appendix A.2

Figure 4 :Figure 5 :

 45 Figure 4: Blue, dashed line: without feedback control, spirals towards a fixed points. Red, plain line: with control, tends to a limit cycle. Common initial condition (x 1 , x 2 , x 3 ) = (0.95, 0.95, 0.1). See parameters in Appendix A.3

RR n 7130

Farcot & Gouzé

7 Acknowledgements J.-L.G. was partly funded by the LSIA-INRIA Colage and the ANR MetaGenoReg projects. E.F. was partly funded by the ANR-BBSRC Flower Model project.

INRIA

We now see that the invariant region in bold is a cycle with no escaping edge. Furthermore, it lies in the region where s -(x 1 , θ 2 1 ) = 1, and it appears from (17) that the system in this region is a negative feedback loop involving the three variables. Hence, from theorem 2, we can conclude that there exists a unique stable limit cycle, which is globally attractive, as deduced from TG ⋆ . Yet, it remains to state the inequalities defining this graph. They follow from the inversion of arrows in contact with some a ∈ A ⋆ = {110, 210, 111, 211, 001, 101, 011}, which leads to    