In systems biology, large biochemical reaction networks can either be represented by bipartite graphs or by systems of ordinary differential equations. Modelers want to determine the existence of reductions between those reaction networks. Because, it is not possible to decide this existence with equation systems, a previous thesis [1] focussed on graph structures. A subgraph epimorphism (SEPI) framework was developed and gave results close to biologists' expectations.

Three main difficulties of the SEPI framework have been identified. First, establishing whether two models are linked through a SEPI is complex and computationally expensive. Second, the number of SEPIs found can be huge, making the analysis of SEPI sets between two given graphs very difficult for biologists. Finally, some existing SEPIs do not have a biological interpretation. This diploma thesis led to three combined ways to improve the framework. One way consisted to redefine the decision problem into an optimisation problem to restrict the set of solutions. A second way was to determine, together with biologists, restrictions on one of the framework's operations in order to filter irrelevant reductions. Lastly, a preprocessing step has been introduced, consisting of rewriting graphs according to subgraph isomorphism relations. The impact of these three combined implementations has been evaluated on models of the BioModels database. Results demonstrated that it contributed to make the SEPI framework more relevant, efficient and functional.

Statement of authorship

I hereby certify that I have authored this Diploma Thesis entitled Graph matching, theory and SAT implementation independently and without undue assistance from third parties. No other than the resources and references indicated in this thesis have been used. I have marked both literal and accordingly adopted quotations as such. There were no additional persons involved in the intellectual preparation of the present thesis. I am aware that violations of this declaration may lead to subsequent withdrawal of the degree. Dresden, 30.09.2019 Orianne Laura Bargain

Introduction

Systems biology relies on computational models. It aims to understand behaviours of complex systems thanks to the analysis and integration of interactions between elementary components. It uses a holistic approach to biological research [2]. For example, systems biology studies molecular interactions in cells and considers the cell as a whole [2]. It uses models to store knowledge. Interactions can be described at different scales, with a variety of quantitative or qualitative approaches. Many formalisms have been developed such as reaction graphs, gene networks and process algebra. Proteins are macromolecules. At their scale, signed graphs, which are graphs with positive or negative labels on each edge, depict the influence of proteins. At a more precise level, signed hypergraphs or signed bipartite graphs represent reaction networks, which indicate which chemical species interact together, by means of a common reaction, to produce other species. With stoichiometry and ordinary differential equations (ODE), simulations can be performed. Simulations give a better understanding of some dynamical behaviour of the system and can highlight existence of stationary states and periodic behaviours.

The Michaelis-Menten enzymatic reaction is a classical example in Systems Biology [3], observed with many enzymes. The reaction can be written in several forms. Michaelis-Menten mechanism is a catalytic mechanism with three reactions and four molecular species in its expanded form: enzyme E, substrate S, complex enzyme-substrate C, also denoted by ES, product P . The mechanism describes kinetics of an enzyme-catalysed reaction acting on a single substrate to irreversibly give a product. Its reaction model and mass action law kinetics is

E + S k d kc C → kp E + P
which gives the following system of ODEs where capital letters denote the concentration of corresponding species:

           . S = -k c • E • S + k d • C . E = -k c • E • S + k d • C + k p • C . C = k c • E • S -k d • C -k p • C . P = k p • C
Michaelis-Menten mechanism can be reduced to a single reaction with three species (E, S, P ) and more complex kinetics:

E + S → km•K•S km+S E + P ,
with system of ODEs:

       . S = -km•K•S km+S . E = 0 . P = km•K•S km+S 1. Introduction
These two mechanisms will sometimes be referred to as expanded and reduced forms of Michaelis-Menten reaction. This reduction of the Michaelis-Menten mechanism is valid for some initial concentrations of molecular species. Figures 1-1 and 1-2 show the dynamic behaviour of species concentrations in expanded and reduced Michaelis-Menten mechanisms. Dynamic behaviours of species of interest (E, S, P ) is almost the same in both mechanisms, which shows the relevance of this reduction. A free open-source repository called BioModels [4] was created for storing, exchanging, and retrieving quantitative models of biological interest, which were described in peer-reviewed scientific literature.

But mathematical biology is not suited for reduction on complex reactions. Similarities of mechanism are better revealed by relationships between network structures. Thus, current approaches focus on reaction graphs by making abstraction of stoichiometry and kinetic.

Figures 1-3 and 1-4 are examples of reaction graphs and they show a same mechanism. Reaction graphs are bipartite graphs that represent a biochemical reaction. Species are depicted by round vertices and reactions are depicted by square vertices. An edge represents a link in a reaction. For example, specie E and specie S are reacting together through reaction R 1 to produce species ES.

Models can have a lot of vertices and can be complex, for example, with more than three hundred vertices on a model representing dynamics of extracellular-signal-regulated kinase [5]. Thus, biologists need automated tools to exhibit reduction relations between models. Especially for big repositories like BioModels [4].

1.1. Related work E R 1 S ES R 2 P R -1

Related work

A few current approaches are reducing biochemical reactions using reaction networks, subgraph matching tools or subgraph epimorphisms.

Subgraph isomorphism

A lot of research has been carried out on the subgraph isomorphism (SISO) problem. Graph-Grep [6] and Grafil [7] for example, are graph matching methods using subgraph isomorphism but it gives limited results on reaction graphs.

A common example for limitations of SISO is Michaelis-Menten reaction graph. Michaelis-Menten kinetic is one of the best-known models of enzyme kinetics in biochemistry [3]. Its reaction is presented in figure 123. Figure 1-4 is a reduced way to represent Michaelis-Menten kinetic. SISO relation can be seen as a sequence of delete operations performed on vertices. This relation is too strong because it doesn't find a relation between both graphs. Here, a way to reach the second reaction graph from the first one is to delete reaction R -1 , merge together reactions R 1 and R 2 and delete specie ES.

Graph minor

Graph minors [START_REF] Robertson | Graph minors. ii. algorithmic aspects of tree-width[END_REF][START_REF] Lovász | Graph minor theory[END_REF] is another well studied graph theory concept. It has a lot of properties but because reaction graphs are bipartite graphs, no solution to extend graph minor definition to biochemical models have been found.

Introduction

Substructure index-based approximate graph alignment

Substructure Index-Based Approximate Graph Alignment (SAGA) tool [START_REF] Tian | Saga: a subgraph matching tool for biological graphs[END_REF] was developed to overcome SISO limitations. SAGA is an approximate subgraph matching tool, which permits node gaps, node mismatches, and graph structural differences. A node gap is a node insertion or deletion, node mismatches exhibit similar functionality of specie nodes and graph structural differences allow flexibility in node connectivity relationships.

But SAGA was implemented to query databases and not to compare two given graphs.

Morphisms of reaction networks

Another approach [START_REF] Cardelli | Morphisms of reaction networks that couple structure to function[END_REF] uses influence networks to represent complex biological systems. Influence networks depict the impact of molecular species on each others. Cardelli studies structural aspects of models to infer properties and makes a connection between network structure and behaviour. Cardelli defines morphisms between reaction networks with structural connections. It achieved to determine statically that a complex network imitate a simpler one. Figure 1-5 to figure 12345678are examples. Simplified influence network reproduce kinetics of complex influence network with same initial conditions and reaction rates.

The study [START_REF] Cardelli | Morphisms of reaction networks that couple structure to function[END_REF] shows that graph morphisms provide a structural reason for kinetic similarity, but it doesn't bring a tool to automatically compare biochemical reactions.

Subgraph epimorphism

SISO problem can be viewed as a sequence of delete operations. To tackle SISO's limitations, a solution is to add a merge operation to the delete operation, which is equivalent to a subgraph epimorphism (SEPI) problem. A previous thesis [1] introduced a SEPI framework 1 for model reductions in systems biology, establishing a formal relation between models.

SEPI framework focuses only on reaction graphs but corresponds to a mathematical reduction, since a vertex deletion can be viewed as neglecting a specie and a fusion of two vertices as an aggregation of two species. To obtain one graph from another, a sequence of vertex deletion and fusion is applied. These two simple rules are accurate enough to capture mathematical model reductions and not too complex. When a model is reachable from another model with this framework, it means that these models might be reduced in some valid way and it give insights as of how much a model is similar to another. SEPI is a powerful framework because, without having to add more information to models, it can organise existing models as hierarchies of refined to reduced. Comparison between reaction graphs can be applied systematically to repositories such as BioModels that contain hundreds of models with no given relation.

Motivations

Three main difficulties of SEPI framework have been identified by biologists when using it. First, establishing whether two models are linked through a SEPI is complex and computationally expensive. Second, the number of SEPIs found can be huge, making analysis of SEPI sets between two given graphs very difficult for biologists. Finally, some existing SEPIs do not have a biological interpretation.

Additionally, there is no guarantee that SEPI pairings are related to correct reductions. Up to now, there is no perfect correspondence between valid model reductions and graph operations, partly because valid model reductions themselves are very difficult to justify mathematically. For example, L. Noethen et al. [START_REF] Noethen | Tikhonov's theorem and quasi-steady state[END_REF] uses Tikhonov [START_REF] Tikhonov | Systems of differential equations containing small parameters in the derivatives[END_REF] and Fenichel [START_REF] Fenichel | Geometric singular perturbation theory for ordinary differential equations[END_REF] theorems.

The goal of this diploma thesis is to study further model reductions and to improve SEPI framework. This work uses and refines the graph theoretical framework of subgraph epimorphisms developed by S. Gay, F. Fages, and S. Soliman in [1,[START_REF] Gay | A graphical method for reducing and relating models in systems biology[END_REF][START_REF] Gay | On the subgraph epimorphism problem[END_REF].

Contributions

This diploma thesis led to three combined ways to improve the framework.

1. One way consists to redefine a decision problem into an optimisation problem to select solutions.

2.

A second way is to determine, together with biologists, restrictions on one of the framework's operations in order to filter irrelevant reductions.

3. Lastly, a preprocessing step is introduced, consisting of rewriting graphs according to subgraph isomorphism relations.

Impacts of these three combined implementations are evaluated on models of BioModels database.

Results demonstrate that it contributed to make SEPI framework more relevant, efficient, and functional.

Outline

Approaches used and improved in this thesis have specific notations and uses specific tools, thus Chapter 2 provides context and explanations to readers not familiar with graph theory or SAT solving. Chapter 3 presents a first approach enhancing SEPI framework: the decision problem is redefined into an optimisation problem by adding an objective of maximisation or minimisation of the vertex deletion number. A second approach determines restrictions on a merge operation in order to filter irrelevant reductions, it is demonstrated in Chapter 4. A more functional variant of this merge restriction is explained and implemented in Chapter 5. Chapter 6 introduces a preprocessing step consisting of rewriting graphs according to subgraph isomorphism relations. Chapter 7 tackles the problem of finding a Greatest Lower Bound (GLB) or a Lowest Upper Bound (LUB) between two graphs that are not linked by a SEPI relation. Bounds correspond to a common reduction or a common refinement of both graphs. Finally, Chapter 8 concludes this thesis and discuss some axes of enhancement.

Preliminaries

This chapter assumes familiarity with basic set theory and provides background knowledge for readers who are not familiar with graph theory and SAT solving.

Graph theory

Basics of graph theory will not be presented in this chapter. Readers can refer to the work of Diestel [START_REF] Diestel | Graph theory, volume 173 of[END_REF] or Bang-Jensen and Gutin [START_REF] Bang-Jensen | Digraphs: theory, algorithms and applications[END_REF].

Definition 2.1.1 (Directed graph) A directed graph G is a pair G = (V, A), where V is a set of vertices and A ⊆ V × V is a set of arcs.
When it is not specified, graph means directed graph.

A graph may have loops and cannot have more than one arc from some vertex to another one. The cardinality of a set S is denoted |S|. When not explicitly defined, G and G denote graphs, with G = (V, A) and G = (V , A ).

Operations

The subgraph epimorphism will be based on a merge and a delete operation.

A delete operation removes a vertex and all arcs connected to it.

V = V \ {u} and A = A ∩ (V × V ). G → d G is written when there exists u ∈ V such that d u (G) = G . G → * d G
is written when G can be obtained from G using zero, one or several deletes.

A merge operation removes both initial vertices and creates a new corresponding vertex that inherits arcs connected to them. Definition 2.1.3 (Merge) For all u, v ∈ V , the result of the merge operation m u,v is graph m u,v (G) = (V , A ) where:

V = V \ {u, v} {uv}, A = A ∩ (V × V ), ∪{(uv, x)|(u, x) ∈ A} ∪{(uv, x)|(v, x) ∈ A}, ∪{(x, uv)|(x, u) ∈ A} ∪{(x, uv)|(x, v) ∈ A}, ∪{(uv, uv)|(u, v) ∈ A} ∪{(uv, uv)|(v, u) ∈ A}. G → m G is written when there exists u and v such that m u,v (G) = G . G → *
m G is written when G can be obtained from G using zero, one, or several merges.

G → md G is written when G → m G or G → d G . G → *
md G is written when G can be obtained from G using zero, one, or several merges and/or deletes.

Merge and delete operations are actions on a graph structure, but they can be biologically construed. Some possible interpretations are:

-Deletion of a specie.

Specie's concentration is constant or its variation is insignificant compared to the absolute value. The specie can be ignored.

-Deletion of a reaction. Reaction's rate is zero or is deemed negligible compared to another.

-Merging of two species. Two species have proportional concentrations and can be confused.

-Merging of two reactions. Two reactions have proportional rates.

Morphisms

From the delete and the merge operation definition, subgraph epimorphism can be defined.

Definition 2.1.4 (Subgraph Epimorphism)

A subgraph epimorphism (SEPI) [1] from G to G is a function µ : V → V ∪ {⊥} such that :

∀(u, v) in µ -1 (V ), (u, v) ∈ A =⇒ (µ(u), µ(v)) ∈ A ,
µ surjective on V and A .

For reaction graphs, a SEPI verifies µ(S) ⊆ S and µ(R) ⊆ R .

Theorem 2.1.1 [1] The existence of a SEPI from G to G is equivalent to the existence of a finite sequence of merge and delete operations that yields a graph isomorphic to G when applied to G.

G SEP I G is written when there exists a SEPI from G to G .

Definition 2.1.5 (SEPI decision problem) Subgraph epimorphism problem is the decision problem:

Instance: Two Graphs G, G .

Question: G SEP I G ?

This chapter assumes familiarities with basics of theoretical computer science. Otherwise readers can refer to the work of Van Leeuwen [START_REF] Van Leeuwen | Handbook of theoretical computer science[END_REF].

Theorem 2.1.2 [1] Subgraph epimorphism problem is NP-complete.

Reaction graph

Partial order

Definition 2.1.6 (Total order) [START_REF] Fuchs | Partially ordered algebraic systems[END_REF] A relation ≤ is a total order on a set S if the following properties hold.

• Reflexivity: a ≤ a for all a ∈ S.

• Antisymmetry: a ≤ b and b ≤ a implies a = b.

• Transitivity: a ≤ b and b ≤ c implies a leqc.

• Comparability: for any a, b ∈ S, either a ≤ b or b ≤ a.

A total order on graphs can be defined.

G ≤ G implies | V |≤| G |.
Definition 2.1.7 (Partial order) [START_REF] Fuchs | Partially ordered algebraic systems[END_REF] A relation ≤ is a partial order on a set S if the following properties hold.

• Reflexivity: a ≤ a for all a ∈ S.

• Antisymmetry: a ≤ b and b ≤ a implies a = b.

• Transitivity: a ≤ b and b ≤ c implies a leqc.

A partial order can be defined on a set of graphs.

G ∩ md G = {H|G → * md H ∧ G → * md H}.
A maximal element of a set X is an element x ∈ X such that for all y ∈ X, y < x. x is the maximum of X if it is unique.

Definition 2.1.10 (Set of greatest lower bounds (glb

)) G ∩ md G is a set of → * md -maximal elements of G ∩ md G .

Definition 2.1.11 (Set of upper bounds) A set of upper bounds can be defined by G∪

md G = {H|H → * md G ∧ H → * md G }.
A minimal element of a set X is an element x ∈ X such that for all y ∈ X, x < y. x is the minimum of X if it is unique.

Definition 2.1.12 (Set of least upper bounds (lub There is a SEPI from G to G by deleting R -1 , merging R 1 with R 2 and merging E with ES.

)) G ∪ md G is a set of → * md -minimal elements of G ∪ md G . E R 1 S ES R 2 P R -1

Boolean satisfiability problem

In previous works [1], subgraph epimorphism problem has been encoded into a SAT problem.

Problem statement

A Boolean variable is a variable that can either be true (represented by a value 1) or false (represented by a value 0). Let X be a set of Boolean variables. A literal l is either x or its negation (denoted x), with x ∈ X. A clause is a disjunction (denoted ∨) of literals. A formula in conjunctive normal form (or CNF formula) is a conjunction (denoted ∧) of clauses. An assignment of a variable is to give a value 0 or 1 to the variable. An interpretation of a CNF formula is an assignment of its variables. A formula is said to be satisfiable if there exists an interpretation that makes the formula true. An interpretation satisfies a formula if it makes the formula true. A formula is said to be unsatisfiable if the formula is false for all possible interpretations.

Maximum satisfiability problem

Definition 2.3.1 (SAT problem) A Boolean satisfiability problem (SAT for short) is the problem of determining whether there exists an interpretation that satisfies a given CNF formula.

Problem resolution

A naive method to solve a SAT problem is to construct a truth table but the number of assignments is exponential in the number of variables of the input formula.

A systematic method with a binary search tree is more efficient. To prove the satisfiability of a CNF formula Φ, a literal l can be assigned, then, satisfiability of l ∨ Φ or ¬l ∨ Φ is recursively proved.

DPLL algorithm [START_REF] Davis | A computing procedure for quantification theory[END_REF][START_REF] Davis | A machine program for theorem-proving[END_REF] is based on this idea and on two other rules (unit and pure) that reduce the search space.

Modern SAT solvers are based Conflict-Driven Clause Learning (CDCL) algorithms [START_REF] Marques-Silva | Grasp: A search algorithm for propositional satisfiability[END_REF] which are inspired by DPLL algorithm, backtracking by conflict analysis with clause recording [START_REF] Silva | Grasp: A new search algorithm for satisfiability[END_REF], and Boolean constraint propagation using watched literals [START_REF] Moskewicz | Chaff: Engineering an efficient sat solver[END_REF].

Those tools are very efficient in practice while the problem is very hard in theory [START_REF] Vardi | Boolean satisfiability: theory and engineering[END_REF]. In worst case scenario, current SAT-solving algorithms must take exponential time to solve a problem but on the other hand solvers are applied successfully on large industrial instances.

Solver input

Input of a SAT solver is a text file in a simplified version of DIMACS format [START_REF] Challenge | Satisfiability: Suggested format[END_REF]. Besides comments, first line of a file should be p cnf nb_variables nb_clauses, where nb_variables is an upper bound on the largest index of a variable appearing in the file or the exact number depending on the implementation, and nb_clauses is the exact number of clauses.

Other lines of a file represent a CNF formula. Each variable corresponds to an integer, its negation is depicted by the inverse of corresponding integer. Each line is a disjunction of literals and stand for a clause. Each line end with a zero. The CNF formula is a conjunction of all clauses.

An example of input file:

c comments p cnf 5 3 1 -5 4 0 -1 5 3 4 0 -3 -4 0

Maximum satisfiability problem

Maximum satisfiability problem [START_REF] Karp | Reducibility among combinatorial problems[END_REF] is a generalisation of Boolean satisfiability problem where a CNF formula can be unsatisfiable.

Previous implementation

Subgraph epimorphism problem has been encoded into instances of SAT by S. Gay [1] and then resolved by a SAT solver. To describe a SAT instance, particular notations will be used.

Notations

Variables in bold font are SAT variables. cl(...) represents a clause. m(u) represents the image of u through morphism m. m u,u = 1 if and only if m(u) = u .

Elements of V ∪ ⊥ are put in a total order

v 0 = ⊥ < v 1 < ... < v n . m < u,u = 1 if and only if m(u) < u . m((u, v)) = (u , v ) represents the image of (u, v) through morphism m. m (u,v),(u ,v ) = 1 if and only if m((u, v)) = (u , v ). deleted (u,v) = 1 if and only if m((u, v)) = ⊥. type(u)
represents the type of a vertex u (specie or reaction).

Partial surjective function coding

First, partial surjective property of SEPI relation need to be encoded.

Left totality clauses are ensuring that each vertex of a starting graph has at least one image in a targeted graph or its image is bottom.

Right totality clauses are ensuring that each vertex of a targeted graph has at least one antecedent in a starting graph.

Functionality clauses are ensuring that each vertex of a starting graph has only one image through the morphism.

III Functionality

. i F f unc-i := ∀(u,u j )∈(S∪R)×(S ∪R ∪⊥) cl(m u,u j =⇒ m < u,u j+1 ), ii F f unc-ii := ∀(u,u j )∈(S∪R)×(S ∪R ∪⊥) cl(m < u,u j =⇒ m < u,u j+1 ), iii F f unc-iii := ∀(u,u j )∈(S∪R)×(S ∪R ∪⊥) cl(m < u,u j =⇒ ¬m u,u j+1 ). Number of clauses: 3 × (|S| + |R|) × ((|S | + |R |) + 1), number of new variables: (|S| + |R|) × (1 + |S | + |R |).

Subgraph epimorphism coding

Subgraph epimorphism property is implemented as follow.

Left totality on arcs clauses are ensuring that each edge of a starting graph has an image in a targeted graph or is deleted.

Right totality on arcs clauses are ensuring that each edge of a targeted graph has an antecedent in a starting graph.

Graph morphism clauses are linking SAT variables for edges and SAT variable for vertex by forcing endpoints of starting graph edges (respectively targeted graph edges) to have the corresponding image (respectively antecedent) through morphism m.

Subgraph morphism clauses are ensuring that if a vertex u is deleted, each edge with u as endpoint is also deleted. They also ensure that if an edge is deleted, at least one of its endpoint is deleted.

Redundant morphism propagation clauses are ensuring that the image of a starting graph edge is an edge in a targeted graph.

Bi-graph constraints clauses are ensuring that the image of a starting graph vertex through a morphism m is a vertex of same type in the targeted graph. Type being species or reaction. The formula is then defined as

Preliminaries

Clauses

III Graph morphism. i F morph-i := ∀(u,v)∈A,(u ,v )∈A cl(m (u,v),(u ,v ) =⇒ m u,u ), ii F morph-ii := ∀(u,v)∈A,(u ,v )∈A cl(m (u,v),(u ,v ) =⇒ m v,v ), iii F morph-iii := ∀(u,v)∈A,(u ,v )∈A cl((m u,u ∧ m v,v ) =⇒ m (u,v),(u ,v ) ). Number of clauses: 3 × |A| × |A |, number of new variables: 0. IV Subgraph morphism. i F sub-morph-i := ∀(u,v)∈A cl(deleted (u,v) ) =⇒ m u,⊥ ∨ m v,⊥ ), ii F sub-morph-ii := ∀(u,v)∈A cl(m u,⊥ =⇒ deleted (u,v) )), iii F sub-morph-iii := ∀(u,v)∈A cl(m v,⊥ =⇒ deleted (u,v) )).
F SEP I := F l-tot ∧F r-tot ∧F f unc-i ∧F f unc-ii ∧F f unc-iii ∧F l-tot-arcs ∧ F r-tot-arcs ∧ F morph-i ∧ F morph-ii ∧ F morph-iii ∧ F sub-morph-i ∧ F sub-morph-ii ∧ F sub-morph-iii ∧ F morph-prop ∧ F bi-graph .

Constraints on the number of deletions

This chapter presents a first approach to reduce the number of SEPI relations between two given graphs. When it is not specified, a SEPI relation, or SEPI, is a sequence of merge and delete operations that relate a graph to another. In a first section, goals are reminded. In a second section, a new approach is presented. An evaluation on BioModels database is shown in a third section. This chapter ends with a brief discussion and a conclusion.

Motivations

A major drawback of SEPI framework is the amount of SEPI found between reaction graphs. Some existing SEPIs do not have a biological interpretation and the set of SEPIs between two given graphs is often too big to be studied and to be understood. One of targeted goals is to reduce the number of distinct SEPIs between two given graphs.

Figure 3-1 is an example of classes of models from BioModels database. Each vertex represents a model and an arrow between two graphs means that there exists a SEPI relation between those two models. When there is a SEPI from a graph G 1 to a graph G 2 another one from G 2 to G 3 and a last one from G 1 to G 3 , arrow corresponding to the last SEPI is not drawn, in other words, transitive closure is not represented in this figure. Numbers on an arrows are numbers of different SEPI existing between two given graphs which is the number of different merge and delete operations that can relate a graph to another. When a number is equal to 200 it means that there is 200 or more different SEPI relations. This size is already too big to be studied and it can exists thousands of relations, all of them were not computed.

Figure 3-1 shows the existence of a lot of SEPI sets with a size greater than 200 and this is not optimal. S(G, G ) will represents a set of SEPI pairings from a graph G to a graph G .

Figure 3-2 and figure 3-3 are an example on small graphs of an explosion of S(G, G ) size. Figure 3-2 is a model created by M. Marhl et al. [START_REF] Marhl | Complex calcium oscillations and the role of mitochondria and cytosolic proteins[END_REF]. It represents complex calcium oscillations and mitochondria and cytosolic proteins roles in cells. Figure 3-3 is a model studied by J. M. Borghans et al. [START_REF] Borghans | Complex intracellular calcium oscillations a theoretical exploration of possible mechanisms[END_REF]. It also represents complex intracellular calcium oscillations. Both graphs are from Calcium Oscillations class.

On both figures, circle vertices represent species and square vertices represent reactions. An arrow is drawn from species to a reaction when they are reactants of the reaction. An arrow is drawn from reactions to a specie when the specie is a product of the reactions. This example also raises a problem of symmetries: parts of a reaction graph can have a same structure, the only difference being specie and reaction names it is also called an automorphism.

A lot of them can be observed in table 3-1. Symmetries in initial graph G or image graph G will have a multiplicative effect on the number of SEPIs: if µ is a SEPI from G to G , then for any σ automorphism of G and σ automorphism of

G , σ • µ • σ is also a SEPI from G to G .
In Michaelis-Menten example as well as in graphs of figure 3-2 and figure 3-3 described previously, a source of combinatorics is vertices that can indifferently be deleted or merged. It motivates a definition of a partial order on sets S(G, G ) that can be described more compactly by exhibiting only maximal or minimal elements. These extremal elements are sufficient to capture relevant information on relation between G and G , for example to get a set of possible images (respectively pre-images) of each species vertex in G (respectively G ).

However, there is no practical way to isolate these extremal elements. Isolating pairings that make extremal the number of deleted vertices (denoted min ⊥ and max ⊥) is another solution. Goals of minimizing or maximizing the number of deleted vertices are heuristics, they highlight SEPI pairings that have a chance to be relevant.

Implementation

To select solutions according to their number of vertex deletions, the SEPI decision problem can be redefined into an optimisation problem. First a few definitions are presented, then clauses are introduced.

Ca PrCyto

R 1 PrCyto R 2 R 7 Ca CytCyo R 4 R 6 R 8 R 5 Ca MMito Ca ErEndo R 3 Figure 3-2.: Calcium Oscillations: initial graph. X R 1 R 2 R 3 Z R 5 Y R 4 R 6 EC Figure 3-3.: Calcium Oscillations: targeted graph. CaPrCyto Ca_cytCyo PrCyto CaMMito CaErEndo R 1 R 2 R 3 R 4 R 5 R R 7 R 8 ⊥ 1 EC Z Z X Y R 5 R 6 R 1 ⊥ R 2 R R 4 R 4 1 2 X Z Z EC Y R 2 R 1 R 6 ⊥ R 5 R R 4 R 4 1 3 EC Z Z Y X R 5 R 6 R 3 R 4 R 4 R R 2 ⊥ 1 4 X Z Z Y EC R 2 R 1 R 3 R 4 R 4 R R 5 ⊥ 1 5 EC Z Z X Y R 5 R 6 R 1 ⊥ R 2 R ⊥ R 4 2 6 X Z Z EC Y R 2 R 1 R 6 ⊥ R 5 R ⊥ R 4 2 7 EC Z Z Y X R 5 R 6 R 3 R 4 ⊥ R R 2 ⊥ 2 8 X Z Z Y EC R 2 R 1 R 3 R 4 ⊥ R R 5 ⊥ 2 9 X Z ⊥ EC Y R 2 R 1 R 6 ⊥ R 5 R R 4 R 4 2 10 EC Z ⊥ X Y R 5 R 6 R 1 ⊥ R 2 R R 4 R 4 2 11 EC Z ⊥ Y X R 5 R 6 R 3 R 4 R 4 R R 2 ⊥ 2 12 X Z ⊥ Y EC R 2 R 1 R 3 R 4 R 4 R R 5 ⊥ 2 13 EC Z ⊥ X Y R 5 R 6 R 1 ⊥ R 2 R ⊥ R 4 3 14 X Z ⊥ EC Y R 2 R 1 R 6 ⊥ R 5 R ⊥ R 4 3 15 EC Z ⊥ Y X R 5 R 6 R 3 R 4 ⊥ R R 2 ⊥ 3 16 X Z ⊥ Y EC R 2 R 1 R 3 R 4 ⊥ R R 5 ⊥ 3

Partial orders

A partial order can be defined for the minimisation and maximisation of deletions. A graph homomorphism is a mapping between two graphs that respects their structure.

Definition 3.2.1 (Partial order ≺)

The partial order ≺ is defined on S(G, G ) by

µ 1 ≺ µ 2 if an d only if ∀v ∈ V, µ 2 (v) ∈ {µ 1 (v)
, ⊥} with µ 1 and µ 2 two graph morphisms (i.e µ 1 ≺ µ 2 if some vertices are deleted by µ 2 but merged by µ 1 ).

If µ 1 ≺ µ 2 then more vertices are deleted by the morphism µ 2 . SEPIs that minimize the number of deleted vertices are minimal for ≺. SEPIs that maximize the number of deleted vertices are maximal for ≺.

However, SEPIs that minimize the number of vertex deletions do not cover all minimal SEPIs for ≺. Conversely, SEPIs that maximize the number of vertex deletions do not cover all maximal SEPIs for ≺.

It would also be possible to minimise the number of species deleted while maximising the number of reaction deleted or vice versa.

Definition 3.2.2 (Partial order ≺ sr )

The partial order ≺ sr is defined on

S(G, G ) by µ 1 ≺ sr µ 2 if an d only if for all v ∈ S(G), µ 2 (v) ∈ {µ 1 (v), ⊥} and for all v ∈ R(G), µ 1 (v) ∈ {µ 2 (v), ⊥}
with µ 1 and µ 2 two graph morphisms (i.e µ 1 ≺ sr µ 2 if some species vertices are deleted by µ 2 while merged by µ 1 or some reaction vertices are deleted by µ 1 while merged by µ 1 ).

Only the first partial order will be implemented and tested. If partial order ≺ sr is interesting for biologists, its implementation would be straightforward.

Minimisation of the number of deletions

In order to minimise the number of deletions, a MAX-SAT solver is used.

|S| + |R| soft clauses are added representing for each vertex of the starting graph the vertex deletion's negation. All the previous clauses are transformed into hard clauses by adding a weight of |S| + |R| + 1 and the solver is used to maximise the number of satisfied clauses.

Maximisation of the number of deletions

Similarly, to maximise the number of deletions, for each vertex of the starting graph a soft clause representing its deletion is added.

No additional variables are needed.

Clauses

Evaluation

Lemma For a MAX-SAT encoding of the minimisation of the number of vertex deletions, with n variables and m = p + q clauses, p corresponding to the number of hard clauses and q to the number of soft clauses, the number of non satisfied clauses r is the minimum number of vertex deletion.

Proof: Each soft clause represents the negation of each vertex deletion. If a soft clause is not satisfied, the corresponding vertex is deleted. A MAX-SAT solver optimises the number of satisfied clauses m -r, thus r is the minimum number of vertex deletion.

Evaluation

To evaluate the implementation, first it was tested on handwritten models. There structure is well known and results can easily be understood. Then real models of the BioModels database [4] were used.

Handwritten models

Combinations of two Michaelis-Menten

Reduction of Michaelis-Menten reaction, presented as the first example is the most simple and classical example of model reduction. Michaelis-Menten reactions are also present in a lot of complex reactions.

Figure 1-3 shows an expanded form of the reaction, where a complex enzyme-substrate ES is present, as well as a reverse reaction R -1 from the complex to initial species substrate and enzyme. Figure 1-4 shows a reduced form of the reaction. Several intermediary forms could be considered, coming examples only take into account intermediary forms with a complex ES and reactions R 1 and R 2 but without reverse reaction R -1 .

This first example is important because Michaelis-Menten-like reactions can be found in more complex models with different forms of reduction. Hence, a basic feature any variant of SEPI framework should ensure is to find Michaelis-Menten-like reductions when they exists. Moreover, there is already three possible SEPIs for the most basic example of Michaelis-Menten reduction.

If two graphs G and G contains Michaelis-Menten patterns, the number of SEPIs between those two graphs is exponential in the number of Michaelis-Menten patterns.

Other small examples are combinations of two Michaelis-Menten reactions. For example, with species aS, aE, aP and bS, bE, bP as substrates, enzymes and products of a first motif a and a second motif b. Reactions are combined if some of their species are shared, for example if aP = bS (in this case it will be represented by aP bS). Combinations are denoted by shared species, for example EE -P S represents graph of figure 34where two enzymes are shared and product of the first reaction plays a role of substrate of the second. In this graph, both Michaelis-Menten reactions are in reduced form. There exist 23 combinations of two Michaelis-Menten patterns: 6 combinations with only one shared specie, 12 combinations with two shared species and 5 combinations with three. Some combinations can seem artificial but some of them are found in biological models, for example in Mitogen-Activated Protein Kinases (MAPK) cascade models: P E when the enzyme is the product of a previous reaction, EE -P S for a two-steps phosphorylation with same enzyme, P S -SP for a Michaelis-Menten reaction that is reversible with another enzyme.

MAPK cascades

More complex combinations of Michaelis-Menten reactions can be used to model MAPK cascades. Figure 3-5 presents a three-levels model corresponding to reduced form of models described by C.Y. Huang et al. [START_REF] Huang | Ultrasensitivity in the mitogen-activated protein kinase cascade[END_REF] and by A. Levchenko [START_REF] Levchenko | Scaffold proteins may biphasically affect the levels of mitogen-activated protein kinase signaling and reduce its threshold properties[END_REF].

Graph if figure 3-5 is a combination of 10 Michaelis-Menten reactions. Between expanded and reduced forms, there are more than 2 10 SEPIs: three SEPIs for each Michaelis-Menten reaction, but when an enzyme is shared it cannot always be deleted and a lot of Michaelis-Menten reactions share an enzyme in MAPK cascades.

Other handwritten models can be considered with one-level or two-level cascades.

Pairs with SEPIs

BioModels

The two cited models from Huang-Ferrel and Levchenko correspond to models 9 and 11 of the publicly-accessible repository BioModels [4]. Previous works on SEPIs [1] and [START_REF] Gay | A graphical method for reducing and relating models in systems biology[END_REF] used a set of 42 models from BioModels for their evaluations. These models have been manually clustered by annotations according to the kind of biological process they represent (MAPK cascade, calcium oscillations, circadian clock, cell cycle). Gathered data on reaction graph sizes are presented in table 3-4. Relevancy of SEPI theory was justified because it permits to find the same clustering a priori: more SEPIs are expected intra-class than inter-class.

This set of models is interesting for an evaluation because different kind of observations can be made. SEPI implementation is too relaxed if a lot of reductions are found between models from different classes but is too strict if it reduces the amount of SEPI intra class. To minimise or maximise the number of vertex deletions, a MAX-SAT solver was used. To select the best one, a benchmark on three different MAX-SAT solvers was made. Those three solvers: RC2 [START_REF] Ignatiev | Rc2: a python-based maxsat solver[END_REF], QMaxSat [START_REF] Zha | Qmaxsat in maxsat evaluation 2018[END_REF] and OpenWBO [START_REF] Joshi | Open-wbo-inc in maxsat evaluation 2018[END_REF] were selected among many other solvers according to results showed in a MAX-SAT competition [START_REF] Bacchus | Maxsat evaluation 2018[END_REF]. Results of each solver were compared on four classes of BioModels.

Number of vertices

Results can be seen in table 3-7. Results are expressed in seconds, timeouts are declared after running the solver for 20 minutes without having results. Additionally table 3-8 presents the time in milliseconds to found first, second and last SEPIs between two given models.

According to these results, the MAX-SAT solver QMaxSat was kept for the implementation in Biocham but after several tests RC2 MAX-SAT solver was better suited for the framework.

After the selection of the MAX-SAT solver, more tests were made on the BioModels database. Results are presented in table 3-5. Additionally, more precise results can be observed in figure A-1 to figure A-12 of the annexes.

On real models, compared to hand-written models, a few observations can be made:

• Globally, the number of timeouts is highly increasing, in particular for maximisation of the number of deletions. Table 3-6 summarises timeouts of each classes. It makes the framework inoperable on big models: there is too few results to interpret them correctly.

• It can still be noted that maximisation of the vertex deletion number filters more reductions.

• Models of different classes behave in distinct ways. Results are better on models representing calcium oscillations than on models depicting circadian clocks.

Drawbacks

However, SEPIs that minimize (respectively maximize) the number of vertex deletion unfortunately do not cover all SEPIs minimal (respectively maximal) for ≺. Next figures are an example of SEPIs without biological interpretation that minimises the number of deletion. This pairing is minimal for ≺ because if same images for non-deleted vertices are kept, these reverse reactions cannot be merged with other ones. For example, if R -1 has to be merged, the only option to preserve an outgoing arc from K is to merge R -1 with R 3 , with image R 3 . But it is not possible because an arc from R -1 to E 1 cannot have an arc from R 3 to E 1 as image.

Constraints on the number of deletions

E1 R -1 R -2 E1-K E1-pK R 1 R 1 R 2 R 2 K pK ppK R 3 R 3 R 4 R 4 E2-pK E2-ppK R -3 R -4 E2 Figure 3-6.: Expected reduction of a one-step kinase phosphorylation. E1 R 1 R 2 K pK ppK R 3 R 4 E2 Figure
However, colours in figure 345678show a pairing where no vertex is deleted. Hence, the minimal number of vertices is 0, which is not achieved with the previous pairing. This example is problematic because the pairing exhibited by figure 3-8 is very unexpected: entry K and output ppK play a role of enzymes, while all enzymes are merged on pK. This unexpected pairing is also valid with accurate merge restriction.

E1

R -1 R -2 E1-K E1-pK R 1 R 1 R 2 R 2 K pK ppK R 3 R 3 R 4 R 4 E2-pK E2-ppK R -3 R -4

E2

Figure 345678.: Unexpected reduction that minimizes bottom.

Conclusion

Adding constraints on the number of vertex deletion gives mixed results. It is reducing a bit the size of the set of SEPIs between two given graphs but the number of timeouts in too big to use this implementation on big models.

After filtration with maximisation or minimisation of the deletion number, SEPIs were more relevant but there were still unwanted SEPIs as the one presented in previous section.

This improvement on its own is not enough but more restrictions on the merge operation can be added. This will be the aim of next sections.

Accurate merge restriction

This chapter will present a restriction of the merge operation, which aims to filter SEPIs without biological interpretation between models, in particular, SEPIs that allow to merge vertices far apart in the reaction graph.

In a first section, motivations will be presented. A formal definition of the restriction will be explained in a second section. Then the implementation will be described in a third section and in the last two sections, an evaluation will be performed and conclusions will be drawn.

Motivations

A major drawback of the SEPI framework is the amount of SEPI found between models. Some existing SEPIs do not have a biological interpretation and the set of SEPIs between two given graphs is often too big to be studied and to be understood according to biologists using the framework. The SEPI framework is not restricting enough.

In the previous chapter extremalisation of the vertex deletion number has been implemented in order to reduce the size of the set of SEPIs found between two given models. But this implementation only reduces the size of a set of SEPIs, it doesn't delete SEPI relations between models.

The SEPI relation is a logical relation defined by the deletion or fusion of vertices. This logical relation between two given graphs doesn't always have a meaning in biology. Therefore, some restrictions to the SEPI relation need to be implemented in order to give more sense. The encoding of previous chapters was still interesting for biologists. It gives insights of structural similarities.

By observing SEPIs drawn between models of opposite classes, it can be seen that sometimes species far apart in the reaction graph are merged together. This is not appropriate because it doesn't have a biological interpretation. It makes more sense to merge species or reactions that are close together. Usually a reduction correspond to the contraction of a chain of reactions or to the reduction of a pattern.

A restriction on the merge operation is thus implemented to allow to merge vertex only when they are close together or when they are close together because their neighbours were also merged together.

Filtering SEPIs is expected to have two main results:

• Reduce the number of SEPIs inter class.

A high number of SEPIs inter class is not wanted because these SEPIs doesn't always have a biological interpretation. In figure 4 • Reduce the combinatoric of SEPIs between two given models.

It is not good to have a high combinatoric of SEPIs between two given models because

Biocham is aimed for biologists that study this SEPIs. Thus having more than 200 SEPIs to study is not practicable. In figure 4-1 each vertex represent a model from BioModels. An edge is drawn from a model to another when there exists a SEPI relation between these two models. Yellow vertices represent models from the class "MAPK", green vertices represent models from the class "Circadian Clock", blue vertices represent models from the class "Cell Cycle" and purple vertices represent models from the class "Calcium Oscillations". And model 73 to model 74 is considered as another pair. A SEPI relation will be searched from model 73 to model 74 but there is no relation. The column "SEPI" represents the amount of SEPI relations found among all pairs of models. The column "No SEPI" counts the absence of SEPI relations between pairs of models. A timeout is declared after running the SAT solver for 20 minutes, it is represented by the column "Timeout".

Pairs

Definitions

A restriction on the merge operation is implemented in this chapter to allow to merge vertices only when they are close to each other or when they are close because their neighbours were also merged together.

For this purpose, the two-neighbours relation will be defined and a notion of good-path will also be explained in this section and implemented in the following section.

As a reminder, a reaction graph G is a triple G = (S ∪ R, A), where S is a set of specie nodes, R is a set of reaction nodes, and A ⊆ S × R ∪ R × S is a set of arcs that describes how species interact through reactions.

To simplify the notation two other sets will be introduced, V and

E. V = S ∪ R is a set of vertices. E ⊆ S × R ∪ R × S is a set of un-oriented arcs: E = {(u, v)|(u, v) ∈ A ∨ (v, u) ∈ A}.
Let G and G be two graphs defined by G = (V, E) and G = (V , E ).

Because reaction graphs are bipartite graphs, graphs whose vertices can be divided into two disjoint and independent sets, only vertices of the same type, specie or reaction, can be merged together. A definition of neighbours cannot be used for a merge restriction because neighbours vertices will always be of different types. A two-neighbours definition is then introduced to explain this neighbour of neighbour notion. It is important to associate a length to each path*, defined by the length of smaller intermediary paths* C i , in order to avoid cyclic situations as described in figure 45. It is the reason why a characterization of paths* of length k is needed for the implementation. In this disconnected graph, a graph where there exist two nodes such that no path has those nodes as endpoints, A and B have the same image, and r1 and r2 as well. The argumentation that there is a path* between A and B with C 1 a path* between r1 and r2, itself defined thanks to A and B, would be infinitely recursive, no base case is defined. Hence, no length could be defined. It shows why a well-defined length is required for paths*.

A new definition of the SEPI relation can then be defined. More formally,

G SEP I * G ⇔ G → * m * d G .
Proof: (⇒) The only if direction will be proved by induction on the number of vertices.

• base case:

|V | = |V |
No vertex is merged or deleted. A SEPI* gives an isomorphism between G and G , hence the empty sequence works.

• • induction step: n = 0

∼ µ : G -→ G x -→ µ(a) = µ(b) if x = ab µ(x) if not ∼ µ is a SEPI*,
Induction hypothesis: the implication is true for any sequence of length smaller or equal to n.

A sequence of length n + 1 can be decomposed in a sequence of length n between G and an intermediary graph G and a single operation o between G and a graph isomorphic to G . The induction hypothesis gives a SEPI* µ from G to G .

-case (i -case (ii Observation: SEPI* is not a well-quasi-order.

): o = d u With µ -1 (u) the antecedent of u through the morphism µ, ∼ µ is defined by: ∼ µ : G -→ G x -→ ⊥ if x ∈ µ -1 (u) µ(x) if not By definition, ∼ µ(G) is isomorphic to G . ∀a, b ∈ V , µ(a) = µ(b) ⇒ ∼ µ(a) =
): o = m * u,v ∼ µ is defined by: ∼ µ : G -→ G x -→ uv if x ∈ {µ -1 (u), µ -1 (v)} µ(x) if not ∼ µ(G) is isomorphic to G . ∀a, b ∈ V \{µ -1 (u), µ -1 (v)}, µ(a) = µ(b) ⇒ ∼ µ(a) =
G G G 1 G 2 G 1 G 3 G 2 G G m * u,v m * uv,w m * v,w m * u,vw d v d u d u d v m * u,v d uv
Proof: SEPI* relation is more restrictive than SEPI. Hence, the same example of infinite antichain given for SEPI by S. Gay et al. [START_REF] Gay | On the subgraph epimorphism problem[END_REF] works.

Lemma:

The SEPI*-decision problem is NP-complete.

This lemma was proved by Eva Philippe.

Proof: The proof of NP-completeness found by C. Solnon [1] only needs a very slight modification to apply to SEPI*: adding a vertex w that makes the merged vertices be two-neighbours.

The idea is to reduce from exactly-k-set-covering problem [START_REF] Karp | Reducibility among combinatorial problems[END_REF], which is NP-complete, to SEPI*. This exactly-k-set-covering problem is defined by:

Instance. A set E,
A = {(u i , r i ) | 1 ≤ i ≤ m} ∪ {(r i , e j ) | e j ∈ u i } ∪ {(r i , w) | 1 ≤ i ≤ m} ∪ {(w, r i ) | 1 ≤ i ≤ m}. G = (S , R , A
) is defined with S = {s 1 , . . . , s m } ∪ {e 1 , . . . , e n }, R = {r }, and

A = {(s i , r ) | 1 ≤ i ≤ k} ∪ {(r , e j ) | 1 ≤ j ≤ n}.
It will be proved that the exactly-k-set-covering problem has a positive answer if and only if there is a reaction graph SEPI* from G to G .

⇒ Hypothesis:

U * = {u * 1 , . . . , u * k } is a solution of (E, U, k). If U \ U * = {u * k+1 , . . . , u * m }.
Then µ defined by u j = u * i → s i , r i → r , e i → e i is a reaction graph SEPI* from G to G (merged vertices all belong to R and have w as a common neighbour).

⇐ Hypothesis: there is a reaction graph SEPI* µ from G to G .

With µ -1 (r ) = {r i 1 , . . . , r i k } and U * = {u i 1 , . . . , u i k }, U * is a solution of the covering problem. First, U * covers E. Indeed, let e ∈ E. Then the arc (r , µ(e)) is covered in G , by some arc (r i , e 1 ) in G. r i ∈ µ -1 (r ), but is e 1 always the vertex e ? Yes: since there are the same number of s-vertices with no exiting arcs in G and G , µ induces a bijection between s-vertices with no exiting arcs of G and G , so e 1 = e. Which proves that e is covered by the subset corresponding to r i . Next, U * has k elements. Vertex types force arcs (u i , r i ) to be the only ones that can cover the (s i , r ). So the preimages of the s i are some u i , and there must be exactly k such u i because of the bijection on s-vertices with no entering arcs. U * is exactly those u i . The coding being polynomial, this concludes the proof of reduction from exactly-k-set-covering to the reaction graph SEPI* decision problem, and the proof of NP-completeness. The graph 4-8(a) may be transformed into graph 4-8(b) by deleting the reaction vertex associated with {a, c, d} and by merging the three other reaction vertices, respectively associated with {a, b, d}, {c, f } and {e, f }. These vertices are two-neighbours thanks to vertex w, which is deleted. This corresponds to the solution of the set covering problem such that the three selected subsets are {a, b, d}, {c, f } and {e, f }.

However, graph 4-8(a) cannot be transformed into graph 4-8(c), as the set covering problem instance has no solution for k = 2.

Implementation

Despite the NP-completeness, the SEPI* problem can be encoded in SAT instances like SEPI. To implement a restriction of the merge operation, constraints are added to the SEPI framework. This new framework is called the SEPI* framework.

For this implementation, new variables need to be introduced and additional clauses need to be written.

New variables

The same notations as the one of section 2.5.1 are used in this section. As a reminder, m is a graph morphism from G to G , m a,y = 1 if and only if m(a) = y.

The following variables will be introduced, in addition to the variables of the SEPI framework.

Clauses

The following clauses will be introduced, in addition to the clauses of the SEPI framework.

Inductive merge definition clauses define the new variable m a,b,y . This variable will be required to remove the conjunction from the clause Good-path of length k.

The definition of a good-path will be made recursively. Good-path of length zero clauses define the base case.

To define the inductive step the variable p a,b,k 1 ,c,d,k 2 need to be introduced to remove a conjunction of literals. Clauses inductive good-path definition define this new variable.

Good

cl(¬p a,b,k 1 ,c,d,k 2 ∨ p a,b,k 1 ), ii F inductive-path-ii := ∀a,b,c,d∈V ∀k 1 ,k 2 ∈ 0,(min(|S|,|R|)-1) cl(¬p a,b,k 1 ,c,d,k 2 ∨ p c,d,k 2 ), iii F inductive-path-iii := ∀a,b,c,d∈V ∀k 1 ,k 2 ∈ 0,(min(|S|,|R|)-1) cl(¬p a,b,k 1 ∨ ¬p c,d,k 2 ∨ p a,b,k 1 ,c,d,k 2 ). IV Good-path of length k: i F path-k-i := ∀a,b∈V ∀k∈ 1,min(|S|,|R|) cl(¬p a,b,k ∨ ( y∈V ∪{⊥} m a,b,y )), ii F path-k-ii := ∀a,b∈V ∀k∈ 1,min(|S|,|R|) cl(¬p a,b,k ∨ c∈V c =a u∈V (a,u)∈E v∈V (c,v)∈E k ∈N k <k p u,v,k ,c,b,k-k -1 ).
V Path*:

F path * := ∀a,b∈V ∀y∈V ∪{⊥} cl(¬m a,b,y ∨ ( k∈ 0,min(|S|,|R|) p a,b,k )).
The formula is then defined as

F SEP I * := F SEP I ∧ F merge-i ∧ F merge-ii ∧ F merge-iii ∧ F path-0-i ∧ F path-0-ii ∧ F inductive-path-i ∧ F inductive-path-ii ∧ F inductive-path-iii ∧ F path-k-i ∧ F path-k-ii ∧ F path * .

Optimisations

To introduce fewer variables and to write fewer clauses, some optimisations can be made. An order on graph G vertices is introduced. Even with the optimisation, the number of clauses to high to compute the SEPI* on big models. To reduce the number of clauses, the length of the good-path can be restricted to a fixed number. It has been implemented by replacing min(|S|, |R|) by a smaller number in the clauses. min(|S|, |R|) is the estimation of the maximal length of the good-path.

Evaluation

This section discusses an evaluation of this new framework and problems it raised.

The problem of reproducibility

Biocham version and subgraph epimorphism implementation changed since the last results of Steven Gay's thesis [1]. The thesis had good results with the number of SEPI connections inter classes but it was not reproducible with the new implementation. Even by using the same old models from 2015 the best current result is an inter class connection of 29% instead of the 9% expected [1].

Tests were performed on the same models, the ones from 2015, with the same timeout of 20 minutes and using the same SAT solver Glucose. Kinases. The number in brackets represents the number of model pairs in the class (or couple of classes). Again, a timeout if declared when the SAT solver didn't give a result after 20 minutes.

It can be observed in table 4-2 that the previous amount of timeouts was much lower. The number of SEPI relations between pairs of models was also higher. It is a bit problematic to not find the same results in 2019.

In table 4-3 a really high number of SEPIs inter class can be noted. It is more than three times more compared to previous results.

Benchmark for further evaluation

Because it is difficult to reproduce results from 2015, benchmarks were computed again on new versions of models from BioModels. These models were curated with the latest version of Biocham and new tests are performed with the latest version of the SEPI framework too.

Accurate merge restriction evaluation

The implementation of this new framework needs to define O(n 6 ) variables p a,b,k 1 ,c,d,k 2 (with n = max(|V |, |V |)). All this variables are required to replace p a,b,k 1 ∧ p c,d,k 2 in the clause characterizing path* of length k. Therefore, for big models even writing the clauses before trying to resolve the SAT problem takes too much time. This implementation is not usable in practice.

Conclusion

The definition of good-path is honourable in a logical point of view. It could achieve the goal of filtering unwanted SEPIs and it brings pleasant properties. But the implementation is very heavy. This implementation was the most optimum found and even with some optimisations and by reducing the maximal length of the good-path between each vertex it is not usable in practice.

Another implementation needs to be found. A lighter one that could be used on big models of the BioModels database. This will be the goal of the next chapter.

Strict two-neighbours restriction

This chapter presents a second merge restriction implemented, called strict two-neighbours merge restriction. A first part explains why a new implementation is necessary, a second part in this chapter presents this implementation. Finally, this chapter ends with an evaluation of this new implementation and a conclusion.

Motivations

Previous implementation was too complex and too slow to be computed on big models of BioModels. For example, the biggest model of MAPK has 334 vertices and 744 edges. This implementation needs to define O(n 6 ) variables p a,b,k 1 ,c,d,k 2 , which are required to replace p a,b,k 1 ∧ p c,d,k 2 in the clause characterizing path* of length k. Hence, for big models even writing clauses requires an unreasonable amount of time. This implementation is not usable in practice for big models.

p a,b,k 1 ,c,d,k 2 ( |S|×(|S|-1) 2 + |S|) × ( |R|×(|R|-1) 2 + |R|) × 2 × (min(|S|, |R|) + 1) 2 Order of magnitude m a,b,y O(n 3 ) p a,b,k O(n 3 ) p a,b,k 1 ,c,d,k 2 O(n 6 )
Unfortunately, a better optimised implementation for the accurate good-path merge restriction was not found. It doesn't seem to be possible to find local properties implied by the existence of a good-path between two vertices: there can be arbitrary size-neighbourhoods of vertices with the same image that have nothing in common. A local definition is regardless implemented. The good-path merge restriction cannot be used on models of BioModels anyway. This implementation is more restrictive and removes SEPIs that had a biological interpretation but was still considered as a good compromise. Strict two-neighbours implementation is equivalent to maxpath = 1 in the previous implementation. Another major drawback of this stronger restriction is the loss of SEPI's transitivity property.

Implementation

Strict two-neighbours merge restriction is a much stronger restriction. With the previous definition of good-path it corresponds to a maximal length of 1. It is a local constraint, two vertices of the initial graph can have the same image in the targeted graph if and only if they are two-neighbours in the initial graph. This restriction doesn't need recursive nor dynamic definitions.

Michaelis-Menten reduction is compatible with this restriction. Then the formula is then defined as F SEP I-strict := F SEP I ∧ F merge .

Clauses

Pseudocode

Before writing clauses of subsection above, the relation two-neighbours between vertices need to be identified.

As a reminder, the two-neighbours definition is the following: a and b are two-neighbours if and only if ∃r ∈ V such that ((a, r) ∈ E) ∧ ((b, r) ∈ E), with E the set of non oriented edges of the initial graph and V the set of vertices of the initial graph.

Pairs

BioModels

An evaluation on models of BioModels was also made. Figure 5-3 shows SEPI relations intra class with the initial SEPI framework. Figure 5-4 shows SEPI relations intra class after the implementation of the strict merge restriction. Numbers on edges represent the amount of SEPIs found between the two models.

When a SEPI relation is found between two models, it can be observed in the figures that the size of the set of SEPIs is decreasing. But a diminution of SEPI relations can also be observed. Table 5-5 will give more precise results.

Figure 5-5 shows SEPI relations inter class with the initial SEPI framework. Figure 56shows SEPI relations inter class after the implementation of the strict merge restriction.

A significant diminution of unwanted SEPIs inter class can be observed in the figures. Quantitative results will be displayed in table 5-6.

Table 5-5 and table 5-6 present more precise results. Each block of two lines present results of a specific class or for the comparison of two specific classes. For each class, or pair of classes, the number in brackets is the number of tested pairs of models. For each class, or pair of classes, two lines of results are presented, the line "Without" presents the number of SEPIs with the initial SEPI framework, the line "Strict" presents the number of SEPIs with the strict merge restriction. For each line, the number of SEPIs, the number of no SEPI relations and the number of Timeouts are given. For SEPIs intra class, the size of the set of SEPIs is also relevant, the number of set with a size bellow 200 is also specified.

Several observations can be made from the two tables:

-The strict merge restriction eliminates 64.49% of unwanted SEPIs inter class.

-The strict merge restriction eliminates 80.23% of timeouts inter class.

-The strict merge restriction increases by a factor of 1.33 the number of set of SEPIs bellow a size of 200 among SEPIs intra class.

-The strict merge restriction eliminates 78.05% of timeouts intra class.

-The strict merge restriction makes 42.05% of SEPI intra class disappear.

Those observations are the targeted results except for the last one. But overall, considering the low complexity of this implementation, those results are encouraging.

SEPI

Execution time

Execution time of the new framework is also an important criterion for biologists.

Table 5-7 shows the time needed to write clauses and to solve specific instances.

To make a benchmark, two models of each class were selected. For the benchmark to be representative, small and big models were selected. Handmade models of MAPK cascade were also used as comparison, they are represented in the table by "mapk1", "mapk2" and "mapk3".

The first column of the table (called Model) represents each problem. For each problem, a SEPI relation is searched between two models of the same class, the corresponding class is specified in between brackets.

Is the first half of the table (called No Restriction), tests were made without any merge restriction. In the second half of the table (called Merge Restriction), the strict two-neighbours merge restriction was used.

The column time specifies the time needed to write the clauses. In column time of the No Restriction's half, it corresponds to the time for all clauses of the SEPI framework. In column time of Merge Restriction's half, it's only the additional time needed to write the clauses corresponding to the merge restriction. All times of the table are in milliseconds.

Then for each half (no restriction or merge restriction) and for each problem (each pair of models) a SEPI relation is searched without any extremalisation of the number of vertex deletions (column Normal) then with minimisation of the vertex deletion number (column Min ⊥) and with maximisation of the vertex deletion number (column Max ⊥).

When a SEPI relation is searched between two models the column nb corresponds to the number of different SEPI relations found. Because there can be thousands of SEPI relation, only the first 200 ones were computed. When the set of SEPI relations exceed 200 it is not humanly possible to look for an interpretation of each of them so the goal is to keep this number below 200. Thus 200 in column nb represents 200 or more SEPI relations.

For each pair of models and for each problem, columns 1st, 2nd and last represent respectively the time needed to find the first the second and the last SEPI. A t represents a timeout. Timeouts are declared after running the SAT solver for 20 minutes without result. A number in brackets represents the time needed to establish that the problem was not satisfiable (no SEPI). Arepresents the absence of result (no time for second and last SEPI when there exists only one SEPI and no number of SEPIs when there is a timeout).

A few observations can be made about the results of this table:

• It takes a relatively low amount of time to write the clauses corresponding to the twoneighbours merge restriction.

• It is more complex for a MAX-SAT solver to maximise the number of deletion than to minimise it.

• The two-neighbours merge restriction makes all timeouts disappear.

• Maximisation of the number of deletions in addition to the two-neighbours merge restriction is very efficient, especially for handwritten models of MAPK.

• The two-neighbours merge restriction makes three SEPI relations disappear (for pair of models 007-168, 021-170 and 011-026).

Conclusion

The strict merge restriction is better than the accurate merge restriction in terms of computational complexity. Results are obtained faster, the strict merge restriction reduces significantly the number of timeouts during tests of models of BioModels.

This restriction is stronger and makes some wanted SEPI relations disappear but it also reduces the amount of unwanted SEPIs inter class. This is a good evolution of the SEPI framework and it is now implemented and accessible for users of Biocham. But results given my the SEPI framework with strict merge restriction need to be taken carefully, a SEPI of a complex form can still exist between big models and not been shown by the updated version of the framework.

Pattern reduction

Motivations

In previous sections, two main improvements were made on the SEPI framework:

. Extremalisation of the vertex deletion number has been implemented (section 3) in order to reduce the size of the set of SEPIs found between two given graphs.

. Restrictions on the merge definition were implemented (section 5) in order to reduce the number of non biologically explainable SEPIs.

These improvements showed encouraging results, for example:

-76.32% of SEPI sets were above a size of 200 units in the class Calcium Oscillations before minimisation of the number of vertex deletions and only 42,11% after.

-27.02% SEPIs were found between models of different classes before implementing the merge restriction, only 9.60% after.

-11.02% of model pairs are too complex for the SEPI framework to compute them in a reasonable time. Only 2.42% are too complex for the SEPI framework with the merge restriction.

But previous approaches still allow some unexpected pairings, as shown in figure 6-1 and figure 6-2.

A SEPI with two-neighbours merge restriction and minimisation of the number of deletions is searched from the initial graph to the image graph. Figure 6-1 is the initial graph: it represents a simple combination of two Michaelis-Menten reactions with complexes substrate enzyme and with reverse reactions R -1 and R -3 . The first product is the substrate of the second reaction and the enzyme is shared. This combination is quite natural because it can be observed in models of MAPK.

E aS R 1 aC R 2 aPbS R 3 bC R 4 bP R -1 R -3
Figure 6-2 is the image graph: it represents a chain of two Michaelis-Menten, both in reduced form, where the first product acts as an enzyme on the same substrate to produce a second product.

A SEPI µ with strict two-neighbours merge restriction and deletion number minimisation is found between the two graphs. µ is represented by the colours in the figures. Vertices of one colour in figure 6-1 are mapped to the vertex of the same colour in figure 6-2:

µ(aP bS) = aE, µ(aC) = S, µ(R 2 ) = µ(R 3 ) = R 0
, mapping those two reactions on the same vertex is possible because aP bS is a common neighbour, µ(bC) = µ(bP ) = aP bE, mapping those two reactions on the same vertex is possible because R 4 is a common neighbour,

µ(R -1 ) = µ(R 4 ) = R 0 , mapping those two reactions on the same vertex is possible because E is a common neighbour, µ(aS) = bP , µ(R 1 ) = µ(R -3 ) = µ(E) = ⊥.
This SEPI does not have a biological interpretation but was not filtered by the strict twoneighbours merge restriction because the shared enzyme makes all reactions two-neighbours.

This pairing is not relevant because, for example, the entry aS of the first expended Michaelis-Menten pattern becomes the output bP of the second reduced Michaelis-Menten pattern and the complex aC of the first expended Michaelis-Menten pattern becomes the entry S of the second reduced Michaelis-Menten pattern.

Some patterns like the Michaelis-Menten pattern were precisely studied. The expected reduction of those patterns are well known. Thus, another way of filtering more SEPIs without biological interpretation is to reduced known expanded patterns.

Reducing patterns before using the SEPI framework between two given graphs will also reduce the size of the set of SEPIs found between those two graphs. Sepi Furthermore, reducing patterns before using the SEPI framework will reduce the computation time of the subgraph isomorphism research as it will reduce the number of vertices. And reduced graphs are also easier to study for biologists.

E S ES P R 1 R -1 R 2 1 E S E P R ⊥ R 2 E S ⊥ P R ⊥ R 3 ⊥ S E P R ⊥ R
The graph rewriting strategy is the following: This pre-processing is expected to reduce both the number of pairs that present a SEPI by eliminating SEPIs without biological interpretation, and to reduce the number of SEPIs for pairs of graphs that have at least one SEPI with interpretation.

Reducing pattern is a variant of the subgraph isomorphism (SISO) problem. But searching SISOs is not sufficient in this case because there is more constraints on a subset of vertices. For example in the Michaelis-Menten pattern, the complex ES in figure 6-3 cannot interfere in any other reactions and reactions R 1 , R -1 and R 2 cannot interfere with other species. This section will present the graph rewriting strategy, also named constrained subgraph isomorphism (cSISO) problem. In part 6.2, a broader framework of cSISO will be presented. This framework allow the search of any pattern with constraints on a subset of vertices. Although part 6.4 will show that the implementation is only made for patterns defined in part 6.3.

Definitions

The Michaelis-Menten pattern will be the main example of this section because it is a well studied reaction but other patterns with known reductions can also be used. The graph in figure 6-3 is an example of reaction graph with expended Michaelis-Menten pattern. The graph in figure 6-4 is the equivalent of the previous reaction graph with reduction of the Michaelis-Menten pattern.

Only searching subgraph isomorphisms is not sufficient in this pattern reduction because there are more constraints :

• the complex ES cannot interfere in other reactions, In this section, theory is given for general graphs (not necessarily bipartite), but it extends immediately to reaction graphs.

• reactions R 1 , R -1 and R 2 cannot use other species. A R 3 E R 1 S ES R 2 P R 4 B R -1
Let G 1 = (V 1 , E 1 ) be a graph with V 1 a set of vertices and E 1 a set of arcs. Let F 1 be a subset of V 1 corresponding to fixed vertices. For example, in the reaction graph with expended Michaelis-Menten pattern of figure 6-3,

F 1 = {ES, R 1 , R -1 , R 2 }.
Let G 2 = (V 2 , E 2 ) be a second graph defined in the same way, with F 2 a subset of V 2 .

A pattern defined by a graph G 1 is searched in a second graph G 2 with constraints defined by F 1 .

If F 1 = ∅, the constrained subgraph isomorphism problem corresponds to searching an induced subgraph isomorphism.

The definition is given with a possible fixed set F 2 to get a partial order on the space of pairs (G, F ) where G is a graph and F a subset of vertices of G. In the case F 1 = F 2 = ∅, the cSISO becomes the induced subgraph isomorphism problem. The induced SISO problem is known to be NP-complete as the search of a k-clique is NP-complete and can be reduced to SAT problem [START_REF] Karp | Reducibility among combinatorial problems[END_REF].

In the case of a pattern reduction in reaction graphs, (G 1 , F 1 ) is fixed and

F 2 = V 2 .
The time and space complexities will therefore be linear in |V (G 2 )|, event to find all the cSISOs from

(G 1 , F 1 ) to (G 2 , F 2 ).
The reduction of patterns is not always commutative, it depends of the pattern.

This definition of cSISO is general but pattern reduction cannot be easily implemented with every patterns. A few characteristics need to be respected:

• for the commutative property, patterns cannot overlap in their extended forms,

• vertices of F 1 form a connected component in G 1 .

For one cSISO µ : (G 1 , F 1 ) → (G 2 , F 2 ) and a reduction r : G 1 → G 1 , let denote r µ (G 2 ) the reduced graph obtained from G 2 by reducing the pattern µ(G 1 ).

Let assume that the reduction r preserves vertices in V (G 1 )\F 1 and modify (delete or merge) vertices in F 1 (in fact, this is how we choose the set F 1 ). Then, relevant patterns G 1 and their reduction r verify : for any G 2 and cSISOs µ, ν :

(G 1 , F 1 ) → (G 2 , V (G 2 )), either r µ (G 2 ) = r ν (G 2 ), either µ(G 1
) is still a subgraph of r ν (G 2 ) and we can apply r µ to r ν (G 2 ). This last requirement is equivalent to ν(F 1 ) ∩ µ(G 1 ) = ∅.

Patterns

In this part will be presented different patterns that will be reduced in the same set of reaction graphs used previously.

A necessary condition for pattern reduction to be commutative is that the vertices of F 1 form a connected component in G 1 . All patterns presented respect this condition.

In all reaction graphs of this section isolated vertices with edge constraints are circled in red in figures and non isolated vertices without constraints are circled in green in figures.

The Michaelis-Menten pattern and its reduction were already shown in figure 6-3 and figure 6-4.

The distributive Michaelis-Menten pattern shown in figure 6-5 will also be reduced. Its reduced form it also the Michaelis-Menten reduced form from figure 6-4. This distributive kinetics used for the dephosphorylation in a MAPK cascade was studied by Markevich in 2004 [START_REF] Markevich | Signaling switches and bistability arising from multisite phosphorylation in protein kinase cascades[END_REF]. It's a model for the dephosphorylation: The Hill pattern with two distinct binding sites shown in figure 6-6 will also be reduced to the Michaelis-Menten reduced form from figure 6-4. This pattern was studied by Moreland et al. [START_REF] Moreland | Studies on the mechanism of herbicidal action of 2-chloro-4, 6-bis (ethylamino)-s-triazine[END_REF] it represents a two-steps enzymatic mechanism with two binding sites.

E + S ⇔ ES → EP ⇔ E + P . E R 1 R 2 S ES R 3 EP P R -1 R -2
E The partial Hill pattern with twice the same binding site shown in figure 6-7 will also be reduced to the Michaelis-Menten reduced form from figure 6-4. This pattern was studied by Good et al. [START_REF] Good | Inhibitors of the hill reaction[END_REF] it is a variant of the Hill reaction.

R 2 1 R 1 1 SE S ES R 3 P R 2 2 R 1 2

SES

The double Michaelis-Menten pattern with two forms of the enzyme shown in figure 6-8 will also be reduced to the Michaelis-Menten reduced form from figure 6-4.

Implementation Strategy

Pattern reduction is a preprocessing step before searching SEPIs between two graphs. E

R 1 1 S ES R 1 2 P R 2 1 ESS R 2 2

Data structure

The data structure used in the previous SEPI framework was not performant enough and the implementation would have been with a high complexity.

With the previous data structure a graph was represented by a number of vertices, a number of species and a list of pairs representing edges:

Evaluation

The first implementation of the pattern reduction function was made only for variants of the Michaelis-Menten pattern.

Figure 6-9 presents all possible Michaelis-Menten patterns. A Michaelis-Menten pattern is always composed of the three non isolated species E, S and P and at least of three isolated vertices, ES, R 1 and R 2 . The reaction R -1 can be present or not. The reaction R -2 is rarely present. The specie EP and the reaction R 3 are always present together but also rarely present. Results are presented in table 6-2. Theses models are coming from BioModels and they are each represented by a number, which correspond to column "Model" in the table. Classes Ca, Circ and MAPK are composed of eleven models, Cell is composed of seven models. Models appearing in this table are only models where at least one variant of Michaelis-Menten pattern was found.

E R 1 R 2 S ES R 3 EP P R -1 R -2
No patterns were found in the other models. Column "Without R-1" represents the smallest variant of the Michaelis-Menten pattern with only R 1 , ES and R 2 as isolated vertices (red and green vertices in figure 6789). Column "With R-1" represents the pattern with the smallest variant of Michaelis-Menten and the reaction R -1 (red, green and pink vertices in figure 6789). Column "With R-1, R-2, EP" represent the biggest variant of the Michaelis-Menten pattern (all vertices in figure 6789). Other variants of the Michaelis-Menten pattern, reaction R -2 without reaction R -1 for example, were not found in the selected models. The other patterns, Hill reaction, partial Hill reaction and double Michaelis-Menten, were also implemented but were not found among models of the four classes. A few observations can be made from these two tables:

Class

Model With R-1 Without R-1 With R-

Evaluation of all methods combined

• The number of timeouts is noticeably decreasing both between models of different classes and between models of the same class.

• The number of SEPIs inter class is reduced considerably. Especially with models of the class Circadian Clock.

• The size of the set of SEPIs between two given models of the same class is also decreasing.

• The number of SEPIs intra class is also reduced by the combined methods. Furthermore, it was noticed in a lot of SEPI sets that species had, most of the time, the same image through SEPI reductions and only reactions had distinct images. However, biologists are more interested by images of species rather than images of reactions.

Another small SEPI filter was implemented to list only SEPIs with distinct species images. Results can be observed in table 6-6. This filter was tested alone (lines Initial f. distinct s.) and combined with the strict two neighbours merge restriction and maximisation of the number of vertex deletion (line Combined m. distinct s.). Previous results are displayed in the table for comparison.

There is a noteworthy diminution of the SEPI number between two given models, especially for models of the class Calcium Oscillations.

Class

Conclusion

The pattern reduction strategy shows good results with both the diminution of the number of timeout for the SEPI framework and a diminution of the size of the set of SEPIs between two given models. But these good results are only in the class MAPK. Other patterns and Michaelis-Menten patterns were not found in other classes.

Implementing the pattern reduction strategy need a lot of knowledge of the patterns and their reduction and the initial goal of the SEPI framework was to be able to search reductions without knowledge of the model and by using only the reaction graph. Maybe other patterns can be found in the models of BioModels but it would require a deeper analysis.

These results are non the less encouraging and the pattern reduction function is implemented in Biocham. It is now available for users to use them on their own models or models of BioModels.

Bounds

Motivations

SEPI can be used to answer a yes/no question "is graph G a reduced version of G?". However, since SEPI is not a total order, some pairs of graphs are incomparable. Instead of using a yes/no question, bounds (GLB or LUB) of G and G can be computed. It corresponds to model an intersection or an union: to find a set of graphs G c that are a common reduction of G and G or that can be reduced both to G and G .

Definitions

As a reminder G = (V, A) and G = (V , A ) are graphs, and o is an operation set among {m * }, {d} and {m * , d} with m * the restricted merge operation. Proof: It is transitive, reflexive and anti-symmetric.

Definition 7.2.1 (Set of lower bounds) A set of lower bounds is defined by G

∩ o G = {H | G → * o H ∧ G → * o H}
A maximal element of a set X is an element x ∈ X such that ∀y ∈ X, y ≤ x. x is the maximum of X if it is unique.

Definition 7.2.2 (Set of greatest lower bounds (GLB))

G ∩ o G is the set of → * o -maximal elements of G ∩ o G .

Definition 7.2.3 (Set of upper bounds) A set of upper bounds is defined by G

∪ o G = {H | H → * o G ∧ H → * o G } A minimal element of a set X is an element x ∈ X such that ∀y ∈ X, x ≤ y. x is the minimum of X if it is unique. Definition 7.2.4 (Set of least upper bounds (LUB)) G ∪ o G is the set of → * o -minimal elements of G ∪ o G . Observation G ∩ m * d G is not empty.
If G and G are connected and have at least an arc, G ∩ m * G is not empty.

The empty graph is in G ∩ m * d G .

Bounds

If G is connected and has at least an arc (s, r) or (r, s), with s a species vertex and r a reaction vertex, then by connectivity all its species vertices can be recursively merged with s and all its reaction vertices can be merged with r. The same goes for G .

Observation G ∪ m * d G is not empty. It contains G ∪ d G , which is not empty (it contains the disjoint union G G ).
However, for G ∪ m * G , the construction of G × G can give disconnected graphs that are not compatible with the merge restriction.

Observation G ∪ m * G can be empty.

G is defined as G = ({s 1 , s 2 }, {r}, {(s 1 , r), (r, s 2 )}) and G as G = ({s }, {r 1 , r 2 }, {(r 1 , s ), (s , r 2 )}). Figure 7-1 represents both graphs.

s 1 r s 2 r 1 s r 2 Figure 7-1.: Graphs G and G such that G ∪ m * G = ∅.
For a vertex u in a graph H, the set of outgoing neighbours of u is defined by N

• (u) = {v ∈ V | (u, v) ∈ A}, and the set of incoming neighbours by • N (u) = {v ∈ V | (v, u) ∈ A}. N • (u) = ∅ and • N (u) = ∅ are invariant for any merge* operation. Hence, if H is a graph such that H → * m * G, it cannot have species vertex u with both N • (u) = ∅ and • N (u) = ∅.
Moreover, H has to be connected, so there must be a reaction vertex t such that N • (t) = ∅ and • N (t) = ∅. This vertex cannot be sent to neither to r 1 nor r 2 in G , so H cannot verify H → * m * G . Additional properties of GLB and LUB can be found in appendix B. Theorem 7.2.2 [1] The LUB decision problem is NP-complete.

Implementation

Implementing GLB problem or LUB problem is more difficult that the SEPI decision problem because the graph size is not known.

The SEPI GLB problem is the problem of finding a maximal G such that G 1 SEP I G and G 2 SEP I G. Maximality should be with respect to the SEPI partial order, but as a proxy maximality will be with respect to the number of vertices. The encoding of SEPI will be changed to find two graph morphisms µ 1 from G 1 to G and µ 2 from G 2 to G. A MAX-SAT solver will be used to maximise the size of G which is bounded by

m = min(| G 1 |, | G 2 |).
For the following implementation G 1 will be the smallest graph and G 2 the biggest graph. S (respectively R ) represents potential specie vertices (respectively reaction vertices) and is an ordered set s 0 = ⊥ < s 1 < . . . < s p (respectively r 0 = ⊥ < r 1 < . . . < r q ), with p =| S 1 | and q =| R 1 |. s 0 (respectively r 0 ) is the image of deleted specie (respectively reaction) vertices.

Morphisms µ 1 and µ 2 define which vertices of S and R are kept as vertices and arcs of G.

As the coding is the same for µ 1 and µ 2 , only the encoding of µ i will be described.

∀s i ∈ S i , ∀s j ∈ S , i cl(m i s i ,s j =⇒ m i,< s i ,s j+1

),

ii cl(m i,< s i ,s j =⇒ m i,< s i ,s j+1

),

iii cl(m i,< s i ,s j =⇒ ¬m i s i ,s j+1

).

∀r i ∈ R i , ∀r j ∈ R , i cl(m i r i ,r j =⇒ m i,< r i ,r j+1

),

ii cl(m i,< r i ,r j =⇒ m i,< r i ,r j+1

),

iii cl(m i,< r i ,r j =⇒ ¬m i r i ,r j+1

).

III Covering of existing vertices.

∀s ∈ S \ {s 0 }, ∀s i ∈ S i , i cl(e s ⇒ s i ∈S i m i s i ,s ),

ii cl(m i s i ,s ⇒ e s ).

∀r ∈ R \ {r 0 }, ∀r i ∈ R i , i cl(e r ⇒ r i ∈R i m i r i ,r ),

ii cl(m i r i ,r ⇒ e r ).

IV Graph morphism.

∀(u i , v i ) ∈ A i , ∀(u , v ) ∈ V 2 , i cl(m i (u i ,v i ),(u ,v ) =⇒ m i u i ,u ), ii cl(m i (u i ,v i ),(u ,v ) =⇒ m i v i ,v ), iii cl((m i u i ,u ∧ m i v i ,v ) =⇒ m i (u i ,v i ),(u ,v ) ).
V Label epimorphism.

i ∀(u i , v i ) ∈ A i , ∀(u , v ) ∈ V 2 , cl(m i (u i ,v i ),(u ,v ) ⇒ a u ,v ).

ii ∀(u , v ) ∈ V 2 , cl(a (u ,v ) ⇒ (u i ,v i )∈A i m i (u i ,v i ),(u ,v ) ).

VI Existing set normalisation.

∀s ∈ S \ {s 0 , s 1 }, i cl(e s ⇒ e s -1 ).

∀r ∈ R \ {r 0 , r 1 }, i cl(e r ⇒ e r -1

∀(r k , r j ) ∈ R 1 × R i cl(B r k+1 ,r j ⇒ B r k ,r j ), ii cl(B r k ,r j ⇒ m < r k ,r j+1
),

iii cl(B r k ,r j ⇒ B r k ,r j+1 ), iv cl(B r k ,r j ⇒ B r k+1 ,r j+1 ), v cl(B r k ,r j ∧ µ < r k+1 ,r j ⇒ B r k+1 ,r j ).

VIII Bound maximisation.

i ∀v ∈ V , cl(e v ).

Conclusion

Greatest Lower Bound definition in honourable in a logical point of view. It could achieve the goal of finding a common reduced graph between two graph not comparable with the SEPI relation. But the implementation is very heavy and it didn't gave exploitable results. GLB decision problem is more complex than SEPI decision problem and even the SEPI framework was to complex to be computed on big models of BioModels database.

GLB implementation is not usable on big graphs. Moreover, implementing LUB is more complex than implementing GLB because, unlike for a lower bounds, no upper bound of a common graph in known. Thus, no implementation of LUB has been made.

Conclusion Summary

The present diploma thesis is aimed to improve the SEPI framework. This SEPI framework is very powerful and already succeeded to compute reductions between biochemical reaction networks. Never less it has three main limits. First, establishing whether two models are linked through a SEPI is complex and computationally expensive. Second, the number of SEPIs found can be huge, making an analysis of SEPI sets between two given graphs very difficult for biologists. Finally, some existing SEPIs do not have a biological interpretation.

This diploma thesis led to three combined ways to improve the existing framework.

One way consisted to redefine the decision problem into an optimisation problem to select solutions by maximising or minimising the number of vertex deletions. This solution reduced the SEPI set's size between two given graphs and made reduction easier to read for biologists.

A second way was to determine, together with biologists, restrictions on the merge operations in order to filter irrelevant reductions. A strict restriction succeeded to delete some SEPIs without biological interpretations.

Lastly, a preprocessing step has been introduced, consisting of rewriting graphs according to subgraph isomorphism relations. The number of timeouts decreased drastically between models containing specific patterns.

The impact of these three combined implementations has been evaluated on models taken from BioModels' database. Results demonstrated that it contributed to make the SEPI framework more relevant, efficient and functional.

SEPI framework is a part of Biochemical Abstract Machine (BIOCHAM) which is a modelling environment for systems biology. All improvements implemented in this project have been integrated into Biocham and are now available for biologists.

Axes of enhancements

Models from BioModels are usually written in Systems Biology Markup Language (SBML). They come from peer reviewed scientific literature and are not always written according to the same conventions. An upstream standardisation work could improve a detection of SEPI reduction between models.

When only a system of ODEs describes a model, reaction graphs can be automatically inferred with Biocham from expressions of kinetic rates. In some cases, inferred graphs are unexpected, due to not well-chosen kinetic expressions. For example, inhibitors are inferred where there should not be. Thus enhancement of Biocham functions could contribute to better results with SEPI framework.

Finally, reaction graphs used for SEPI are making abstraction of each species' roles in reactions. Inhibition and enzymes are equally treated. Labelling edges accordingly could improve consistency.

B. Additional properties of GLB and LUB

These properties were studied by Eva Philippe. In G (respectively G ), a pre-image of a blue/green part is necessarily a blue/green part in H n . The only possible images for orange vertices in H n are orange vertices in G and G . There can be no merge* between blue and green parts because it would create a double arrow from an orange vertex. But there cannot be either a merge* inside a blue (respectively green) part, because then it could not be matched to G (respectively G ). Thus there can be no sequence of merge*-operations from H n to H, also in G ∪ m * G .

C. Implementation

This diploma thesis was aimed to improve a SEPI framework. This framework is a part of Biochemical Abstract Machine (BIOCHAM) which is a modelling environment for systems biology. It is developed by the Lifeware research team of Inria Saclay-Ile de France.

Biocham is a free software and is publicly available at http://lifeware.inria.fr/biocham4. Models created for the evaluation and models curated from BioModels are available at gitlab.inria.fr/evphilip/sepi-benchmark.

If needed, sources of this report and sources of presentations given during internal seminars are available at gitlab.inria.fr/obargain/diploma-thesis.

Figure 1 - 1 .:

 11 Figure 1-1.: Evolution of concentrations: model with three reactions.

Figure 1 - 2 .:

 12 Figure 1-2.: Evolution of concentrations: model with one reaction.

Figure 1 - 3 .:Figure 1 - 4 .

 1314 Figure 1-3.: Complete Michaelis-Menten reaction graph. E

Figures 1 -

 1 Figures 1-5 and 1-6 are examples of influence networks. They are an abstraction for more detailed biochemical interactions. Each node represents an influence species. A ball-head represents an activation influence, a bar-head represents an inhibition influence, and a simple edge represents an outgoing influence to another node.

Figure 1 - 5 .:

 15 Figure 1-5.: Complex influence network.Figure 1-6.: Simplified influence network.

Figure 1 - 6 .:

 16 Figure 1-5.: Complex influence network.Figure 1-6.: Simplified influence network.

Figure 1 - 7 .:

 17 Figure 1-7.: Concentration of species over time of complex network.

Figure 1 - 8 .:

 18 Figure 1-8.: Concentration of species over time of simplified network.

Figure 1 -

 1 Figure 1-6 is a reduction of figure 1-5.Figure 1-8 and figure 1-7 are corresponding concentration graphs.

Definition 2 . 1 . 2 (

 212 Deletion) Let u ∈ V . The result of the delete operation d u is graph d u (G) = (V , A ) where:

Figure 2 Figure 2 - 2 .Definition 2 . 2 . 1 (

 222221 Figure 2-1.: Complete Michaelis-Menten reaction graph E

2. 5 .

 5 Previous implementation Clauses I Left totality. F l-tot := ∀u∈S∪R cl( u ∈S ∪R ∪⊥ m u,u ). Number of clauses: |S| + |R|, number of new variables: (|S| + |R|) × (1 + |S | + |R |). II Right totality. F r-tot := ∀u ∈S ∪R cl( u∈S∪R m u,u ). Number of clauses: |S | + |R |, number of new variables: 0.

I

  Left totality on arcs. F l-tot-arcs := ∀a∈A cl(deleted a ∨ a ∈A m a,a ). Number of clauses: |A|, number of new variables: |A| × (|A | + 1). II Right totality on arcs. F r-tot-arcs := ∀a ∈A cl( a∈A m a,a ). Number of clauses: |A |, number of new variables: 0.

Number of clauses: 3 ×

 3 |A|, number of new variables: 0.V Redundant morphism propagation.

F

  morph-prop := ∀(u,v)∈A, (u ,v )∈((S ∪R )×(S ∪R ))\A cl(¬m u,u ∨ ¬m v,v ). Number of clauses: |A| × ((|S | + |R |) 2 -|A |), number of new variables: 0. VI Bi-graph constraint. F bi-graph := ∀(u,u )∈(S∪R)×(S ∪R ), type(u) =type(u ) cl(¬m u,u )). Number of clauses: |R| 2 , number of new variables: 0.

Figure 3 - 1 .:

 31 Figure 3-1.: Number of SEPIs between models in Calcium Oscillations class.

I

  Minimisation of the number of deletions. ∀u ∈ S ∪ R, cl(¬m u,⊥ ). Number of clauses: |S| + |R|, number of new variables: 0. II Maximisation of the number of deletions. ∀u ∈ S ∪ R, cl(m u,⊥ ). Number of clauses: |S| + |R|, number of new variables: 0.

Figure 3 - 4 .

 34 Figure 3-4.: EE-PS: combination of two Michaelis-Menten, with shared species.

Figure 3 - 5 .:

 35 Figure 3-5.: Completely reduced three-levels MAPK cascade.

Figure 4 -

 4 Figure 4-1.: All SEPIs inter and intra classes.

Definition 4 . 2 . 1 (Figure 4 - 2 .:Definition 4 . 2 . 2 (BFigure 4 - 3 .:

 4214242243 Figure 4-2.: Example of two-neighbours vertices.

Figure 4 -Figure 4 - 4 .:

 444 Figure4-4 is an illustration of the good-path definition. Dotted edges represent a good-path between two vertices. There exists a good-path between vertices a and b because there exists three vertices u, v and c such that a is related to the reaction u, c is related to the reaction v and there exists a good-path between reactions u and v and there exists a good-path between species c and b.

Figure 4 -

 4 Figure 4-5 is an example of cyclic configuration. In this disconnected graph, a graph where there exist two nodes such that no path has those nodes as endpoints, A and B have the same image, and r1 and r2 as well. The argumentation that there is a path* between A and B with C 1 a path* between r1 and r2, itself defined thanks to A and B, would be infinitely recursive, no base case is defined. Hence, no length could be defined. It shows why a well-defined length is required for paths*.

Definition 4 . 2 . 3 AFigure 4 - 5 .:

 42345 Figure 4-5.: Example of a cyclic configuration.

Figure 4 -

 4 Figure 4-3 was an example of good-path reduction because vertices with the same image through the morphism are linked by a good-path.

Figure 4 -BFigure 4 - 6 .:Theorem 4 . 2 . 1

 446421 Figure 4-6 is an example of invalid good-path reduction. Species A and B can't have the same image in the targeted graph because reactions r1 and r2 are not merged together.

  induction step: |V | > |V | Induction hypothesis: the implication is true for any G such that |V | = |V | -1. -case (i): ∃c such that µ(c) = ⊥ G will be defined by G := d c (G) and ∼ µ by:∼µ : G -→ G x -→ µ(x) ∼ µ is a SEPI* from G to Gand there is a sequence from G to G by induction hypothesis. µ is a SEPI* from G to G and G → d G → * m * d G . Hence there is a sequence from G to G . -case (ii): ∃a, b ∈ V such that µ(a) = µ(b) = ⊥ Distance d(a, b) will be defined as the length of the smallest path* between a and b in µ. a and b will be selected such that d(a, b) is the smallest distance in G with associated path* {x 0 = a, . . . , x n = b}. The minimality of d(a, b) implies that n = 1 and length(C) = 0, i.e a and b are two-neighbours. G is defined by G := m * a,b (G) and ∼ µ by:

  with paths* given by µ where vertices a and b are replaced by ab.∼ µ is a SEPI* from G to G and there is a sequence from G to G by induction hypothesis. µ is a SEPI* from G to G and G → m * G → * m * d G .Hence there is a sequence from G to G .(⇐) The if direction will be proven by induction on n, the length of the sequence of merge and delete operation between G and G . • base case: n = 0 If the sequence is empty, by definition G and G are isomorphic, thus an isomorphism µ between them is a SEPI* (since µ(a) = µ(b) ⇔ a = b for any vertices a, b of G).

∼

  µ(b), hence the same path* works for µ and ∼ µ. Then ∼ µ is a SEPI*.

∼

  µ(b) hence the same path* works for µ and ∼ µ for a and b. For a, b ∈ V such that ∼ µ(a) = ∼ µ(b) = uv. If µ(a) = µ(b), the same paths* work for µ and ∼ µ for a and b. u and v are defined by ∼ µ(a) = u and ∼ µ(b) = v. The possibility to apply merge* operations to u and v implies that there exists w ∈ V such that {u, w} ∈ E , and {v, w} ∈ E . Let {a , c 1 } and {b , c 2 } be antecedents of edges {u, w} and {v, w}, with µ(a ) = u, µ(b ) = v and µ(c 1 ) = µ(c 2 ) = w. The concatenation of paths* between a and a , c 1 and c 2 , and b and b given by µ gives a path* from a to b for ∼ µ. Then ∼ µ is a SEPI*. This conclude the proof that G SEP I * G ⇔ G → * m * d G . Property 4.2.1 Merge* and delete operations enjoy the commutation and association properties of figure 4-7, where the solid arrow represents universal quantification and the dashed arrow represents existential quantification. Merge* and delete operation are commutative because changing the order of the operands does not change the result.

Figure 4 - 7 .:

 47 Figure 4-7.: Properties of merge* and delete operations.

  subsets U ⊆ P(E), and an integer k ≤ |U |. Question. Is there a family U * ⊆ U such that U * = k and s∈U * s = E? (E, U, k) is an instance of exactly-k-set-covering, E = {e 1 , . . . , e n } and U = {u 1 . . . , u m }. G = (S, R, A) is defined with S = {u 1 , . . . , u m } ∪ {e 1 , . . . , e n } ∪ {w}, R = {r 1 , . . . , r m }, and

Figure 4 -

 4 Figure 4-8 is, for example, an instance of the set covering problem (E, U, k) such that E = {a, b, c, d, e, f } and U = {{a, c, d}, {a, b, d}, {c, f }, {e, f }}.Figure 4-8(a) displays the corresponding source graph, figure 4-8(b) the target graph corresponding to k = 3, and figure 4-8(c) the target graph corresponding to k = 2.

Figure 4 -

 4 Figure 4-8 is, for example, an instance of the set covering problem (E, U, k) such that E = {a, b, c, d, e, f } and U = {{a, c, d}, {a, b, d}, {c, f }, {e, f }}.Figure 4-8(a) displays the corresponding source graph, figure 4-8(b) the target graph corresponding to k = 3, and figure 4-8(c) the target graph corresponding to k = 2.

Figure 4 - 8 .:

 48 Figure 4-8.: Instance of the set covering problem (E, U, k).

For

  all a, b ∈ V , ∀k ∈ 0, min(|S|, |R|) , p a,b,k represents the existence of a good-path of length k between a and b. p a,b,k = 1 if and only if there exists a good-path of length k between a and b. m a,b,y = 1 if and only if m(a) = y and m(b) = y. Vertices a and b have the same image y through the morphism m. p a,b,k 1 ,c,d,k 2 = 1 if and only if p a,b,k 1 = 1 and p c,d,k 2 = 1. Variables m a,b,y and p a,b,k 1 ,c,d,k 2 are introduced mainly to remove conjunctions from disjunctions of literals.

  -path of length k clauses are the recursive definition's inductive step. Path* clauses are linking together the merge definition and the good-path definition, restricting the merge definition. I Inductive merge definition: i F merge-i := ∀a,b∈V ∀y∈V ∪{⊥} cl(¬m a,b,y ∨ m a,y ), ii F merge-ii := ∀a,b∈V ∀y∈V ∪{⊥} cl(¬m a,b,y ∨ m b,y ), iii F merge-iii := ∀a,b∈V ∀y∈V ∪{⊥} cl(¬m a,y ∨ ¬m b,y ∨ m a,b,y ). II Good-path of length zero: i F path-0-i := ∀a,b∈V cl(p a,a,0 ), ii F path-0-ii := ∀a,b∈V cl(¬p a,b,0 ). III Inductive good-path definition: i F inductive-path-i := ∀a,b,c,d∈V ∀k 1 ,k 2 ∈ 0,(min(|S|,|R|)-1)

I 2 × 2 × 2 + 2 + 2 + 2 + 2 × 2 × 2 + 2 + 2 + 2 + 2 + 2 ) 2 × 2 ×

 2222222222222222 Inductive merge definition. Clauses are only written when a < b and type(a) = type(b) = type(y). II Good-path of length zero. Clauses are only written when a ≤ b and type(a) = type(b). III Inductive good-path definition. Clauses are only written when a ≤ b, c ≤ d, type(a) = type(b), type(c) = type(d) and type(a) = type(c). IV Good-path of length k. Clauses are only written when a < b and type(a) = type(b). V Path*. Clauses are only written when a < b and type(a) = type(b) = type(y). Number of variables Optimisation reduced the number of variables. For m a,b,y , |S|×(|S|-1) (|S | + 1) + |R|×(|R|-1) (|R | + 1) variables are introduced. For p a,b,k , ( |S|×(|S|-1) |S| + |R|×(|R|-1) |R|) × (min(|S|, |R|) + 1) variables are introduced. For p a,b,k 1 ,c,d,k 2 , ( |S|×(|S|-1) |S|) × ( |R|×(|R|-1) |R|) × 2 × (min(|S|, |R|) + 1) 2 variables are introduced.Number of clausesThe optimisation also reduced the number of clauses added for the definition of SEPI*:I Inductive merge definition: 3 × ( |S|×(|S|-1) (|S | + 1) + |R|×(|R|-1) (|R | + 1)) clauses.II Good-path of length zero: |S|×(|S|-1) |R|×(|R|-1) |V | clauses. III Inductive good-path definition: 3 × ( |S|×(|S|-1) |S|) × ( |R|×(|R|-1) |R|) × 2 × (min(|S|, |R|) + 1) 2 clauses. IV Good-path of length k: 2 × ( |S|×(|S|-1) |R|×(|R|-1) × min(|S|, |R|) clauses. V Path*: |S|×(|S|-1) (|S | + 1) + |R|×(|R|-1) (|R | + 1) clauses.

2 × 2 × 2 + 2 +

 2222 As a reminder, |S| (respectively |S'|) is the number of species of the initial graph (respectively targeted graph). |R| (respectively |R'|) is the number of reactions of the initial graph (respectively targeted graph). |V| is the number of vertices, |V | = |S| + |R| and |V | = |S | + |R |. Orders of magnitude are given with n = max(|V |, |V |). New variables Exact number m a,b,y |S|×(|S|-1) (|S | + 1) + |R|×(|R|-1) (|R | + 1) p a,b,k ( |S|×(|S|-1) |S| + |R|×(|R|-1) |R|) × (min(|S|, |R|) + 1)

2 × 2 × 2 + 2 + 2 + 2 + 2 Path* of length k 2 × 2 + 2 ) 2 × 2 × (|R | + 1 )

 2222222222221 merge def. 3 × ( |S|×(|S|-1) (|S | + 1) + |R|×(|R|-1) (|R | + 1)) Path* of length zero |S|×(|S|-1) |R|×(|R|-1) |V | Inductive path* def. 3 × ( |S|×(|S|-1) |S|) × ( |R|×(|R|-1) |R|) × 2 × (min(|S|, |R|) + 1) ( |S|×(|S|-1) |R|×(|R|-1) × min(|S|, |R|) (|S | + 1) + |R|×(|R|-1) Order of magnitude Inductive merge def. O(n 3 ) Path* of length zero O(n 2 ) Inductive path* def. O(n 6 ) Path* of length k O(n 3 ) Path* O(n 3 )Table 5-2.: Number of clauses added by the accurate merge restriction.

Figure 5 -Figure 5 - 1 .:

 551 Figure 5-1 is an example of this impossibility. It would make sense to merge together reactions r 1 and r 2 and merge together reactions r 3 and r 4 . And then merge together A, C and E. These species are far apart from each other and species B and D don't have the same image in the targeted graph. Thus, merging A, C and E can't be defined locally.

Figure 5 -Figure 5 - 2 .:

 552 Figure 5-2 is an example of simple not allowed vertex merge that would have been nice to keep

  These clauses are added to the already existing SEPI implementation's clauses. With this stronger restriction, no new variables need to be defined and it is introducing less than |V | 2 × |V | new clauses. neigh(a, b) is written when a and b are two-neighbours. Clauses are defined by F merge := ∀a,b∈V ¬neigh(a,b) ∀y∈V cl(¬m a,y ∨ ¬m b,y ).

Figures 5 -

 5 Figures 5-3 to 5-6 are a graphical presentation of the results. In all four figures, yellow vertices represent MAPK models, green vertices represent models of the Circadian Clock class, blue vertices represent Cell Cycle models and purple vertices represent models of the Calcium Oscillation class.

Figure 5 - 3 .:

 53 Figure 5-3.: SEPIs intra class without merge restriction.

Figure 5 - 4 .:

 54 Figure 5-4.: SEPIs intra class with strict merge restriction.

Figure 5 - 5 .:

 55 Figure 5-5.: SEPIs inter class without merge restriction.

Figure 5 - 6 .:

 56 Figure 5-6.: SEPIs inter class with strict merge restriction.

Figure 6 Figure 6 - 2 .

 662 Figure 6-1.: Unexpected SEPI: initial graph. aE

Figure 6 - 3 .:Figure 6 - 4 .:

 6364 Figure 6-3.: Reaction graph with expended Michaelis-Menten pattern. A R 3 E

Figure 6 - 5 .:

 65 Figure 6-5.: Distributive Michaelis-Menten reaction graph.

Figure 6 - 6 .:

 66 Figure 6-6.: Hill reaction graph.

2 ESSFigure 6 - 7 .:

 267 Figure 6-7.: Partial Hill reaction graph.

Figure 6 - 8 .: 2 , 2 .

 6822 Figure 6-8.: Double Michaelis-Menten reaction graph.

  Graph = [|V|, |S|, list_edges] With V a set of vertices and |V| it's cardinality and S a set of species vertices and |S| it's cardinality. Vertices were represented by integers and ordered. The set of specie vertices is 0, |S| -1 . The set of reaction vertices is |S|, |V | -1 .

Figure 6 - 9 .:

 69 Figure 6-9.: All variants of the Michaelis-Menten pattern.

Figure 6 -

 6 Figure 6-10.: SEPIs on MAPK without pattern reduction.

Figure 6 -

 6 Figure 6-11.: SEPIs on MAPK with pattern reduction.

Property 7 .

 7 2.1 (→ * o partial order) The arrow → *o defines a partial order on the space of graphs.

Figure 7 -Theorem 7 . 2 . 1 Definition 7 . 2 . 6 (

 7721726 Figure 7-2 gives examples of graph pairs such that G ∩ o G and G ∪ o G are not singletons. Graphs G, G , H and H show there is no unicity of GLB and LUB for the operation o = m * . Relations → * m * are indicated with dashed arrows, the fusions can be recognized with the names of H and H vertices. In this example, H, H ∈ G ∩ m * G and G, G ∈ H ∪ m * H .

Figure 7 - 2 .:

 72 Figure 7-2.: Illustration no unicity.

Figure A- 1 .

 1 Figure A-1.: All SEPIs between models of Cell Cycle.

Figure A- 2 .:

 2 Figure A-2.: Cell Cycle -Minimisation of the number of deletions.

Figure A- 3 .:

 3 Figure A-3.: Cell Cycle -Maximisation of the number of deletions.

Figure A- 4 .:

 4 Figure A-4.: All SEPIs between models of Calcium Oscillations.

Figure A- 5 .

 5 Figure A-5.: Calcium Oscillations -Minimisation of the number of deletions.

Figure A- 6 .

 6 Figure A-6.: Calcium Oscillations -Maximisation of the number of deletions.

Figure A- 7 .

 7 Figure A-7.: All SEPIs between models of Circadian Clock.

Figure A- 8 .

 8 Figure A-8.: Circadian Clock -Minimisation of the number of deletions.

Figure A- 9 .

 9 Figure A-9.: Circadian Clock -Maximisation of the number of deletions.

Figure A- 10

 10 Figure A-10.: All SEPIs between models of MAPK.

Figure A- 11

 11 Figure A-11.: MAPK -Minimisation of the number of deletions.

Figure A- 12

 12 Figure A-12.: MAPK -Maximisation of the number of deletions.

Property B.0. 1 G

 1 ∩ o G is finite. Proof: Delete and merge reduce the number of vertices, and graphs with less than min(| G |, | G |) vertices are in finite number. Property B.0.2 G ∪ m * G is not finite. It is unknown either G ∪ m * d G (i.e SEPI*) is finite or not. Proof: Figure B-1 and figure B-2 give an example of two graphs G and G such that G ∪ m * G is infinite.

Figure B- 3

 3 Figure B-3 illustrate how an infinity of graphs can be build. ∀n ∈ N\{0}, by concatenating n times a certain loop pattern, there is an infinity of graphs H n in G ∪ m * G . A sequence of merge* operations that shows that H n → * m * G (respectively G ) is given by folding the light-green parts by merging vertices a, a and b, b (respectively the light-blue parts), then merging the orange vertices, which have become 2-neighbours. Hence, H n ∈ G ∪ m * G .

1 a 1 a

 11 Figure B-1.: Graph G.

Figure B- 3 .:

 3 Figure B-3.: Infinite set illustration.

Definition 2.4.1 (MAX-SAT problem) A

  A weighted clause is a pair (C i , w i ), where C i is a clause and w i a positive integer representing its weight. A weighted CNF formula, WCNF formula for short, is a set of weighted clauses. The weighted version of MAX-SAT problem is to find an interpretation that maximise the combined weight of the satisfied clauses. Hard clauses are clauses that must be satisfied, soft clauses are non-mandatory clauses. A cost of an interpretation is the sum of weights of clauses that are unsatisfied.

	Input of a weighted MAX-SAT solver is also a text file in a simplified version of DIMACS
	format. Besides comments, first line of a file should be p cnf nb_variables nb_clauses top, where
	nb_variables is an upper bound on the largest index of a variable appearing in the file or the
	exact number, nb_clauses is the exact number of clauses and top is an integer representing the
	maximum weight. Hard clauses have weight top and soft clauses have a weight smaller than
	top.
	An example of MAX-SAT solver input file:
	c comments
	p wcnf 4 5 16
	16 1 -2 4 0
	16 -1 -2 3 0
	8 -2 -4 0
	4 -3 2 0
	3 1 3 0
	Maximum Satisfiability problem, MAX-SAT prob-
	lem for short, is the problem of determining a maximum number of clauses for which there exists
	an interpretation that satisfies a given CNF formula.

Table 3

 3

	-1 presents all possible SEPIs from initial graph of figure 3-2 to targeted graph of figure
	3-3. Each line represents a SEPI. Each column represents an initial graph's vertex and for each
	SEPI it's image in the targeted graph is presented. ⊥ represents a vertex deletion.
	It can be observed in table 3-1 that, even on a small example, there are 16 solutions but only
	4 that minimise the number of deletions. It is interesting for biologists to see only SEPIs that
	minimise or maximise the number of vertex deletion. Minimising indicates what has to be deleted
	imperatively and where vertices can be merged. Maximising indicates all vertices that can be
	deleted and makes a SEPI more readable.

Table 3 -

 3

1.: All reductions between models.

Table 3 -

 3 2 shows the same effects of reducing the number of pairings with Min ⊥ and Max ⊥ for combinations of two Michaelis-Menten patterns. It decreases the number of pairs with more than 200 SEPIs and decrease the mean and median of the number of pairings.

		Pairs with SEPIs	Number of SEPIs	No SEPI
	Type of pairing	≥ 200	< 200	mean	median
	Normal	99 (2.1%) 508 (10.8%) 20.1	2	4085 (87%)
	Min ⊥	6 (0.3%) 601 (12.8%) 12.4	2	4085 (87%)
	Max ⊥	6 (0.3%) 601 (12.8%)	8.4	2	4085 (87%)

Table 3 -2.:

 3 Combinations of two Michaelis-Menten.

Table 3 -

 3

	Number of SEPIs No SEPI Timeouts

3.:

Statistics on the sets of pairings for 1, 2 and 3-levels MAPK cascades.

Table 3 -

 3 3 presents statistics on sets of pairings searched among all one, two and three-levels MAPK cascades. Again, Min ⊥ and Max ⊥ filter reductions but Max ⊥ causes many timeouts.

	E1
	KKK	pKKK
	E2
	KK	pKK	ppKK
		KKpase
		K	pK	ppK
			Kpase

Table 3 -4.:

 3 Reaction graph characteristics for each class.

	Number of arcs

Table 3 -5.:

 3 Number of SEPIs in each class.

	Graphs 3-6, 3-7 and 3-8 used in this example could represent a one-step kinase double phospho-
	rylation, with all dephosphorylations. It is composed of four Michaelis-Menten reactions. Figure
	3-7 presents the reduced form, Michaelis-Menten reactions are recognizable with double-arrows
	towards enzymes.

A SEPI µ with minimisation of the number of vertex deletions is expected between both graphs. µ is represented by colours in the figures. Vertices of one colour in figure

3

-6 are mapped to the vertex of the same colour in figure

3-7

. Colours represent an expected SEPI pairing from expanded form to reduced one, where the four reverse reactions from complexes to enzyme and substrate are deleted, like in Michaelis-Menten reduction.

Table 3 -6.:

 3 Percentage of timeout in each class.

Table 3 -

 3

	Model	Normal Reductions					Bottom Minimisation
				Glucose				RC2		OPEN-WBO		QMAXSAT
		nb SEPI	1st 2nd	last nb SEPI	1st	2nd	last	1st	2nd	last	1st	2nd	last
	043-045	6	13	7	7	2	173	82	(75)	95	14	(11)	61		(4)
	145-098	200	7	6	196	200	72	69	96	10	7	139	2		3
	144-008	200	104	97	248	200	691	726	826	577	702 29773	426		706
	021-171	2	326	808 (1728)	2	383	378 (4191)	603	1871 (2229)	421	(3430)
	083-084	200	202	87	247	200	806	783	769	608	829 38823	305		275
	029-027	10	12	12	8	6	103	105	104	19	18	21	3		3
	011-026	200 4381 4501 61875	200 11745 11690 28507 39617 39614 94653 30641 62250

7.: Solver comparison on Bottom minimisation and maximisation (time in s).

Table 3 -8.:

 3 Time to compute the first second and last SEPI for each solver (time in ms).

Table 4 -1.: Number

 4 of SEPIs inter and intra class.

			SEPI	No SEPI Timeouts
	Intra class	372	88	243	41
			23,66%	65,32%	11,02%
	Inter class 1188	321	776	91
			27.02%	65.32%	7.66%

Table 4 -

 4 1 represents the number of SEPIs inter and intra class. The column "Pairs" gives the amount of ordered model pairs for each problem. For example, model 74 to model 73 is considered as one pair, a SEPI relation will be searched from model 74 to model 73 and there exists one.

Table 4 -

 4 2 compares intra class results. Table 4-3 shows inter class results from 2019. Precise results from 2015 on inter class performance are not available, the number of 9% for the number of SEPIs inter class is the only indication.

Table 4 -

 4 2 and table 4-3 present results on old models. These models are also ordered in four classes, Ca represents the class Calcium Oscillations, Cell represents the class Cell Cycle, Circ represents the class Circadian Clock and MAPK represents the class Mitogen-Activated Protein

				2015 results				2019 results
			SEPI	No SEPI Timeouts	SEPI	No SEPI Timeouts
	Ca	(110)	38	72	0		38	72	0
			34.55% 65.45%	0%		34.55% 65.45%	0%
	Cell	(72)	12	51	9		12	49	11
			16.67% 70.83%	12.50%		16.67% 68.06%	15.27%
	Circ	(110)	37	73	0		25	61	24
			33.64% 66.36%	0%		22.73% 55.45%	21.82%
	MAPK (110)	38	63	9		40	60	10
			34.55% 57.27%	8.18%		36.36% 54.55%	9.09%
	Total	(402)	125	259	18		115	242	45
			31.09% 64.43%	4.48%		28.61% 60.20%	11.19%
	Table 4-2.: Comparison with results of 2015: number of SEPI relations intra class.
						2019 results
					SEPI	No SEPI Timeouts
		Ca	Cell	(198)	79		114	5
				39.90% 57.57%	2.53%
		Ca	Circ	(242)	100		133	9
				41.32% 54.96%	3.72%
		Ca	MAPK (242)	82		155	5
				33.88% 64.05%	2.07%
		Cell Circ	(198)	31		129	38
				15.66% 65.15%	19.19%
		Cell MAPK (198)	48		115	35
				24.24% 58.08%	17.68%
		Circ MAPK (242)	55		147	40
				22.73% 60.74%	16.53%
			Total	(1320)	395		793	132
				29.92% 60.08%	10.00%

Table 4 -3.:

 4 Comparison with results of 2015: number of SEPI relations inter class.

Table 4 -

 4 4 presents characteristics of the recent version of models. It can be noticed that the number of edges is higher compared to old models (same characteristics on previous models were displayed in table[3][4]. And two models of the Cell Cycle class are missing, they were too big to be curated by the latest version of Biocham.

			Number of vertices	Number of arcs
	Class	Nb models Min Max Mean Min Max Mean
	Ca Oscillations	11	8	44	16	11	72	28
	Cell Cycle	7	20	189	70	38	364	155
	Circadian Clock	11	25	82	55	38	130	86
	MAPK	11	10	334	63	29	744	137

Table 4 -4.:

 4 Reaction graph characteristics for each class. Results presented in table 4-5 and table 4-6 will be used as benchmarks for further tests. The important result of table 4-5 is the number of timeouts (11.02%) which need to be decreased. Important results of table 4-6 is the number of SEPIs inter class (27.02%) and the number of timeouts (7.66%). Both also need to be decreased.To be able to observe a diminution of the number of SEPIs between two given models, it is also specified in table 4-5 the number of SEPI sets that are below a size of 200. Results on the new set of models are similar to results on previous models with the same SEPI implementation. With just a lower amount of timeouts and SEPI relations inter class and intra class.

			SEPI	No SEPI	Timeout
	Class	≥ 200 < 200	total
	Ca (110)	20	6	26 (23.64%) 83 (75.45%)	1 (00.91%)
	Cell (42)	8	0	8 (19.05%)	28 (66.67%)	6 (14.28%)
	Circ (110)	16	3	19 (17.27%) 70 (63.64%) 21 (19.09%)
	MAPK (110)	35	0	35 (31,82%) 62 (56,36%) 13 (11,82%)
	Total (372)	79	9	88 (23,66%) 243 (65,32%) 41 (11,02%)

Table 4 -5.:

 4 SEPIs intra class without restriction.

	Results inter class	SEPIs	No SEPIs	Timeouts
	Ca	Cell	(154)	43 (27.92%) 109 (70.78%) 2 (01.30%)
	Ca	Circ	(242)	100 (41.32%) 136 (56.20%) 6 (02.48%)
	Ca	MAPK (242)	74 (30.58%) 162 (66.94%) 6 (02.48%)
	Cell Circ	(154)	22 (14.29%) 109 (70.78%) 23 (14.93%)
	Cell MAPK (154)	30 (19.48%) 105 (68.18%) 19 (12.34%)
	Circ MAPK (242)	52 (21.49%) 155 (64.05%) 35 (14.46%)
		Total	(1188) 321 (27.02%) 776 (65.32%) 91 (7.66%)

Table 4 -

 4

6.: SEPIs inter class without restriction.

Table 5 -

 5 1 and table 5-2 summarise the complexity of the previous implementation. This implementation is based on the SEPI framework, new clauses and new variables are added to already existing clauses and variables of the SEPI framework to restrict the merge operation.

Table 5 -1.:

 5 Number of variables added by the accurate merge restriction.

Table 5 -3.:

 5 SEPIs for combinations of Michaelis-Menten.

				with SEPIs	Number of SEPIs	No SEPI
	Restriction		≥ 200	< 200	mean	median
	Expanded and intermediary
	Without		0	158 (07.60%) 1.65	2	1912 (92.40%)
	Strict		0	23 (01.00%)	1.6	2	2047 (99.00%)
	Intermediary and Reduced
	Without		7 (00.30%) 225 (10.90%) 17.5	4	1838 (88.80%)
	Strict		0	76 (03.70%)	10.3	6	1994 (96.30%)
	expanded and Reduced
	Without		92 (04.40%) 260 (12.60%)	25	4	1718 (83.00%)
	Strict		19 (00.90%) 173 (08.40%)	29	18	1878 (90.70%)
	Total
	Without		99 (01.60%) 643 (10.40%) 16.4	2	5468 (88.00%)
	Strict		19 (00.30%) 272 (04.40%) 21.4	10	5919 (95.30%)
		Pairs with SEPIs	Number of SEPIs	No SEPI	Timeouts
	Reduction	≥ 200	< 200	mean	median
	Without	64 (30.50%) 9 (04.00%)	25	14	127 (60.50%) 10 (05.00%)
	Strict	54 (26.00%) 21 (10.00%)	60	39	135 (64.00%)	0

Table 5 - 4

 54

.: SEPIs for MAPK cascades.

Table 5 -

 5 5.: SEPI relations intra class with and without strict merge restriction.

	No SEPI	Timeout

Table 5 -

 5 6.: SEPI relations inter class with and without strict merge restriction.

Table 6 -

 6 1 is a good example. It represents all possible reductions between the complete Michaelis-Menten reaction graph represented in figure 6-3 and the reduced Michaelis-Menten reaction graph represented in figure6-4. One line represents a SEPI. One column represents one vertex in the initial graph and it's values in the image graph through each SEPIs. It can already be observed that there is three different SEPIs just for the reduction of a Michaelis-Menten reaction.

Table 6 -1.:

 6 All possible reductions between complete and reduced Michaelis-Menten reaction graphs.

  1. Search known patterns in G 1 and G 2 , 2. Rewrite graphs to obtain graphs G 1 and G 2 by replacing expanded patterns by their reduced form, 3. Search SEPIs between G 1 and G 2 .

  (G 2 , F 2 ) is written when there exists a cSISO from (G 1 , F 1 ) to (G 2 , F 2 ). The partial order given by cSISO is not a well quasi-order.Proof: The set {(C n , ∅), n ≥ 3}, where C n is the cycle graph with n vertices, is an infinite antichain.

	Definition 6.2.1 (Constrained subgraph isomorphism) A constrained subgraph isomor-phism (cSISO) from G 1 to G 2 constrained by F 1 and F 2 is a function µ : G 1 → G 2 such that: 1. µ(G 1 ) is an induced subgraph of G 2 , 2. µ is a graph isomorphism from G 1 to µ(G 1 ), This definition extends easily to bipartite directed graphs. (G 1 , F 1 ) Property 6.2.1 Definition 6.2.2 (Constrained subgraph isomorphism problem.) The constrained sub-graph isomorphism problem is the decision problem: Instance: Two graphs and associated subsets of vertices (G 1 , F 1 ), (G 2 , F 2 ) Question: (G 1 , F 1 ) cSISO (G 2 , F 2 )? 3. µ(F 6.3. Patterns Theorem 6.2.1 The constrained subgraph isomorphism problem is NP-complete.

1 ) ⊂ F 2 , and

4. for all u ∈ F 1 , for all z ∈ V 2 , (µ(u), z) ∈ E 2 =⇒ z ∈ µ(G 1 ) (all neighbours of µ(u) are images of neighbours of u) cSISO

Table 6 -2.:

 6 Number of Michaelis-Menten patterns.It can be observed in table 6-2 that only three models outside of the MAPK class were concerned by Michaelis-Menten patterns. Thus after reducing all variants of the Michaelis-Menten pattern, tests with the SEPI framework were only conducted on the MAPK class. Tests were performed with the initial SEPI framework, without merge restrictions and without extremalisation of the number of vertex deletion. Results can be observed in table6-3.

	1, R-2, EP

Table 6 -3.:

 6 SEPIs in MAPK class.

Table 6 -

 6 3 presents results for the 110 pairs of models of the MAPK class. Figures 6-10 and 6-11 also represent the results.Three observations can be made from the table and the figures:-The number of set of SEPIs which size is bellow 200 goes from zero to six. Reducing this size was the goal of pattern reduction.-The number of timeout of considerably decreasing which is a very encouraging result.-The number of total SEPIs is decreasing, which can be concerning. In fact, models 27, 29 and 31 of the MAPK class are variant of respectively models 26, 28 and 30. Markevich et al. also tried to reduce patterns but didn't use the same methodology, that's why, after reducing Michaelis-Menten patterns, some SEPIs with models 27, 29 and 31 cannot be found any more.

Table 6 -

 6 4 and table 6-5 present results of the updated SEPI framework with the three improvements: maximisation of the vertex deletion number, strict two-neighbours merge restriction and reduction of Michaelis-Menten patterns.

Table 6 -

 6 4.: SEPI relations intra class with combined restrictions.

	6.6. Evaluation of all methods combined
					SEPI		No SEPI	Timeout
	Class		Restriction ≥ 200 < 200	total		(20 min)
	Ca (110)		Without	20	6	26 (23.64%) 83 (75.45%)	1 (00.91%)
			Combined	7	10	17 (15.45%) 93 (84.55%)	0
	Cell (42)		Without	8	0	8 (19.05%)	28 (66.67%)	6 (14.28%)
			Combined	1	0	1 (02.38%)	35 (83.33%)	6 (14.29%)
	Circ (110)	Without	16	3	19 (17.27%) 70 (63.64%) (19.09%)
			Combined	0	1	1 (00.91%) 105 (95.45%) 4 (03.64%)
	MAPK (110) Without	35	0	35 (31,82%) 62 (56,36%) (11,82%)
			Combined	1	10	11 (10.00%) 89 (80.91%) (09.09%)
	Total (372)	Without	79	9	88 (23,66%) 243 (65,32%) (11,02%)
			Combined	9	21	30 (08.06%) 322 (86.56%) (05.38%)
	Pair of classes	Restriction		SEPI		No SEPI	Timeout
	Ca	Cell	(154)	Without	43 (27.92%)		109 (70.78%)	2 (01.30%)
				Combined	16 (10.49%)		126 (81.82%) (07.69%)
	Ca	Circ	(242)	Without	100 (41.32%) 136 (56.20%)	6 (02.48%)
				Combined	12 (04,96%)		230 (95.04%)	0
	Ca	MAPK (242)	Without	74 (30.58%)		162 (66.94%)	6 (02.48%)
				Combined	5 (02.07%)		226 (93.39%) (04.54%)
	Cell Circ	(154)	Without	22 (14.29%)		109 (70.78%) (14.93%)
				Combined		0		144 (93,51%) (06,49%)
	Cell MAPK (154)	Without	30 (19.48%)		105 (68.18%) (12.34%)
				Combined	4 (02.60%)		133 (86.36%) (11.04%)
	Circ MAPK (242)	Without	52 (21.49%)		155 (64.05%) (14.46%)
				Combined		0		231 (95.45%) (04.55%)
	Total	(1188) Without	321 (27.02%) 776 (65.32%) (07.66%)
				Combined	37 (3.11%)	1090 (91.75%) 61 (5.14%)

Table 6 -5.:

 6 SEPI relations inter class with combined restrictions.

Table 6 - 6

 66

		Filter	≥ 200 < 200 mean
	Ca	Initial framework	20	3	58.00
		Initial f. distinct s.	8	15	7.33
	Ca	Combined methods	7	10	19.20
		Combined m. distinct s.	1	16	7.19
	Cell	Initial framework	8	0	-
		Initial f. distinct s.	7	1	1.00
	Cell	Combined methods	1	0	-
		Combined m. distinct s.	0	1	15.00
	Circ	Initial framework	16	3	13.67
		Initial f. distinct s.	16	3	1.67
	Circ	Combined methods	0	1	6.00
		Combined m. distinct s.	0	1	2.00
	MAPK Initial framework	35	0	-
		Initial f. distinct s.	34	1	12
	MAPK Combined methods	1	10	29.60
		Combined m. distinct s.	0	11	18.36
	Total	Initial framework	79	6	35.83
		Initial f. distinct s.	66	19	6.39
	Total	Combined methods	9	21	23.52
		Combined m. distinct s.	1	29	11.52

.: SEPI relations intra class.

  ). VII Symmetry breaking.∀(s k , s j ) ∈ S 1 × S i cl(B s k+1 ,s j ⇒ B s k ,s j ), ii cl(B s k ,s j ⇒ m < s k ,s j+1 ), iii cl(B s k ,s j ⇒ B s k ,s j+1 ), iv cl(B s k ,s j ⇒ B s k+1 ,s j+1 ), v cl(B s k ,s j ∧ µ < s k+1 ,s j ⇒ B s k+1 ,s j ).

  All improvements implemented in this project are integrated into Biocham and are available from version 4.2.6 at gitlab.inria.fr/lifeware/biocham. Documentation on how to execute the developed functions are in section 7.3. "Detecting model reductions" and section 7.4. "Pattern reduction" of the documentation created once Biocham is installed. Good-path merge restriction described in chapter 4 was not efficient enough to be integrated in the latest release. It can only be viewed in the branch feature/sepi_refinement of Biocham's repository. Bounds presented in chapter 7 are also in a separated branch called feature/sepi_bounds.

SEPI framework was implemented in an open-source software modelling platform called Biocham, available at https://lifeware.inria.fr/biocham4.

Acknowledgements

Foremost, I would like to thank Dr. François Fages, who offered me this diploma thesis subject and who supervised me during this thesis. I'm grateful to Dr. Johannes Fichte for his support during this project.

I would like to thank Dr. Sylvain Soliman, for sharing with me his expertness of Prolog.

I would like to give credit to Eva Philippe who also worked on subgraph epimorphisms.

Then the function to identify two-neighbours vertices is: for all( A in V ) { for all ( (A,B) in E ) { for all ( (B,C) in E ) { assert(neighbours(A,C)) } } } With the definition of two-neighbours, clauses can now be written.

∀a, b ∈ V such that a and b are not two-neighbours, ∀y ∈ V , cl(¬m a,y ∨ ¬m b,y ).

These clauses can easily be translated with the pseudocode: for all( A in V and B in V ) { if ( not(neighbours(A,B)) ) { for all( Y in V' ) { write(-m(a,y) -m(b,y)) } } } Strict two-neighbours restriction is much faster and simpler than good-path merge restriction. An evaluation on BioModels will be performed in next section.

Evaluation

The strict merge restriction was first evaluated on handmade models and then tested on bigger models from BioModels.

Handmade models

Table 5-3 and table 5-4 compare the number of SEPIs with the initial SEPI framework (line "Without") and with the strict merge restriction (line "Strict"). Table 5-3 presents results on combinations of two Michaelis-Menten patterns. It can be observed that strict merge restriction permits to reduce the number of unwanted reductions. With only three models for three-levels MAPK cascade in the three reduction forms, no difference is observed with or without the strict merge restriction when the number of SEPIs is recorded until 200.

With the set of 15 models corresponding to MAPK cascades with 1, 2 or 3 levels and the three different forms of reduction, merge restriction permits to eliminate the timeouts and reduce the number of reductions. With this data structure, the Michaelis-Menten pattern reduction problem would have a complexity of O(|V | 7 ): the existence of a Michaelis-Menten pattern need to be checked for each vertex of the graph.

A better solution was to change the data structure by constructing two dictionaries, ingoing_edges and outgoing_edges.

With this new data structure, the complexity for every pattern is only O(|V |), as it will be shown in the pseudocode subsection.

Pseudocode

This strategy does not need a SAT solver, all the encoding is made in Prolog C.

The Michaelis-Menten pattern will be taken as example for the implementation. Checking if a vertex is the core vertex of a Michaelis-Menten pattern has a constant cost. The core vertex of the pattern is the complex ES. To check the pattern, it is sufficient to check all ongoing and outgoing edges of isolated vertices. Is the Michaelis-Menten pattern, ES, R 1 , R -1 and R 2 are fixed vertices. The function michaelis_menten_pattern checks ongoing and outgoing edges and identifies species E, S and P that will be needed for the graph rewriting part.

michaelis_menten_pattern(+ES, -E, -S, -P):-ingoing_edges[ES] = [R1], outgoing_edges[ES] = [R_1, R2], ingoing_edges[R1] = [E, S], outgoing_edges[R1] = [ES], ingoing_edges[R_1] = [ES], outgoing_edges[R_1] = [E, S], ingoing_edges[R2] = [ES], outgoing_edges[R2] = [E, P].

To identify all patterns of a graph, it is then sufficient to review all specie vertices to check if they are the core of the pattern. for all( Vertex in V ) { if ( michaelis_menten_pattern(+Vertex, -E, -S, -P) ): rewrite_reaction("E+S=>E+P.") } If they are the core of the pattern, reactions R 1 , R -1 and R 2 can all be rewritten by a simple reaction R: E + S => E + P .

After reviewing all vertices of the initial graph, all reactions that were not in the pattern can be added to the new graph.

Bounds

Notations

• For u ∈ S i (resp. R i ) and u ∈ S (resp. R ), m i,< u,u = 1 if and only if µ i (u) < u . Left totality clauses are ensuring that each vertices of the starting graph have at least one image (or s 0 or r 0 for ⊥) in the targeted graph.

Functionality clauses are ensuring that each vertices of the starting graph have only one image through the morphism.

• For v ∈ S ∪ R , e v = 1 if and only if v is in the image of µ i (V i ), which means that vertex v exists in the glb graph G.

Covering of existing vertices clauses put together the two morphisms µ 1 and µ 2 .

Bound maximisation clauses allow to maximize the size of the search graph G, with MAX-SAT.

Graph morphism clauses ensure that an arc has an image through the morphism µ i if both its vertices have an image through µ i .

which means that the arc (u , v ) exists in the glb graph G.

Label epimorphism clauses link together the arcs through the two morphisms µ 1 and µ 2 .

• Existing set normalisation clauses force existing vertices to be the smallest in S and R .

The species (resp. reaction) vertices of the smallest initial graph (supposed to be G 1 ) are put in a total order s 1 1 < . . . < s 1 p (resp. r 1 1 < . . . < r 1 q ). B S v,v = 1 if and only if ∀u ≤ v, µ(u) ≤ v . Symmetry breaking clauses force the tuple of the image vertices of the smallest initial graph (µ(v 1 ), µ(v 2 ), . . . , µ(v n )) to be written in the smallest lexicographical order, with respect to all possible permutations of the image set. To do this, we use a precedence constraint (see [START_REF] Walsh | Symmetry breaking using value precedence[END_REF] and [1]) with the bounding functions B S (resp. B R ) such that B S : S → S is defined by B(v) = max{µ(u), u ≤ v}, or equivalently by induction : B(v 1 ) = µ 1 (v 1 ) and B(v i+1 ) = max(B(v i ), µ(v i+1 )). The condition on the lexicographical minimality of (µ(v 1 ), µ(v 2 ), . . . , µ(v n )) is equivalent to impose that B(v i+1 ) -B(v i ) ≤ 1∀i ∈ 0, n -1 .

Clauses

I Left totality.

∀s i ∈ S i , i cl( s ∈S m i s i ,s ).

Appendices

List of Figures

List of Tables