Anne Siegel

Ovidiu Radulescu

Michel Le Borgne

Philippe Veber

Julien Ouy

Sandrine Lagarrigue

P V Eber

Qualitative analysis of the relation between DNA microarray data and behavioral models of regulation networks

Keywords: systems biology, qualitative analysis, steady state shift, sign algebra

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Systems biology: models and data. The eld of systems biology appeared as a response to increasing need for analytical approaches in molecular biology. Its goals include modeling interactions, understanding the behaviour of a system from the interplay of its components, confronting the prediction of the model to data, and inferring models from data. Solutions to these challenges are often interdisciplinary.

The dynamical framework includes simulations and prediction of behaviours models can be either qualitative or quantitative, as reviewed in [START_REF] De | Modeling and simulation of genetic regulatory systems: A literature review[END_REF][START_REF] Chaves | Robustness and fragility of boolean models for genetic regulatory networks[END_REF][START_REF] King | On the use of qualitative reasoning to simulate and identify metabolic pathways[END_REF]. A rst approach m a k es use of continuous models: the concentrations of products are modeled by continuous functions of time, governed by di erential equations. This framework allows one to state biological properties of networks, eventually by using simulation software [START_REF] Bakker | Glycolysis in bloodstream from trypanasoma brucei can be understood in terms of the kinetics of the glycotic enzymes[END_REF][START_REF] Eisenthal | Propsects for antiparasitic drugs: the case of trypanasoma brucei, the causative agent of african sleeping sickness[END_REF][START_REF] Tyson | Sni ers, buzzers, toggles and blinkers: dynamics of regulatory and signaling pathways in the cell[END_REF]Mendes, 1997 Tomita et al., E-mail address: anne.siegel@irisa.fr Tel: (33)299847448 Fax: (33)299847171. 1999). The properties of continuous models can be studied with convex analysis, linear and nonlinear control techniques [START_REF] Fell | Understanding the Control of Metabolism[END_REF]Heinrich a n d Schuster, 1996[START_REF] Papin | Comparison of network-based pathway analysis methods[END_REF][START_REF] Angeli | Detection of multi-stability, bifurcations, and hysteresis in a large class of biological positive-feedback systems[END_REF]. Stochastic models transform reaction rates into probabilities and concentrations into numbers of molecules, allowing to understand how noise in uences a system [START_REF] Rao | Control exploitation and tolerance of intracellular noise[END_REF][START_REF] Kaern | Stochasticity in gene expression: from theories to phenotypes[END_REF]. Finally, in the discrete models, each component is assumed to have a small number of qualitative states, and the regulatory interactions are described by discrete functions. Relevant discrete frameworks can be boolean (Kauman, 1993 Sanchez and[START_REF] Sanchez | A logical analysis of the drosophila gap-gene system[END_REF], logical [START_REF] Karp | Eco-cyc: Encyclopedia of escerichia coli genes and metabolism[END_REF][START_REF] Reiser | Developing a logical model of yeast metabolism[END_REF], or Petri networks [START_REF] Matsuno | Hybrid petri net representation of gene regulatory network[END_REF]Chaouiya e t a l . , 2004). The bridge between continuous and discrete models is made by piecewise linear di erential models (de Jong et al., 2004 Ghoshn and[START_REF] Ghoshn | Symbolic reachable set computation of piecewise a ne hybrid automata and its application to biological modelling: Delta-notch protein signalling[END_REF].

Each of these methods addresses in complementary ways dynamical properties such as the existence of attractors (limit cycles or steady states) and the behavior of these with respect to changes in the parameters (Thomas, 1973 Soul e, 1 Steady state shift experiments and microarray data. Qualitative data such as DNA microarrays data cannot be easily used in most of the frameworks described above for two main reasons. First, the model-based identi cation approach has di culties to take into account the errors and the variability that commonly a ect measured expression leve l s i n D N A microarrays.

Secondly, time series data is not easily available and in many situations (for instance disease studies on clinical tissues) microarrays provide static data, meaning that they inform more on steady state shifts under perturbations than on the dynamics of the system.

In this paper, we develop a mathematical framework that allows us to check the compatibility between gene expression and metabolite concentrations di erential data and a graphical model for the interactions among the measured products (metabolites and genes). If incompatibility is found we propose corrections to data or to the model. Our mathematical results connect network topology and the response to steady state shift experiments. Steady state shift experiments are useful tools in chemistry allowing in principle to recover the reaction mechanisms [START_REF] Chevalier | Toward a systematic determination of complex reaction mechanisms[END_REF]. We argue that similar approaches are well adapted to di erential microarray experiments which compare gene expressions between two di erent states.

In our approach DNA microarrays are interpreted as qualitative data: we only consider the sign of the di erence in the expression levels among the two situations. Additionally, we can consider qualitative data on metabolites, provided by biochemical measurements.

Biological literature is rich i n k n o wledge about interactions between molecules. This information is scattered in many publications and it is not globally compared to the results of large scale experiments. A reason is the di culty to translate the biological knowledge into models: this requires further information on chemical kinetic parameters and mechanisms.

In our models, we do not need such details. We simply use the biological information to build a labeled oriented graph (interaction graph). The nodes of this graph represent molecules (such a s mRNA, proteins or metabolites). The edges represent interactions, labeled by their sign. Edge signs can be interpreted dynamically: a \+" on an edge between A and B means that an increase of the concentration of A increases the rate of production of B. A \;" on an edge stands for a decreasing e ect. The necessary information to build such a graph can be collected either from interaction databases or can be manually extracted from the literature [START_REF] Calvano | A network-based analysis of systemic in ammation in humans[END_REF].

Restricting ourselves to a qualitative framework allows us to compose the interaction rules, but it also tolerates redundancies. Particularly, we can include in the model information on indirect interactions, that is, the modeler can add an edge for A to B as soon as he knows that variations of A a ect the variations of B, even if the underlying mechanism is not known.

Qualitative analysis We present a formal analysis of interaction graphs based on a mathematical model, that allows to confront the interaction graph and qualitative d a t a .

We suppose that the dynamics of the system can be described by a system of di erential equations. An experiment is modeled as a steady state shift of the di erential dynamical system. We perform transformations on the di erential system that allow to connect the variations of products in a linear system. The quantitative system generalizes the discrete Laplace equation on graphs. Like for the Laplace equation, we solve the linear Dirichlet problem [START_REF] Chung | Discrete green's functions[END_REF][START_REF] Soardi | Potential theory on in nite networks[END_REF][START_REF] Campanino | Random walks in randomly oriented lattices[END_REF]. We also solve the non-linear Dirichlet problem: this means relating the interior values to the values on the entrance boundary of the oriented graph even when the products variations are large.

We further suppose that experimental information on products is qualitative, that is, we can nd out if the concentration of a product has increased, it has decreased, or it has not changed signi cantly. In order to exploit this information, qualitative models are derived from the interaction graph and from the di erential dynamics: we transpose the quantitative equations into a linear qualitative system in the sign algebra (Kuipers, 1994 Trav e-Massuy es andDague, 2003). Provided that the signs of the interactions are constant for states within the range of experimental variations, the linear qualitative system applies to both small and large variations of the products.

We have developed several complementary methods to analyze the system of qualitative equations and to test whether experimental data are solutions of the system.

The graph valuation algorithm is an automated method to partially or totally solve the qualitative system of equations. When there is no solution, incompatibilities between the model and the data are detected, as well as the place where these incompatibilities occur. When a solution exists, the algorithm predicts the signs of variations of products that have not been measured.

The essential balance c omputation allows us to detect the variables in the system that are in uenced by competitive pathways and to nd which pathway wins the competition. Eventually, this computation detects the presence of errors in the data.

The modeling and the algorithms will be illustrated on a simpli ed model of regulated lipogenesis in liver. Then, we apply the algorithm to an extended model of lipogenesis. With these algorithms, we are able to validate or correct data and models. Also, they serve to emphasize products and pathways that have biological dynamic interest. Further improvements of the algorithms such as coding of qualitative equations over Galois elds and applying them to models containing hundreds of products will be presented elsewhere [START_REF] Veber | Complex qualitative models in biology: a new approach[END_REF]. This paper is organized as follows. In Section 1 we describe our working example which is the regulation of lipogenesis in liver. Section 2 is devoted to modeling assumptions. We introduce the mathematical concepts and we s h o w how the qualitative equations can be obtained. We detail in Section 3 two methods to solve the system of qualitative equations derived in Section 2. Section 4 is devoted to an independent method for analyzing incompatibilities and competitions. Section 5 presents the application of the algorithm to an extended model of lipogenesis.

Working example: regulation of the synthesis of fatty acids

Gene regulation associated with fatty acid synthesis in liver is our working example. The corresponding interactions are intricate and involve hundreds of molecules. By way of illustration we have kept as nodes of our illustrative model only the most important biological molecules and their intermediates. The choice of the observed molecules was not optimized we shall see how our analysis suggests improvements on it.

Two w ays of production of fatty acids coexist in liver. Saturated and mono-unsaturated fatty acids are produced from citrates due to the presence of a metabolic pathway composed of four enzymes, namely ACL (ATP citrate liase), ACC (acetyl-Coenzyme A carboxylase), FAS (fatty acid synthase) and SCD1 (Stearoyl-CoA desaturase 1). Polyunsaturated fatty acids (PUFA) such as arachidonic acid and docosahexaenoic acid are synthesized from essential fatty acids provided by n utrition D5D (Delta-5 Desaturase) and D6D (Delta-6 Desaturase) catalyze the key steps of the synthesis of PUFA.

PUFAs play pivotal roles in many biological functions. They regulate the expression of genes that impact on lipid, carbohydrate, and protein metabolism. The e ects of PUFA are mediated directly or indirectly as shown in the following examples. Direct e ects are due to bindings leading to changes in the trans-activating activity o f n uclear receptors (PPAR { p e r o xisome proliferator activated receptors, LXR { Liver-X-Receptor , HNF-4 ). Indirect e ects result in changes in the abundance of regulatory transcription factors (SREBP-1c { sterol regulatory element binding-protein{, ChREBP, etc.) [START_REF] Jump | Fatty acid regulation of gene transcription[END_REF].

Variables in the model

We have considered in our model the transcription factors PPAR , LXR and SREBP-1c (denoted by P P AR, LXR, SREBP), as they are synthesized from the corresponding genes. We have included the transactivating active forms of these nuclear receptors: LXR-a (denoting a complex LXR :RXR ), PPAR-a (complex PPAR :RXR ) and SREBP-a (cleaved form of SREBP-1c). We h a ve a l s o c o nsidered SCAP (SREBP cleavage activating pro-tein), a key enzyme involved in the cleavage of SREBP-1c. SCAP interacts with another family of proteins called INSIG, showing the complexity of molecular mechanism. We have included PUFA to symbolize metabolites. Finally, w e h a ve considered the enzymes ACL, ACC, FAS, SCD1 (implied in the fatty a c i d s y n thesis from citrate) and D5D, D6D (implied in PUFA s y n thesis).

Interactions in the model Relations between the variables are as following. SREBP-a is an activator of the transcription of ACL, ACC, FAS, SCD1, D5D and D6D [START_REF] Nara | The e-box like sterol regulatory element mediates the suppression of human delta-6 desaturase gene by highly unsaturated fatty acids[END_REF][START_REF] Jump | Fatty acid regulation of gene transcription[END_REF]. LXR-a is an activator of the transcription of SREBP and FAS. It also indirectly activates ACL, ACC and SCD1 (Ste ensen and [START_REF] Ste Ensen | Putative metabolic e ects of the liver x receptor (lxr)[END_REF]. These indirect actions are kept in the model because we don't know whether they are only SREBP-mediated.

PUFA activates the formation of PPAR-a from PPAR. It inhibits the formation of LXR-a from LXR. It also inhibits the formation of SREBP-a (by inducing the degradation of mRNA and inhibiting the cleavage) [START_REF] Jump | Fatty acid regulation of gene transcription[END_REF]. SCAP represents the activators of the formation of SREBP-a from SREBP. SCAP is inhibited by P U F A.

PPAR directly activates the transcription of SCD1, D5D, D6D [START_REF] Miller | Peroxisome proliferators induce mouse liver stearoyl-coa desaturase 1 gene expression[END_REF][START_REF] Tang | Regulation of human delta-6 desaturase gene transcription: identi cation of a functional direct repeat-1 element[END_REF][START_REF] Matsuzaka | Dual regulation of mouse delta(5)-and delta(6)-desaturase gene expression by srebp-1 and pparalpha[END_REF].

The activation of SCD1, D5D and D6D by both SREBP and PPAR is paradoxical because fatty acid synthesis (partially governed by SREBP) and oxidation (partially induced by PPAR) are antagonistic in liver. Nevertheless, PUFA h a ve a regulatory role in oxidation. Hence, the induction of D5D and D6D gene by PPAR could be a compensatory response to the increased PUFA demand caused by induction of fatty a c i d o xidation.

Working data set: virtual fasting protocol

Our goal is to test whether data can t properly with the full knowledge on the system. Such d a t a can be an experimental set that have to be analyzed. It could also simply be a set of hypothetical data whose coherence have to be tested. Here, our working data set is inspired by the results of fasting extracted from literature as described hereafter, modi ed so that a potential error is in-troduced.

A compilation of recent literature on lipogenesis regulation indicates that SREBP, A CL, ACC, FAS and SCD1 decline in liver during the fasted state in rodents [START_REF] Liang | Diminished hepatic response to fasting/refeeding and liver x receptor agonists in mice with selective de ciency of sterol regulatory element-binding protein-1c[END_REF] this state is characterized by an inhibition of fatty acid synthesis and an activation of the fatty acid oxidation. However, Tobin et al [START_REF] Tobin | Crosstalk between fatty a c i d a n d c holesterol metabolism mediated by l i v er x receptor-alpha[END_REF] showed that fasting rats for 24h increased the hepatic LXR mRNA and [START_REF] Matsuzaka | Dual regulation of mouse delta(5)-and delta(6)-desaturase gene expression by srebp-1 and pparalpha[END_REF] observe no di erence in either the hepatic D5D or D6D mRNA level between fasted and reefed mouse livers. Moreover, PUFA l e v els can be considered to be increased in liver following starvation because of the important lipolysis from adipose tissue as shown by Lee et al in mice after 72h fasting [START_REF] Lee | Requirement of pparalpha in maintaining phospholipid and triacylglycerol homeostasis during energy deprivation[END_REF].

So we de ne our working data as results from a virtual experimentation on fasted animals compared to fed state and on which the hepatic mRNA of di erent genes could be quanti ed for example by DNA microarray analysis and the variation of hepatic PUFA could be measured by b i o c hemical analysis. Hence, we assume that ACL, LXR, PPAR and PUFA increase while SREBP, ACC, FAS and SCD1 are supposed to decrease, in response to the fasting. We also assume that D5D and D6D are unchanged. By doing this, we v oluntarily introduce a data that does not t with results appearing in the literature: ACL should decline instead of inclining. This error is introduced to test whether the analyze is able to notice it and possibly correct it or not.

Steady state shift: qualitative description

We justify here the approach and impose some applicability limits. The basic mathematical object is a qualitative graph called the interaction graph. This graph represents the present k n o wledge on interactions for the biological process under study. Experimental data give variations of the concentrations of some molecules between two external conditions. The corresponding nodes of the interaction graph will be called observed nodes.

Di erential model The model describes the in-teractions between a set of molecules indexed by a set I = f1 : : : n g. The concentration of the molecule indexed by i is X i . The quantity dXi dt denotes the speed of variation of X i .

We assume that the various concentrations involved in the model evolve according to a di erential dynamics. The set of control parameters of the dynamics is denoted by P. The state vector X = ( X 1 : : : X n ) satis es a di erential equation dX dt = F(X P): 2.1. Interaction graph Qualitative graph The set of nodes I f E g gathers all the molecules I in the model together with an extra generic node E to represent the exterior world.

The set of edges A (I f E g ) (I f E g ) contains oriented edges joining j to i whenever @F i @X j 6 = 0 . This description implies that the interaction graph contains no multi-arcs. Additionally, the set A contains edges from the exterior world E to I specifying which are the entrance nodes of the system. If there is an edge from E to i 2 I, the node i is potentially a ected by n o d e s that are not included in the model. If there is no edge between E and i, then the node i is supposed to be in uenced only by I. Specifying these edges is an important modeling hypothesis.

The interaction graph M = (I f E g A s ) is signed by the sign function s(j i) = sign( @F i @X j ) 2 f+ ;g and weighted by the interaction coe cient a(j i) = @F i @X j . Notice that the sign function s(j i)

is state dependent (it depends on the variables X and P). We suppose that s(j i) The only di erence between E and P is dynamical. Indeed, E is rigorously a dynamical variable, while P is xed during dynamics. From a static perspective E and the parameters P represent variables that may shift the steady state of the model. Throughout this paper we shall consider that there is no direct e ect of P on the nodes of G: all in uences arrive v i a the exterior E. This modeling hypothesis applies to the analyzed example. It can be lifted with minor changes. If the modeler judges this assumption to be not reasonable, he can just add a (non-dynamical) variables P and connect them to the nodes that are directly in uenced.

Extracting an interaction graph from biological facts described in the literature Intuitively, an arc j ! i describes the alteration of the speed of variation of the product i when the concentration of j is modi ed. We build the interaction graph associated to a given biological question from the biological literature. More precisely, w e consider an edge from j to i in the graph when experimentation in a paper shows that the production of i is modi ed after a change in the concentration of j.

The interaction graph is supplied with negative self-interaction on all nodes of the graph, implying s(i i) = ; for any i (in the absence of additional positive self-regulation). These result for instance from degradation processes or from growth induced dilution. It will be shown (Theorem 2.6 of section 2) that negative self-interaction ensures the existence of a steady state.

Such a n i n teraction graph might be incomplete: it is not necessary to gather all the biological facts in a model. Actually, we intend to understand whether experimental data ts with the theoretical model or not.

Example: fasting protocol (Fig. 1) For the biological problem used as an illustration, we consider a model M = (I f E g A s ) with vertices I = fPUFA, LXR, LXR-a, PPAR, PPARa, SREBP, SREBP-a, SCAP, ACL, ACC, FAS, SCD1, D5D, D6Dg. The interaction graph M is shown in Fig. 1. It is obtained by translating the biological information given in Sec. 1. As an example, since the active form of LXR (that is, LXR-a) is an activator of the transcription of SREBP, w e consider an edge from LXRa to SREBP, labeled by \+". We also consider a positive edge from SREBP to its active form SREBP-a since increasing the concentration of SREBP should increase the concentration of its active f o r m . There is a positive edge from LXRa to ACL because of the activation of ACL by the active form of LXR. Notice that this action is taken into account e v en if it is indirect: no detail on this action exists in the literature a possible mechanism for this action involves SREBP as an intermediate, as described by the path LXR-a ! SREBP ! SREBP-a ! ACL. However, since we don't know whether the action of LXR on ACL is only SREBP-mediated, the modeler has decided to keep an additional edge from LXR to ACL, redundant with the path containing SREBP.

Paths and loops on graph A path j i is a sequence of nodes fi 1 i 2 : : : i p g such that (i k i k+1 ) 2 A, i 1 = j and i p = i. The nodes of a path are visited just once. Let the path inuence a j i denote the product of interaction coe cients along a path a j i = Q k=p;1 k=1 a(i k i k+1 ): The sign of the path is s(j i) = sign(a j i ).

A l o o p i s a p a t h w i t h i d e n tical terminal nodes

i 1 = i p . In a loop all nodes are visited just once, excepting the terminal nodes.

Analyzed subgraph, observed nodes The subset of the nodes in the model on which we perform the analysis is denoted by G I. The nodes in G inherit the connections from M.

The set of observed nodes will be denoted as O. One can de ne the nodes variations X i i 2 O, given by experimental data on O. The sign of these variations is denoted by s( X i ).

In The entrance boundary of G, denoted in G is formed by all the nodes of G that are supposed to have e n tering connections from the outside:

PPAR [+] PPAR-a PUFA [+] ACL [+] LXR-a SREBP-a SCAP D5D [0] D6D [0] ACC [-] FAS [-] SCD1 [-] LXR [+] EXT SREBP [-]
in G = fi 2 G j 9 j 2 (I f E g ) n G (j i) 2 Ag: The set of interior nodes is denoted G = Gn in G.
Notice that when G = I, the entrance boundary is entirely made of those nodes that have connections from E.

In the model M of our working example, the analyzed subgraph is I. Its entrance boundary is in I = fPPAR, LXR,PUFAg. Predecessors We denote by pred(i) the set of predecessors of i: pred(i) = fj 2 I f E g j 6 = i (j i) 2 Ag. 2.2. Qualitative linear equations Steady state shift In our approach w e suppose that the waiting time after a change of the parameters is su cient for the biological system to reach steady state.

Steady states are characterized by: F(X P) = 0 :

(1) The steady state shift is the result of a change in the control parameters P. In our example fast-ing is the only control parameter. Experimental data are interpreted as variations, between two steady states, of the variables associated with observed nodes. This point of view, which is assumed in the rest of the paper, is not correct in all situations. It cannot be assumed for time series where observations are made at relatively close instants. In such a situation, it is not true that the system has reached steady state at each instant w h e n i t is observed. For our example, this hypothesis is justi ed since data come from metabolism of animals that where kept under normal diet or fasted during a long time before the measurements. Furthermore, we should check carefully if the system is not capable of autonomous oscillations (limit cycles) in which case steady state again can not bereached.

Small variation and self-susceptivity We d i fferentiate Eq.1 in order to understand how the steady state changes for small changes in the parameters. Fo r a n o d e i 2 I one gets the equation: @F i @X i X i + P j6 =i @F i @X j X j + P k @F i @P k P k = 0 :

We can notice that if @F i @X j 6 = 0 and j 6 = i if and only if there is an edge from j to i in the interaction graph, that is, j 2 pred(i). Let us consider a node i such that @F i @P k = 0: the production or the consumption of X i does not depend directly on the parameters.

Property 2.1

Let us consider a system whose steady state depends on the parameters. Let i be a node that satis es the following conditions: the node i is not directly in uenced by t h e parameters (all in uences come via its predecessors), i.e. @F i @P k = 0 , @F i @X i 6 = 0 , then the variation in i c a n b e e n tirely calculated from the variations in the predecessors nodes:

X i = ; @F i @X i ;1 X k2pred(i) a(k i) X k : (2)
Using a mechanical analogy, i = ; @F i @X i can be called self-susceptivity: it is the ratio between the force produced by the predecessors of i and the variation in i. As discussed, there are always e ects that produce self-susceptivity (degradation, growth induced dilution, self-regulation).

The variation X i can be calculated from the force exerted by the predecessors only if the self-susceptivity i is non-zero. A zero selfsusceptivity would represent a non-generic case when the sum of all these e ects cancel exactly.

Sign algebra If we only know the signs of the variations and of the interaction coe cients, we intend to nd the relations between the signs of the variations and the signs of the interaction coe cients. To d o s o , we transcribe Eq.2 into the sign algebra. Remind that s( X i ) s t a n d s f o r t h e sign of the variation of X i .

The sign algebra is de ned as the set S of subsets of f0 + ;g, that is, S = ff0g f+g f;g f0 +g f0 ;g f+ ;g f0 + ;gg:

We d e n o t e b y ? = f0 + ;g the subset standing

for an undetermined sign.

The set S is then provided with the sum and product rules de ned on singletons as follows:

f+g + f;g =? f+g + f+g = f+g f+g + f0g = f+g f;g + f;g = f;g f;g + f0g = f;g f+g f;g = f;g f+g f +g = f+g f;g f;g = f+g f+g f 0g = f0g f;g f0g = f0g:

These rules are extended to S by the union rule:

(A B) + C = ( A + C) (B + C).
An equality in an equation means that the sign corresponding to the l.h.s. (left hand side) and the sign of the r.h.s. (right hand side) have a non empty intersection (Kuipers, 1994 Trav e-Massuy es andDague, 2003).

Generalization to large variations Eq.2 represents a di erential constraint that should be satis ed by (small) variations of the concentrations X i . Naturally, w e m a y ask whether there is also a constraint among large variations. This is equivalent to asking whether there is any globally de ned function relating the value X i to the concentrations of its predecessors that we d e n o t e b y the vector X(i) . If there is no direct in uence of the parameters P on X i , stationarity o f X i reads:

F i (X i X(i) ) = 0 (3)
We h a ve the following:

Property 2.2 If @F i @X i < ;C C > 0, if there is no direct in u- ence of P on X i and if F i (0 X(i) ) > 0 then there is a function i : R ni + ! R (n i being the num- ber of predecessors of i) s u c h that X i = i ( X(i) ).
The di erential of i is d i = ;( @F i @X i

) ;1 X j2pred(i) @F i @X j dX j :

Now, let us suppose that the steady state shifts from the vector X 1 to the vector X 2 and that the projection of the steady state onto the predecessors of i shifts from X(i 1) to X(i 2) . In order to obtain the variation of X i it is enough to integrate the di erential of i on a segment C 1 2 connecting X(i 1) X(i 2) :

X i = Z C 1 2 ;( @F i @X i ) ;1 X j2pred(i)
@F i @X j dX j :

(4)

Qualitative system of equations Transcribing Eq. 4 in the sign algebra provides the following result.

Theorem 2.3

Supposing that the conditions of Property 2.2 are satis ed and supposing that the signs of @F i @X j ( i ( X(i) ) X(i) ) are constant on the segment C 1 2 , then the signs of nite variations of steady state variables are connected by: s

( X i ) X k2pred(i) s(k i)s( X k ):
(5)

Remark Although the ends of the segment C 1 2 are steady states there is no reason for the interior points to ful ll stationarity equations (they fulll only one of the stationarity equations, namely Eq.3). Nevertheless, supposing that the absolute values of the derivatives @F i @X j are bounded (which is a natural modeling assumption) it follows from Eq.4 that the excursion of i and therefore of X i is bounded when X(i) describes the segment C 1 2 . This means that the segment C 1 2 corresponds to a range of states reasonably close to the experimental range, thus justifying the modeling assumption concerning the signs. This reasoning warns against the naive i n terpretation of interaction signs that assumes constancy of these signs over arbitrary ranges of the concentrations. The problem is not trivial and should be considered carefully.

Example Theorem 2.3 provides the following equations for the working example M (for ease of notations, s( A) is denoted as A). One can apply Theorem 2.3 because the self-susceptivities are positive in the absence of positive self-regulation.

Nodes i that satisfy the hypothesis of Theo- rem 2.3 are those which are not connected to the exterior world E (implying that the rst condition is satis ed). In order to analyze data, we need to transform the qualitative system derived in Theorem 2.3 to obtain a system where the variables are functions of a xed set of variables. In general, this kind of system transformation in the sign algebra is a NP-complete problem [START_REF] Dormoy | Controlling qualitative resolution[END_REF].

In the following paragraphs we show that quantitatively this operation is equivalent to solving the Dirichlet problem for a discrete operator on the interaction graph. A solution to the Dirichlet problem can be found under mild conditions. We use this result for deriving a large class of qualitative systems.

Transmission of in uences on graphs One can read Eq.2 as the propagation of in uence on the interaction graph: the concentration variation of the molecule i is an average of the in uences from its predecessors. Let us de ne as follows the discrete operator L : R I ! R I :

(L X) i = X i + @F i @X i ;1 X k2pred(i) a(k i) X k : (6)
Then Eq.2 reads L X = 0 which is analogous to the Laplace equation on graphs. Let us remind that the Laplace equation ful lls an analogous property: the value of an harmonic function (satisfying the Laplace equation) in a node is the average of the values of the function on the predecessors of this node. The analogy can be pushed further as one may try to de ne the in uence of the boundary on an interior point (question known as the Dirichlet problem) [START_REF] Chung | Discrete green's functions[END_REF].

Laplace equations on graphs were intensively studied in connection to electrical networks and random walks on undirected graphs [START_REF] Soardi | Potential theory on in nite networks[END_REF]. Let us notice that our model describes propagation of in uence on oriented graphs. This could be related to random walks on oriented graphs. Nevertheless, the propagation opera-tors resulting from random walks have negative non-diagonal elements [START_REF] Campanino | Random walks in randomly oriented lattices[END_REF], while in our case there are no such constraints on the signs of the interactions. Little is known on the properties of random walks on oriented graphs. For instance it is not known whether the solution to the Dirichlet problem has a representation in terms of stopped random walks or not. The connection between propagation of in uences and random walks on oriented graphs could thus be pro table both for biology and for applied mathematics.

In uence of the boundary on the interior From now, we x a subnetwork G I of the model. We have already de ned the entrance boundary in G, formed by all the nodes of G that are supposed to have entering connections from the outside. The set of interior nodes was denoted G = G n in G.

By di erentiating the components of Eq.1 corresponding only to the interior nodes of G one gets for every i 2 G: P j2 G @F i @X j X j + P k2 in G @F i @X k X k + P l2InG @F i @X l X l + @F i @P P= 0 :

Noticing that @F i @X l = 0 for all i 2 G , l 2 I n G and assuming that @F i @P = 0 for all i 2 G it follows: 8i 2 G X j2 G @F i @X j X j = ; X k2 in G @F i @X k X k :

(7) Eq.7 is the starting point for the proof of the following theorem (details can be found in the Appendix):

Theorem 2.4 (Linear Dirichlet solution)

Let G be a subgraph of I. Let A be the restriction of the jacobian of F to the interior nodes: A ij = @F i @X j i j 2 G. Let us suppose that: det( A) 6 = 0 , there is no direct in uence of the parameters on the interior nodes of G, i.e. @X i @P k = 0 , 8i 2 G.

Then, the response of i 2 G to small changes on the boundary of G is given by:

X i = X j2 in G X j i2P G a j i C j i X j (8)
with the following notations: P G denotes the set of paths included in G, starting on the boundary and that do not return to the boundary: P G = fj i = (i 1 = j : : : i k : : : i p = i), i k = 2 in G, 8k > 1g a j i = Q p;1 k=1 a(i k i k+1 ) C j i = (;1) l k(j) i +1 det A det A fk(j) ig c denotes the \path modulus" of the internal path k(j) i where { k(j) 2 G is the second node after j of the path j i { l k(j) i is the length of the path k(j) i, i . e . the number of arcs

Moduli and sign algebra

Path modulus Eq.8 describes the in uence on an interior node of variations on the boundary. j acts with a force a j i X j on the node i, along the path j i. This force is bigger when the product of interaction coe cients @F i k+1 @X i k along the path is bigger.

C k(j) i is the ratio force/response and therefore can be called path modulus. According to Theorem 2.4, the contribution to the modulus comes only from the internal part of the path: k(j) i. If the internal part reduces to a point (k(j) = i, meaning that the path from j to i is a direct arc), then the contribution to the modulus comes from the node i only: the path modulus becomes node modulus. A large path (or node) modulus implies a small response at the end of the path, even if the force is big. Therefore, the modulus can be related to sensitivity.

Signs of moduli One way to test the compatibility b e t ween micro-array data and interaction models is to compare predicted and observed signs of the variations. In order to do so, we n e e d to know the signs of the path moduli. The following property i s v ery useful: Property 2.5

If the interior of the subgraph G contains no positive loops, then all internal path moduli C k(j) i and node moduli C i in Theorem 2.4 are strictly positive.

If G has no positive self-loops, than all the selfsusceptivities i are positive.

Following a reasoning similar to the one leading to Property 2.2 and to Theorem 2.3, and using Property 2.5 and Theorem 2.4, we can write down (under some conditions) linear qualitative equations valid for large variations.

Let X G , XG be the set of variables internal to G, and on the entrance boundary, respectively. If the parameters P have direct in uence on the interior of G, then the stationarity equations for G (non-linear Dirichlet problem) read:

F G ( X G XG ) = 0 (9)
where F G is the restriction of F to G.

Let us suppose that steady state shifts from X 1 to X 2 and that the projection of the steady state onto the boundary variables shifts from X1 to X2 .

Then, we obtain (see Appendix) the following:

Theorem 2.6

Suppose that the parameter variations have no direct inuence on the interior of G, i.e. @F i @P k = 0 , 8i 2 G the restriction of the Jacobian of F to G is non-singular the interaction graph has no positive l o o p s in the interior of G the restriction of F to G has the form

F G = G ( X G XG ) ; G X G ( 10
)
where G is bounded and satis es G i (: : : X i = 0 : : : XG ) > 0 8i 2 G (11) and G is a diagonal matrix with positive entries Then there is a function G such that X G = G ( XG ) gives the unique solution of the nonlinear Dirichlet problem (Eq.9).

If furthermore the signs s(j i) of the path in uences from the boundary to the interior point i are constant for X G on a segment connecting X1 G X2 G and X G = G ( XG ), then the signs of the nite variations of the internal variables X i satisfy the following relation for every i 2 G: s(

X i ) = X j2 in G X j i2P G s(j i)s( X j ): (12)
Remark The decomposition in Eq.10 corresponds to a regulation term G and to a linear dissipation term G X G . The condition Eq.11 naturally con nes the dynamics to positive c o ncentrations. A remark similar to the one for Theorem 2.3 holds as to the boundedness of the variations of X G when the boundary variables describe the segment from X1 to X2 .

Model and data compatibility

In this section we introduce two methods of qualitative analysis of biological systems and experimental data. These methods use the systems of qualitative equations generated by T h e orem 2.3. Nevertheless, Theorem 2.6 can be used to generate other systems of qualitative equations. Indeed, it is enough to consider any c o vering of the network I by subgraphs fG j g j=1 p . Then, the Eq. 12 written for all nodes i 2 G j and for all subgraphs G j provide another system of qualitative equations.

Hence, even if algorithms for incompatibility are presented for system 1 in this section, they can be easily adapted for another system of qualitative equations as well.

Solving qualitative systems

Let us recall the qualitative equations provided by Theorem 2.3 System 1

(1) PPAR-a = PPAR + PUFA (2) LXR-a = -PUFA + LXR (3) SREBP = LXR-a (4) SREBP-a= SREBP + SCAP -PUFA (5) ACL= LXR-a + SREBP-a -PUFA (6) ACC= LXR-a + SREBP-a -PUFA (7) FAS= LXR-a + SREBP-a -PUFA (8) SCD1 = LXR-a + SREBP-a -PUFA + PPAR-a (9) SCAP= -PUFA (10) D5D= PPAR-a + SREBP-a -PUFA (11) D6D = PPAR-a + SREBP-a -PUFA PPAR=+, PUFA=+, LXR=+, SREBP=-, SCD1=-, FAS=-, ACL=+, ACC=-, D5D=0, D6D=0

Compatible system A valuation of the qualitative system is called compatible if it satis es:

1. all nodes are determined (have values different f r o m ? ) 2. all the observed nodes are given the observed values 3. the qualitative Eqs. 5 are ful lled. Solving the qualitative system means nding all compatible valuations. If there is no compatible valuation, the system is said to be incompatible.

The interest of nding compatible valuations lies in the following property: any solution of the quantitative system provides a compatible valuation. Consequently, the incompatibility of the qualitative system implies that the quantitative system is incompatible as well. This immediately indicates a contradiction between the model and the observations. Furthermore, the signs corresponding to a solution of the quantitative system should be given by one of the solutions of the qualitative system.

Let us notice that when the system is compatible, the number of compatible valuations can be big and this could bring very limited information on the biological system. In general we would like to know which kind of information one can extract from compatible valuations. In order to do this an useful concept is the following notion of competition.

Competitions We say t h a t there is a c o m p e t ition or a balance on a node i if the right side of the Eq.5 is equal to \?". The biological interpretation of this indetermination is that the node i is submitted to competitive actions of di erent signs. An example is provided by the node LXR-a in the working example, that satis es PUFA = + LXR = + LXR-a = ; PUFA + LXR = ?:

Determining the sign of the balanced node says which paths win the competition. If the variation of the balanced node is zero we s a y that the competition is neutral. If this variation is negative ( p o s i t i v e) we s a y that the negative (positive) paths win the competition, or that the competition leans toward the negative (positive) direction. The information on competitions has biological interest, because it replaces undetermined actions by determined ones and facilitates the understanding of the interaction graph.

In the next subsections we use simple rules of qualitative algebra to solve System 1 and extract information on competitions. Then, we propose an algorithmic approach t o s o l v e the system obtained from Theorem 2.3 by using these rules systematically. Algorithmic methods related to competitions are detailed in Section 4.

Hand solving of qualitative equations

L e t u s s o l v e System 1 by using the rules of the sign algebra and the following computation rule: Property 3.1

If A = X + B and if the sign of A is opposite to the sign of B, then the sign of X is determined and it is equal to the sign of A.

The result of the application of these rules to System 1 is summarized as following: Property 3.2 System 1 is not compatible with the observed variations.

Simple change of the sign of ACL (ACL=;) renders it compatible and determined. The values of non-observed variables are PPAR-a=+, LXR-a=;, SCAP=;, SREBP-a=;.

In this case there are competitions on LXRa, SCD1 that are observed to lean towards the negative directions. There are also competitions on D5D, D6D that are observed to be neutral. There is no competition on FAS, ACC and ACL.

Proof From Eq.1 PPAR-a=+. From Eq.3 LXR-a=SREBP=;. From Eq.9 SCAP=;. From Eq.4 SREBP-a=;.

We have already some biologically interesting results: by Eq.2, LXR-a is supposed to be undetermined, that is, the result of a competition between the positive action of LXR and the negative action of PUFA. Using Eq.3 we h a ve s h o wn that the competition on LXR-a should lean toward the negative direction. From Eq.8 it follows that there is a competition on SCD1. Data on SCD1 implies that this competition leans toward the negative direction.

From Eq.5 and Property 3.1, SREBP-a=+. This is inconsistent with the previously found value, which p r o ves that the system is incompatible. There is necessarily an error, either in the interactions or in the data. Let us assume that the interactions are correct. Then, considering that data on ACL is wrong is the simplest way t o make the system compatible (because Eqs.6 and 7 do not force signs of FAS and ACC).

With the correction ACL=;, the unknown values are determined in the same way as before. From Eq.1 PPAR-a=+. From Eq.3 LXR-a=SREBP=;. From Eq.9 SCAP=;. From Eq.4 SREBP-a=;. Now Eq.5 is satis ed with no competition. The same is true for Eqs.6,7: there are no competitions on ACC and FAS. From Eqs.2,8 we nd that there are competitions on LXRa,SCD1 which leans towards the negative direction. From Eqs.10,11 we nd that there are competitions on D5D, D6D. The observed variations of D5D, D6D are vanishing, therefore these competitions are neutral.

Manual resolution of a qualitative system is only possible for small examples. In the next section we describe an algorithm that allows the systematic application of these rules to large systems.

Graph valuation algorithm

Basic ideas underlying the algorithm A natural way to reduce the number of unknowns in a system of equations is to proceed by elimination. This means to compute the value of one unknown from an equation and then substitute this value in all the other equations. Unfortunately, in the sign algebra, the usual elimination rules do not apply (qualitative equality of undetermined quantities is not transitive).

In the proposed algorithm, we will rely on the rules used in Section 3.2. Let us formulate them more precisely.

Eq. 5 implies that the sign of X i is determined if s(k i)s( X k ) are determined and have the same value for all the predecessors k of the node i. This rst elementary rule is used in the forward propagation step of the algorithm. The second rule follows from Property 3 . 3 . Suppose that

X i = s(k i)s( X k ) + P l6 =k s(l i)s( X l ):
If both the signs of P l6 =k s(l i)s( X l ) a n d s( X i ) are determined and if they are opposite then s(k i)s( X k ) = s( X i ). The second rule is used in the backward propagation step of the algorithm. In the case of a sign algebra with 0 value, the rst rule is unchanged while the second one needs a straightforward adaptation.

The graph valuation algorithm 1. Initialization. All the nodes are given a value in f+ ; ?g. Observed nodes are ini- tialized with the value corresponding to the observed variations of the concentration of the molecules. The other nodes are initialized with the value ?. The value associated with a node i is denoted val(i). Let us consider two sets of nodes U and T.

The nodes with values in f+ ;g are put in T. The remaining nodes (with value \?") are put in U.

Forward p r opagation:

Create a set F that contains every node of U that is also a successor of a node in T.

Repeat: Find a node i in F which has at least a predecessor and such that for all predecessors j of i, s(j i)val(j) = +, or for all predecessors j of i, s(j i)val(j) = ;. Then set val(i) = s(j i)val(j), remove i from F, remove the successors of i which a r e i n U and put them in F put i in T.

until no node found 3. Backward p r opagation.

Create a set B that contains every node of U that is also a predecessor of a node in T.

Repeat: Find a node j in B which has a successor i in T, s u c h that val(i) P k2pred(i) k6 =j s(k i)val(k) = ; or P k2pred(i) k6 =j s(k i)val(k) = 0 . Then val(j) = s(j i)val(i), remove j from B, remove the predecessors of j which a r e in U and put them in B, put j in T.

until no node found 4. Alternate propagations. Alternate forward and backward propagations until no new value is added to T.

Notice that each time that a node is valuated, it is moved to T and it never returns to F or B.

T is an increasing subset of a nite set, so that the algorithm stops.

Application to the detection of incompatibilities All the valuations introduced by t h e a l g orithm are unique sign choices compatible with the equations provided by Theorem 2.3. By doing so, no compatible valuation is lost. More precisely:

Property 3.3

The graph valuation algorithm provides a new graph such that any compatible valuation of the nodes in the initial graph is a compatible valuation of the new graph. Nevertheless, nothing guarantees that a compatible valuation exists, i.e. that all Eqs. 5 are satis ed. The resulting graph may contain incompatibilities. These incompatibilities are detected simply by c hecking Eq. 5 o n e a c h node.

Criterion to compare model and data 1

If the graph output by the algorithm described below c o n tains a node where Eq. 5 is not satised, then experimental data and the model are incompatible.

Working example

In the example M the initialization is PUFA=+, LXR=+, LXR-a=?, PPAR=+, PPAR-a=?, SREBP=-, SREBP-a=?, ACL=+, ACC=+, FAS=-, SCD1=-, SCAP=?, D5D=?, D6D=?. We assume that the variations of D5D and D6D are not known.

The forward propagation begins with The backward propagation begins with B=f LXR-a g and ends with B = , LXR-a=-(from SREBP=-), or LXR-a=+ (from ACL=+, SREBP-a=-and PUFA=+). The second forward propagation begins with F=fD6D, D5Dg and stops with F = but no new value for D5D or D6D.

After running the algorithm, if LXR-a=-equation ( 5) is not satis ed at node ACL and if LXR-a=+ equation ( 3) is not satis ed at node SREBP. This shows a contradiction between the model and the values assigned to some nodes by experimental data or hypothesis.

The graph valuation algorithm is not complete for the resolution of linear qualitative systems. As a consequence, the algorithm concludes only when it detects an incompatibility or when all the values are determined.

Assessment of competitions

The graph valuation algorithm aims at detecting incompatibilities or producing compatible valuations. Nevertheless, even when these tasks are successful the result may be di cult to use. If for instance data is wrong for a node that happens to be very connected to the rest of the graph, then many incompatibilities appear and it is difcult to localize and to correct the error. Also, if there are not enough observations then compatible valuations are produced in huge number. It is di cult to extract useful biological information from this situation.

In this section, we develop another approach, that exploits the global qualitative information on the propagation of in uences on the graph. This approach is complementary to the graph valuation algorithm since it identi es competition processes and the places where these occur.

Competitions As already explained, we s a y t h a t there is a competition or a balance on a node i if the right side of the Eq.5 is equal to \?". The biological interpretation of this indetermination is that the node i is submitted to competitive actions of di erent signs. This de nition can be extended to any equation connecting i to other nodes in the graph instead of Eq.5 (related to the predecessors of a node i).

Such equations are provided by Theorem 2.6 that can be applied to every subgraph G containing i.

Then, we say that there is a competition on a node i when there exists an equation whose left side is s( X i ) and whose right side is equal to \?".

Set of equations describing the sign of variation of a variable X i0 Let M = ( I f E g A s )

beanetwork with no positive l o o p i n t h e i n terior of I and such that the sign function s is constant within the entire range of parameter variations.

Let us focus on a variable X i0 with index i 0 in the model M. From Theorem 2.6, as soon as X i0 is an internal variable of a subgraph G of M, the sign of variation of X i0 can be expressed as a combination of the signs of paths starting on the boundary of G and never coming back t o i t :

8G I X i0 2 G s( X i0 ) X j2 in G j i 0 2P G s(j i 0 )s( X j ): (13)
The set of equations describing the variable X i0 is the set of equations ( 13) for all subgraphs G such t h a t X i0 is an internal variable of G.

Example of SREBP In the working example of lipogenesis, let us consider the variable X = SREBP. W e w ant to determine all the subgraphs G that contain X as an internal variable. A necessary and su cient condition is that G contains all the predecessors of X, that is, LXR-a. In each subgraph G, SREBP is determined by all the paths that start from the boundary of G and never go back to it. Notice that any path arriving in SREBP contains LXR-a. If it is longer than the path LXR-a ! SREBP, i t m ust contain PUFA or LXR, that both belong to the boundary of G (since they are connected to the exterior). Hence, there are three paths to SREBP that start from the boundary of such subgraphs G and never come back t o i t . They are given in Table 4, together with the set of equations associated to SREBP.

Paths appearing in an equation (13) Let P(M) contain all paths of M that appear in an equation ( 13). Hence, a path j i belongs to P(M) if there exists a subgraph G I such t h a t j belongs to the boundary of G and j i is an internal path of G. The reasoning detailed in the example allows to characterize this set. A full proof is given in the Appendix. Property 4.1 A path j i I belongstoP(M) i pred(j) 6 = and no node of j i except possibly j is connected to the exterior E.

If so, there exists a subgraph G such that the variations of X i satisfy:

s( X i ) = s( X j )s(j i)

+ P k 2 in G k i 2 P G k i 6 = j i s(k i)s( X k ):
Example of SREBP-a Let us x X =SREBPa. As an application of Property 4.1 seven paths start from the boundary of a graph G such t h a t X 2 G and the path never come back to the boundary. They are given in Table 2 together with the equations describing SREBP-a. Notice that Eq. 19, 20 and 21 correspond to eq. 16, 17 and 18 in which s(P 7 )s( PUFA) has been replaced by s(P 10 )s( SCAP). This corresponds to the particular case of SCAP, whose only predecessor is PUFA, so that SCAP can be considered either on the boundary or in the interior of a subgraph G, p r o viding two di erent equations.

In uences Property 4.1 suggests the following de nition. For a path j i 2 P (M) w e s a y t h a t s( X j )s(j i) i s the in uence of j i on the variable X i .

Let us also de ne the sets of positive and negative in uences acting on X i , denoted by I + (i) and I ; (i). I + (i) = fj i 2 P (M) j 2 O s( X j )s(j i) = + g I ; (i) = fj i 2 P (M) j 2 O s( X j )s(j i) = ;g The set of all in uences is denoted by I(i).

I(i) = I + (i) I ; (i):

PUFA [+] LXR-a LXR [+] EXT SREBP [-] (P1) LXR-a!SREBP, (P2) PUFA ! LXR-a ! SREBP, (P3) LXR! LXR-a!SREBP. s( SREBP) (;) = s(P1 )s( LXR-a) (14) s( SREBP)
(;) = s(P2 )s( PUFA) (;) + s(P3)s( LXR): (+) (15)
Table 1 System of qualitative equations describing the variations of SREBP sign of in uences associated to each path. The graph is restricted to nodes having an in uence on SREBP such that this in uence possibly starts from the boundary (PPAR, PUFA, LXR) but does not come back t o i t .

Examples of SREBP and SREBP-a. The observed variations are PUFA= +, SREBP = + and LXR= +, so that the positive and the negative in uences are I + (SREBP) = fP 3 g I ; (SREBP) = fP 2 g I + (SREBP-a) = fP 9 g I ; (SREBP-a) = fP 4 P 5 P 6 P 7 g:

The in uence of paths starting from LXR-a and SCAP is not taken into account since the variations of these molecules are not observed.

Counterbalanced in uences In Table 2, the right sides of Eqs. 15, 16 and 18 are undetermined since they are the sums of at least two terms of di erent signs. This means that there are competitions on SREBP and SREBP-a. The observation of SREBP shows that the competition between P 2 and P 3 leans toward the negative direction: we say that the in uence of P 3 is counterbalanced by the in uence of P 2 . Since SREBP-a is not observed, we can say that the inuence of P 9 is counterbalanced by the in uence of P 4 P 5 P 6 P 7 , but we d o n o t k n o w the inclination of the competition.

This example justi es the following:

De nition 4.2

The set of counterbalanced in uences on i is denoted by C(i). It is de ned as follows:

if i is observed, the in uence of j i is counterbalanced if it is di erent from the variation of X i , i.e. if s( X j )s(j i) 6 = s( X i ). Then if s( X i ) = + then C(i) = I ; (i) if s( X i ) = ; then C(i) = I + (i)
if s( X i ) = 0 then C(i) = I(i): If i is not observed, the in uence of j i is counterbalanced if there exists another path k i whose in uence on X i has a di erent sign, i.e. if I + (i) and I ; (i) are both nonempty. Then

{ if i 6 2 O , I + (i) 6 = and I + (i) 6 = then C(i) = I(i), { if i 6 2 O , I + (i) 6 = or I ; (i) = then C(i) = .
A direct consequence of the de nition is that an in uence which is not counterbalanced absorbs all the in uences in the equations it appears in, so that it does not introduce indetermination's: Property 4.3 If a p a t h j i is not counterbalanced then the sign of its in uence on the node i is the same as the sign of X i . Table 2 System of qualitative equations describing the variations of SREBP-a sign of in uences associated to each path.

To summarize, in uences are counterbalanced when they compete with other in uences. When we can observe w h i c h is the result of this competition, the counterbalanced in uence is the dominated one. In uences are not counterbalanced either if they are non-competing, or if they dominate their competitors.

Examples

In the example of SREBP-a, positive in uences I + (SREBP-a) = fP 9 g compete with negative ones I ; (SREBP-a) = fP 4 P 5 P 6 P 7 g for determining the unknown sign of the variation of SREBP-a. All these in uences are considered to be counterbalanced: C(SREBP-a) = I(SREBP-a)

= fP 9 P 4 P 5 P 6 P 7 g:

The observed negative variation of SREBP implies that positive in uences are dominated, hence counterbalanced: C(SREBP) = I + (SREBP) = fP 3 g: Redundant and essential in uences Let us consider the in uences of P 1 and P 2 on SREBP. The observed sign of SREBP implies that P 1 is dominated (thus counterbalanced) by P 2 . With no loss, we can skip P 1 in our analysis of in uences. We can not skip P 2 . We say that P 2 is essential, while P 1 is redundant. Let us now consider the in uences of P 9 and P 6 on SREBP-a.

These two in uences are competing, but we do not know w h i c h one wins. We decide to consider that both are redundant.

Generally, we call essential those in uences that can not be skipped without loss. The discussion of the examples justi es the following:

De nition 4.4

A path is said to have an essential in uence on a node i if it belongs to the set E(i) de ned as follows: E(i) = fj i = ( j 0 = j j 1 : : : j l = i) 2 I (i) 80 < k < l j j k 6 2 C (k)g:

Essentially balanced nodes Let us consider LXR-a. It collects two e s s e n tial in uences of different signs (from PUFA and LXR). We s a y t h a t LXR-a is essentially balanced. Let us also consider SREBP. This node collects two competing in uences: P 1 and P 2 . Neverthe-less, we know that P 1 is dominated, hence redundant. We s a y that SREBP is not essentially balanced.

More generally, w e h a ve the following:

De nition 4.5

An observed variable X i is said to be essentially balanced if it collects essential in uences that are also counterbalanced, i.e. if E(i) \ C (i) 6 = .

A non observed variable X i is said to be essentially balanced if it collects essential in uences with opposite signs, that is E(i) \ I + (i) 6 = and E(i) \ I ; (i) 6 = .

An essential balance on a variable X i gathers competing in uences such that no information allows to say that one of the competing in uences has been previously absorbed in a previous competition. This implies that the essential balance is localized precisely at the node i. The variation of X i gives the inclination of the competition.

From the above discussion it follows:

Criterion to compare model and data 2

Let us consider a model of molecular in uences.

Let us consider that the interaction graph contains no positive loop made of nodes that are not connected to the exterior. Let X i be a variable that collects essential in uences, i.e. E(i) 6 .

Then several situations are possible:

1. X i is essentially balanced. Then there is a non-redundant competition process localized at i. If X i is observed, we c a n s a y which one of the competing in uences is dominating.

2. X i is not essentially balanced and all the collected essential in uences have t h e s a m e sign.

If X i is observed then:

(a) if the sign of X i coincides with the one of the essential in uences then we can say nothing.

(b) if the sign of X i is di erent from the one of the essential in uences then there is probably an error in data.

Complementarity with the graph valuation algorithm As already stated, the assessment o f competitions is not devoted to proving incompatibilities. For instance the case 2b) of the Criterion 2 suggests, but does not prove the presence of an error. Its utility is not obvious in our simple working example, but should appear from its application to more complex networks. Case 2b) of Criterion 2 is a hint for the location of the error.

Case 1) of Criterion 2 is obviously useful as it allows to localize non-redundant competitions. These non-redundant competitions are usually important regulation checkpoints.

Working example and biological discussion

Computation on the graph related to the working examples M of regulation of lipid synthesis showed that there are competitions on 9 nodes among 11: LXR-a, SREBP, SREBP-a, ACL, ACC, FAS, SCD1, D5D, D6D. However, some of these competitions are redundant. Among the 32 counterbalanced paths, 20 are essential. Essentially balanced nodes are LXR-a, ACL and SCD1, D5D and D6D.

LXR-a collects non-redundant competitive paths. ACL collects four essential negative in uences and no essential positive in uence. This suggests that the data on ACL is false.

The system has to be studied more precisely, in order to check whether there is a real incompatibility b e t ween the model and the data. A w ay to do this is by using the methods of the previous section. An alternative w ay is a re ned assessment of competitions (see the paragraph improvements). SCD1, D5D and D6D are observed nodes that collect essential negative and positive in uences. This suggests that a nonredundant competition occurs on each of these nodes.

The observed negative variation on SCD1 suggests that the negative path through SREBP dominates the positive path through PPAR.

Consequently, t wo p o i n ts can be extracted from such an analysis: at the SCD1 level, the negative pathway containing SREBP dominates the positive pathway containing PPAR. At t h e A CL level, the data on ACL should be checked carefully or paths in the model leading to ACL should be studied in detail.

Improvements Our discussion of competitions aims at identi ng and localizing essential competitions. Let us describe modi cations of the de nitions that allow us to re ne the analysis.

De nitions 4.2, 4.4 declare redundant in uences that compete with others, even if the result of the competition is unknown. We h a ve did so in order to localize essential competitions: once they occur somewhere, competitions are transported along the graph as non-essential.

With altered de nitions we can identify and localize incompatibilities. Let us change the De nition 4.2 of counterbalanced in uences by skipping its second part. This means deciding that in uences of di erent signs on an unobserved node are not counterbalanced. Then, according to De nition 4.4, the redundant in uences are only those that are dominated somewhere. Hence the case 2b) of the Criterion 2 identi es an incompatibility: it is impossible to act with non-dominated in uences of even sign on a node and to obtain a di erent sign on the node. This is actually the case of the node ACL in our working example. This node collects non-dominated negative in uences and one dominated positive in uence.

Thus, by altering De nitions 4.2, 4.4 we c a n d etect incompatilities, but we can no longer localize competitions. In order to satisfy both purposes we can introduce two types of redundancy, one suitable for incompatibility detection, the other for localizing competitions.

Application to an extended model of the synthesis of fatty acids

We applied the graph valuation algorithm and the essential balance computation to an extended model of the synthesis of fatty acids shown in Fig. 2. The model contains 59 products and 110 interactions implied in the synthesis metabolism. The product list includes ACL, ACC, D5D, D6D, FAS, LXR, PPAR, PUFA, SCAP, SCD1 that were introduced in the working example given in Fig. 1. Compared to the working example, this new model includes neither the active forms of SREBP and LXR, nor SCAP. The model was obtained as a compilation of interactions extracted from the litterature. Notice that in this model, only the nodes without predecessors (that is, INS and retinoid) are connected to the exterior. Hence, we suppose implicitely that the model is complete enough to explain all the possible actions among products.

We have used a theoretical experimental set of measurements for 16 products, built as a result of a virtual experimentation on fasted animals. Hence, we suppose that PUFA, PPAR, LXR, ACL increase while ACC, D6D, FAS and SREBP decrease.

Backward-forward algorithm The backwardforward algorithm states that the system of qualitative equations is incompatible with the set of data. The qualitative equation that cannot be satis ed by the data is the equation corresponding to ACL.

Essential balances

The result of the computation of essential balances is shown in Table 3. Six products are essentially balanced, among which D5D, D6D, SCD1 and SREBP. In the working example, an essential balance occured on the active form of LXR, denoted LXR-a it cannot occur in this model since LXR-a is not considered any more. Hence, the essential balance on LXR-a on the working example is shifted in the new model as a balance on the successor of LXR-a, that is, SREBP.

16 products collect only positive pure in uences. Four of them were supposed to be positive during the experimentation, which is con rmed by the analysis. The other 12 products collect only positive in uences, meaning that either they are balanced, either they incline in the positive direction.

19 products collect only negative pure in uences. Among them, ACL, FACE and FDFT1 were supposed to incline in the positive direction during the experimentation. We suspect that there is an error in the data.

Nothing can be said about the 17 products that do not collect eseential in uence.

Backward-forward algorithm on a corrected set of data The backward-forward algorithm on the set of data stated that the set of data was not compatible with the model. The computation of essential balances suggested that the data on ACL, FACE and FDFT1 were false. Hence, we propose as a new set of data the following set: +] PUFA, CPT1A, LXR, RAET1C, PPAR -] D6D, SREBF1, ACC, FAS, HMGCR, LDLR, ACL, FACE, FDFT1

Then the backward-forward algorithm terminates so that this set of data is compatible.

Finally, the backward-forward algorithm states that the inclinaison of nonobserved values suggested by the essential balanced algorithm is compatible with the interaction graph. +] PUFA, CPT1A, LXR, RAET1C, PPAR, ACS, ABCA1, D5D, SCD1, ABCD3, APOE, CD5L, ADAM11, CYP7A1, FBXO3, FABP1, HMGCS2, LPL, MTP, UCP2, VLDL TG -] D6D, SREBF1, ACC, FAS, HMGCR, LDLR, ACL, FACE, FDFT1

Remarks and conclusion

In this paper we have discussed how a qualitative theoretical model of a mixed gene and metabolic network can be confronted to DNAmicroarray and metabolic data.

Using a di erential description of dynamics we have obtained quantitative equations connecting the theoretical variations of genes and metabolic products to the changes of the external conditions. In their linear version, that connects small variations of genes and products, these equations are analogous to the discrete linear Laplace equations on graphs that have been studied in the context of Markov c hains [START_REF] Chung | Discrete green's functions[END_REF][START_REF] Soardi | Potential theory on in nite networks[END_REF][START_REF] Campanino | Random walks in randomly oriented lattices[END_REF] 3 Computation of essential in uences on the extended model of the synthesis of fatty acids shown in Fig. 2. The computation suggests that competitions occur in D6D, SREBP, A CS, ABCA1, D5D and SCD1. Moreover, the competition should incline in the negative direction in D6D and SREBP. The computation also suggest that the data about ACL, FACE and FDFT1 have to be corrected. Finally, the computation suggests that ACS, ABCA1, D5D and SCD1 either incline in the positive direction, or collect a competition process. Similarly, ABCD3, APOE, CD5L, ADAM11, CYP7A1, FBXO3, FABP1, HMGCS2, LPL, MTP, UCP2 and VLDL TG either incline in the negative direction, or collect a competition process. order to calculate the variations we h a ve t o c o mpute the matrix inverse of the Jacobian of the di erential dynamical system. This matrix inverse is in fact the discrete Green function of the problem. Contrary to [START_REF] Chung | Discrete green's functions[END_REF] we h a ve not searched for explicit formulas for the Green function (di cult to obtain for irregular graphs), but we have emphasized how one node collects in uences from the other nodes. A similar path representation of the Green function is used in statistical mechanics [START_REF] Campanino | Random walks in randomly oriented lattices[END_REF].

We have extended the quantitative equations to the case of large variations of the products and we h a ve established general linear qualitative equations that are valid both for small and large variations. The validity of the approach r e l i e s o n uniform signs of the interaction coe cients (elements of the Jacobian) within the range of the experimental conditions.

We h a ve i n troduced the graph valuation algorithm that can reduce or solve the qualitative system. The computation rules are not those of usual linear spaces, but belong to the sign algebra.

When the qualitative system is incompatible, it is generally di cult to propose corrections. When it is compatible, we m a y obtain many solutions. For this reason we h a ve a l s o i n troduced a complementary approach which is the study of competitions. This approach treats all nodes on equal footing. It does not prove incompatibilities, but it identi es biologically important competition processes and the places where they occur.

The analysis of the examples taken from the regulation of lipogenesis suggests that this type of methods could be used to optimize experimental protocols. The optimization criteria are not formalized here. Nevertheless, as a general rule, if a node is observed we w ould like t o h a ve observations as complete as possible on its predecessors. In the working example, SCAP is not observed. Instead, we h a ve predicted the sign of its variation. The experimental knowledge of SCAP variation would further facilitate the localization of the discovered incompatibility between interaction and data.

Using Eqs.23,24 it follows X i = ; P k2 G det A fig c det A @F i @X k X k ; P j2 G j6 =i P j i P k2 G (;1) lj i a j i det A fj ig c det A @F j @X k X k i :

We can identify in the sum the contributions of direct (one arc) and indirect (at least two arcs) paths from the boundary to the interior. Thus,

X i P k2 in G h P indirect k i a k i C k(j) i + P direct k i a(k i) Ci i X k
from which i t f o l l o ws Eq.8.

Proof of Theorem 2.6 The existence of a solution of the nonlinear Dirichlet problem Eq.9 is a rather standard mathematical result based on Eqs.10,11 [START_REF] Radulescu | Topology and linear response of interaction networks in molecular biology[END_REF]. The uniqueness of such a solution follows from the absence of positive loops in the interior of G (Thomas's conjecture proven in [START_REF] Soul E | Graphic requirements for multistationarity[END_REF]). Thus, the unique solution of Eq.9 can be represented as a function X G = G ( XG ). The di erential of G is the solution of the linear Dirichlet problem, given by Eq.8. Integrating this di erential on the segment connecting X1 X2 and passing to sign algebra we nd Eq.12. We also need the fact that all internal moduli are positive which b y Property 2.5 follows from the absence of positive l o o p s .

Proof of Property 2.5 If there are no positive loops, then it can be shown (see Soul e (Soul e, 2003)) that the signs of the principal minors of the Jacobian det A are alternating. Thus, sign(det A fjg c ) = (;1) # G;1 sign(det A (j i) c ) = (;1) # G;lj i ;1 sign( A jj ) = sign(a jj ) = ;1 from which it follows that all the moduli are positive.

Proof of Property 4.1 If P = j i = (j j 1 : : : i ) cont a i n s a n i n termediate node j k that is connected to the exterior, then for every subgraph G, j k belongs to the boundary of G, hence P = 2 P G .

Conversely, i f P contains no intermediate node connected to the exterior E let G = fj 0 = j j 1 : : : j n = ig 0<k n pred(j k )g. Then P 2 P G , j 2 in G, hence P 2 P (M).

In order to prove the second part of the Property, Theorem 2.6 can be applied to G. 3.

  our example, we consider G = I. Only some nodes are observed. The variations of observed nodes are shown in Fig.1. O = fPPAR LXR ACL ACC SREBP FAS SCD1 PUFA D5D D6Dg: Entrance boundary The analyzed subgraph G is naturally partitioned into two s u b s e t s : the entrance boundary of G and the interior nodes.

Figure 1 .

 1 Figure1. Interaction graph A for a model M = ( I A s) of regulation of the synthesis of fatty acids. The node EXT stands for E and represents the exterior world. Self-regulation loops on nodes are omitted for sake of clarity. Observed variations are s( PPAR) = +, s( PUFA) = +, s( LXR) = +, s( SREBP) = ;, s( SCD1) = ;, s( FAS) = ;, s( AC L ) = + , s( ACC) = ;, s( D5D) = 0, s( D6D) = 0.

  PPAR-a = PPAR + PUFA (2) LXR-a = -PUFA + LXR (3) SREBP = LXR-a (4) SREBP-a= SREBP + SCAP -PUFA (5) ACL= LXR-a + SREBP-a -PUFA (6) ACC= LXR-a + SREBP-a -PUFA (7) FAS= LXR-a + SREBP-a -PUFA (8) SCD1 = LXR-a + SREBP-a -PUFA + PPAR-a (9) SCAP= -PUFA (10) D5D= PPAR-a + SREBP-a -PUFA (11) D6D = PPAR-a + SREBP-a -PUFA Known values are the following: Values 1 PPAR=+, PUFA=+, LXR=+, SREBP=-, SCD1=-, FAS=-, ACL=+, ACC=-, D5D=0, D6D=0 2.3. In uences and their transmission across a boundary

{

  det A fk(j) ig c is the principal minor obtained by eliminating all the lines and all the columns whose indices are in the path k(j) i from the jacobian A. Conventionally, if the resulting mi- nor is empty w e c hoose it equal to one: det A G c = 1 . If k(j) = i the path modulus becomes the node modulus C i = ; det A det A fig c . If we consider G = fig pred(i) t o b e t h e s e t of a predecessors of a node i we obtain the Property 2.1 as a corollary.

F

  = fPPAR-a LXR-a SREBP-a SCAP D6D D5Dg and ends with F = , PPAR-a =+ (from PUFA =+ and PPAR=+), SREBP-a=-(from SREBP=-) and SCAP=-(from PUFA=+). D5D and D6D are removed from F and keep their value ?.

  PUFA!SREBP-a, (P5) SREBP!SREBP-a, (P6) PUFA! LXR-a!SREBP!SREBP-a, (P7) PUFA!SCAP!SREBP-a, (P8) LXR-a!SREBP!SREBP-a, (P9) LXR! LXR-a!SREBP!SREBP-a, (P10) SCAP!SREBP-a.

Figure 2 .

 2 Figure 2. Interaction graph for an extended model of regulation of the synthesis of fatty a c i d s . The model contains 59 nodes and 110 interactions. 16 products are supposed to be observed. The list of products is given in Table3.

Acknowledgments This research was supported by the ACI IMPBio, French Ministry for Research program on interdisciplinarity. We would like to thank M. Crouzeix, I. Gruais, Y. Guivarch, D. Petritis, and C. Soul e f o r m a n y inspiring discussions.

Appendix: Proofs Proof of Property 2.2 From @F i @X i < ;C it follows that for any xed X(i) , F i (X i X(i) ) depends monotonically on X i . Furthermore, F i becomes negative for large enough X i . As it is positive for X i = 0 , b y the intermediate value theorem it vanishes somewhere. Eq.3 has an unique solution X i = i ( X(i) ) for any X(i) . The di erential of i follows from the implicit function theorem.

Proof of Theorem 2.4 From Eq. 7 it follows:

where B = A ;1 .

Let us consider the following well-known property of the inverse matrix [START_REF] Bloom | Linear algebra and geometry[END_REF]:

where M ij is the minor obtained by deleting line i and column j in matrix A.

Furthermore, it is a simple exercise [START_REF] Bloom | Linear algebra and geometry[END_REF] to show that:

where j i is any path leading from j to i with no loop, l j i is the length of the path, a j i is the product of elements of A along this path, A (j i) c

is the principal minor de ned by the set of indices complementary to those in the path. Conventionally, det A G c = 1 (determinant of the empty ma- trix).