Nathann Cohen

Frédéric Havet

William Lochet

Nicolas Nisse

Subdivisions of oriented cycles in digraphs with large chromatic number

published or not. The documents may come   L'archive ouverte pluridisciplinaire

Subdivisions of oriented cycles in digraphs with large chromatic number

Nathann Cohen 1 , Frédéric Havet 2 William Lochet 2 and Nicolas Nisse 3

1 CNRS, LRI, Univ. Paris Sud, Orsay, France 2 Univ. Côte d'Azur, CNRS, Inria, I3S, France 3 Univ. Côte d'Azur, Inria, CNRS, I3S, France

Extended Abstract

The chromatic number χ(D) of a digraph D is the chromatic number of its underlying graph.

The chromatic number of a class of digraphs D, denoted by χ(D), is the smallest k such that χ(D) ≤ k for all D ∈ D, or +∞ if no such k exists. If χ(D) = +∞, we say that D has bounded chromatic number.

We are interested in the following question : which are the digraph classes D such that every digraph with sufficiently large chromatic number contains an element of D ? Let us denote by Forb(H) (resp. Forb(H)) the class of digraphs that do not contain H (resp. any element of H) as a subdigraph. The above question can be restated as follows :

Problem 1. Which are the classes of digraphs D such that χ(Forb(D)) < +∞ ?

An oriented graph is an orientation of a (simple) graph. An oriented path (resp., an oriented cycle) is said directed if all vertices have in-degree and out-degree at most 1.

Observe that if D is an orientation of a graph G and Forb(D) has bounded chromatic number, then Forb(G) has also bounded chromatic number. A classical result by Erdős implies that G must be a tree. Burr proved that every (k -1) 2 -chromatic digraph contains every oriented tree of order k and conjectured Burr [3] that it could be further improved to (2k -2)-chromatic digraphs.

For special oriented trees T , better bounds on the chromatic number of Forb(T ) are known. The most famous one, known as Gallai-Hasse-Roy-Vitaver Theorem [START_REF] Gallai | On directed paths and circuits[END_REF] states that χ(Forb(P + (k))) = k, where P + (k) is the directed path of length k (a directed path is an oriented path in which all arcs are in the same direction).

The chromatic number of the class of digraphs not containing a prescribed oriented path P on n vertices with two blocks (blocks are maximal directed subpaths) has been determined by Addario-Berry et al. [START_REF] Addario-Berry | Paths with two blocks in n-chromatic digraphs[END_REF] : Theorem 2 (Addario-Berry et al. [START_REF] Addario-Berry | Paths with two blocks in n-chromatic digraphs[END_REF]). Let P be an oriented path with two blocks on n ≥ 4 vertices, then χ(Forb(P )) = n -1.

In this paper, we are interested in the chromatic number of Forb(H) when H is an infinite family of oriented cycles. Let us denote by S-Forb(D) (resp. S-Forb(D)) the class of digraphs that contain no subdivision of D (resp. any element of D) as a subdigraph. We are particularly interested in the chromatic number of S-Forb(C), where C is a family of oriented cycles.

Let us denote by C k the directed cycle of length k. For all k, χ(S-Forb( C k )) = +∞ because transitive tournaments have no directed cycle. Let us denote by C(k, ) the oriented cycle with two blocks, one of length k and the other of length . Observe that the oriented cycles with two blocks are the subdivisions of C(1, 1). As pointed by Gyárfás and Thomassen (see [START_REF] Addario-Berry | Paths with two blocks in n-chromatic digraphs[END_REF]), there are acyclic oriented graphs with arbitrarily large chromatic number and no oriented cycles with two blocks. Therefore χ(S-Forb(C(k, ))) = +∞. We first generalise this result to every oriented cycle. Theorem 3. For any oriented cycle C, χ(S-Forb(C)) = +∞.

In fact, we show the following stronger theorem. We need a construction due to Erdős and Lovász [START_REF] Erdős | Problems and results on 3-chromatic hypergraphs and some related questions[END_REF] of hypergraphs with high girth and large chromatic number. Theorem 5. [5, Theorem 1'] For k, g, c ∈ N, there exists a k-uniform hypergraph with girth larger than g and weak chromatic number larger than c.

We assume g is being fixed, the following construction allow us to find D c+1 from D c . Let p be the number of proper c-colourings of D c , and let those colourings be denoted by col 1 c , ..., col p c . By Theorem 5 there exists a c × p-uniform hypergraph H with weak chromatic number > p and girth > g/2. Let X = {x 1 , . . . , x n } be the ground set of H.

We construct D c+1 from n disjoint copies D 1 c , ..., D n c of D c as follows. For each hyperedge S ∈ H, we do the following :

• We partition S into p sets S 1 , . . . , S p of cardinality c.

• For each set S i = {x k1 , . . . , x kc }, we choose vertices v k1 ∈ D k1 c , . . . , v kc ∈ D kc c such that col i c (v k1 ) = 1, . . . , col i c (v kc ) = c
, and add a new vertex w S,i with v k1 , . . . , v kc as in-neighbours.

On the other hand, considering strongly connected (strong for short) digraphs may lead to dramatically different result. An example is provided by the following celebrated result due to Bondy [START_REF] Bondy | Disconnected orientations and a conjecture of Las Vergnas[END_REF], which can be rephrased as follows when denoting the class of strong digraphs by S.

Theorem 6 (Bondy [START_REF] Bondy | Disconnected orientations and a conjecture of Las Vergnas[END_REF]).

χ(S-Forb( C k ) ∩ S) = k -1.
Inspired by this theorem, Addario-Berry et al. [START_REF] Addario-Berry | Paths with two blocks in n-chromatic digraphs[END_REF] posed the following problem. A consequence of the previous lemma is that, if we partition the arc set of D into set A 1 • • • A k , then bounding the chromatic number of all digraphs induced by the A i implies that D has bounded chromatic number.

Proof. Let D be a strong digraph without any copy of C(k, ), we exhibit a colouring of D using a bounded number of colours. The proof heavily relies on the technique of levelling. Let u be a vertex of D. The level of a vertex x, noted lvl(x) is the length of the shortest dipath from u to x. L(i) is the set of vertices at level i.

Since D is strongly connected, it has an out-generator u. Let T be a BFS-tree with root u. We define the following sets of arcs.

A 0 = {xy ∈ A(D) | lvl(x) = lvl(y)}; A 1 = {xy ∈ A(D) | 0 < | lvl(x) -lvl(y)| < k + -3; A = {xy ∈ A(D) | lvl(x) -lvl(y) ≥ k + -3}.
Since k + -3 > 0 and there is no arc xy with lvl(y) > lvl(x) + 1, (A 0 , A 1 , A ) is a partition of A(D). Observe moreover that A(T ) ⊆ A 1 . We further partition A into two sets A 2 and A 3 , where A 2 = {xy ∈ A | y is an ancestor of x in T } and

A 3 = A \ A 2 . Then (A 0 , A 1 , A 2 , A 3 ) is a partition of A(D). Let D j = (V (D), A j ) for all j ∈ {0, 1, 2, 3}. Claim 10. χ(D 0 ) ≤ k + -2. Proof. Observe that D 0 is the disjoint union of the D[L i ] where L i = {v | dist D (u, v) = i}.
Therefore it suffices to prove that χ(D[L i ]) ≤ k + -2 for all non-negative integer i.

L 0 = {u} so the result holds trivially for i = 0.

Assume now i ≥ 1. Suppose for a contradiction χ(D[L i ]) ≥ k + -1. Since k ≥ 3, by Theorem 2, D[L i ] contains a copy Q of P + (k -1, - 1
), the path on two blocks of length k -1 and -1 with one vertex of indegree 2. Let v 1 and v 2 be the initial and terminal vertices of Q, and let x be the least common ancestor of v 1 and v 2 . By definition, for j ∈ {1, 2}, there exists a dipath P j from x to v j in T . By definition of least common ancestor, V (P 1 ) ∩ V (P 2 ) = {x}, V (P j ) ∩ L i = {v j }, j = 1, 2, and both P 1 and P 2 have length at least 1. Consequently,

P 1 ∪ P 2 ∪ Q is a subdivision of C(k, ), a contradiction. Claim 11. χ(D 1 ) ≤ k + -3.
Proof. Let φ 1 be the colouring of D 1 defined by φ 1 (x) = lvl(x) (mod k + -3). By definition of D 1 , this is clearly a proper colouring of D 1 .

The following two claims are more complicated, we refer the reader to [START_REF] Cohen | Subdivisions of oriented cycles in digraphs with large chromatic number[END_REF] for the complete proofs. Proof. The general idea is the same as in the proof of Theorem 8.

Suppose that D admits an out-generator u and let T be an BFS-tree with root u. We partition A(D) into three sets according to the levels of u.

A 0 = {(x, y) ∈ A(D) | lvl(x) = lvl(y)}; A 1 = {(x, y) ∈ A(D) | | lvl(x) -lvl(y)| = 1}; A 2 = {(x, y) ∈ A(D) | lvl(y) ≤ lvl(x) -2}. For i = 0, 1, 2, let D i = (V (D), A i ).
Claim 15. χ(D 0 ) ≤ 3.

Proof. Suppose for a contradiction that χ(D) ≥ 4. By Theorem 2, it contains a P -(1, 1) (y 1 , y, y 2 ), that is (y, y 1 ) and (y, y 2 ) are in A(D 0 ). Let x be the least common ancestor of y 1 and y 2 in T . The union of T [x, y 1 ], (y, y 1 ), (y, y 2 ), and T [x, y 2 ] is a subdivision of Ĉ4 , a contradiction.

Claim 16. χ(D 1 ) ≤ 2.

Proof. Since the arc are between consecutive levels, then the colouring φ 1 defined by φ 1 (x) = lvl(x) mod 2 is a proper 2-colouring of D 1 .

Let y ∈ V i we denote by N (y) the out-degree of y in 0≤j≤i-1 V j . Let D = (V, A ) with A = ∪ x∈V {(x, y), y ∈ N (x)} and D x = (V, A x ) where A x is the set of arc inside the level and from V i to V i+1 for all i. Note that A = A ∪ A x and

Theorem 4 .

 4 For any positive integers b, c, there exists an acyclic digraph D c with χ(D c ) ≥ c in which all oriented cycles have more than b blocks.

Problem 7 .

 7 Let k and be two positive integers then χ(S-Forb(C(k, ) ∩ S)) < k + l.We give evidence for this problem by showing the following weaker statement. Theorem 8. Let k and be two positive integers such that k ≥ max{ , 3}, and let D be a digraph in S-Forb(C(k, )) ∩ S. Then, χ(D) ≤ (k + -2)(k + -3)(2 + 2)(k + + 1). We need the following lemma. The union of two digraphs D 1 and D 2 is the digraph D 1 ∪D 2 with vertex set V (D 1 )∪V (D 2 ) and arc set A(D 1 ) ∪ A(D 2 ). Lemma 9. Let D 1 and D 2 be two digraphs. χ(D 1 ∪ D 2 ) ≤ χ(D 1 ) × χ(D 2 ).

  Claim 12. χ(D 2 ) ≤ 2 + 2. Claim 13. χ(D 3 ) ≤ k + + 1. Claims 10, 11, 12, and 13, together with Lemma 9 yield the result. More generally, one may wonder what happens for other oriented cycles. Our next result generalises Theorem 8 for Ĉ4 the cycle with 4 blocks. Theorem 14. Let D be a digraph in S-Forb( Ĉ4 ). If D admits an out-generator, then χ(D) ≤ 24.

Claim 17. χ(D 2 ) ≤ 4.

Proof. We refer to [4] for the proof of this statement.

Claims 15, 16, 17, and Lemma 9 implies χ(D) ≤ 24.


Frédéric Cazals

Dorian Mazauric

Optimal transportation problems with connectivity constraints

Keywords: Optimal transportation, connectivity constraints, NP-completeness, graph algorithms, clustering analysis, molecular simulation Transport optimal, contraintes de connectivité, NP-complétude, graphes et algorithmes, clustering, simulation moléculaire

The earth mover distance (EMD) or the Mallows distance are example optimal transportation (OT) problems reducing to linear programs. In this work, we study a generalization of these problems when the supply and demand nodes are the vertices of two graphs called the supply and the demand graphs. The novel problems embed connectivity constraints in the transport plans computed, using a Lipschitz-like condition involving distances between certain subgraphs of the supply graph and certain subgraphs of the demand graph. More precisely, we make three contributions. First, we formally introduce two optimal transportation problems generalizing EMD, namely Minimum-cost under flow, transport size, and connectivity constraints problem (problem EMD-FCC) and Maximum-flow under cost, transport size, and connectivity constraints problem (problem EMD-CCC). We prove that problems EMD-CCC and EMD-FCC are NP-complete, and that EMD-FCC is hard to approximate within any given constant. Second, we develop a greedy heuristic algorithm returning admissible solutions, of time complexity O(n 3 m 2 ) with n and m the numbers of vertices of the supply and demand graphs, respectively. Third, on the experimental side, we apply our novel OT algorithms for two applications, namely the comparison of clusterings, and the analysis of so-called potential energy landscapes in molecular science. These experiments show that optimizing the transport plan and respecting connectivity constraint can be competing objectives. Implementations of our algorithms are available in the Structural Bioinformatics Library at http://sbl.inria.fr.

Introduction

Optimal transportation (OT) problems have a long standing history in mathematics and computer science, originating with the works of Monge on earth moving ( la théorie des déblais et des remblais ) [START_REF] Monge | Mémoire sur la théorie des déblais et des remblais[END_REF]. Such problems were later rephrased in terms of Riemannian geometry and measure theory [START_REF] Villani | Topics in optimal transportation[END_REF], one key concept being the distance between two distributions, namely is the minimal amount of work that must be performed to transform one distribution into the other by moving distribution mass around. Various applications were developed across all sciences, one of the early ones in computer science being the earth mover distance (EMD), used to compare two images using their color histograms [START_REF] Rubner | The earth mover's distance as a metric for image retrieval[END_REF]. The EMD is also related to the Mallows distance used in statistics [START_REF] Levina | The earth mover's distance is the mallows distance: Some insights from statistics[END_REF]. Both distances are of special interest since they can be phrased as linear programs, and hence, are amenable to efficient algorithms. Beyond linear programming, OT problems are generally hard to solve, which motivated the development of stochastic optimization techniques [START_REF] Genevay | Stochastic optimization for large-scale optimal transport[END_REF], or the exploitation of specific (geometry) properties of the functional studied [START_REF] Solomon | Convolutional wasserstein distances: Efficient optimal transportation on geometric domains[END_REF].

In this work, we explore a new dimension of OT problems, namely when the supply and demand nodes are the vertices of two graphs called the supply and the demand graphs. In the presence of connectivity information provided by these graphs, it is natural to expect transport plans to comply with this connectivity information, while still minimizing the transport cost. For example, a transport plan may be termed valid provided that any connected component of the supply graph induces (via edges carrying flow) a connected component of the demand graph-a condition called strict connectivity constraints in the sequel. Naturally, depending on the nature of the constraints, one may end up with more costly transport plans, or worse, may face unfeasible problems. This paper precisely addresses these problems, the notion of connectivity constraints being modeled by a Lipschitz-like condition involving distances between certain subgraphs of the supply graph and certain subgraphs of the demand graph.

Practically, we apply our novel OT algorithms for two applications, namely the comparison of clusters obtained by mode seeking algorithms [START_REF] Cheng | Mean shift, mode seeking, and clustering[END_REF][START_REF] Chazal | Persistence-based clustering in riemannian manifolds[END_REF], and also the analysis of potential energy landscapes (PEL) of molecular systems [START_REF] Wales | Energy Landscapes[END_REF][START_REF] Cazals | Conformational ensembles and sampled energy landscapes: Analysis and comparison[END_REF][START_REF] Carr | Energy landscapes and persistent minima[END_REF].

Paper overview and contributions. Section 2 introduces OT problems with connectivity constraints. In Section 3, we prove that the problems are NP-complete and not in APX for very simple classes of instances, and we prove stronger NP-completeness results for strict connectivity constraints. In Section 4, we develop exact, approximation and heuristic polynomial time algorithms, and we investigate the case of strict connectivity constraints. Finally, Section 5 present experiments on the comparison of clusterings and on the analysis of molecular potential energy landscapes.

Due to the lack of space, all proofs are presented in appendix.

Problem formulation and models

We first define notations, and proceed with various constraints (connectivity constraints, transport size), from which our optimal transportation problems are defined.

Pre-requisites

Notations. Consider two graphs: a supply graph G = (V, E) and a demand graph G = (V , E ). The set V = {v 1 , . . . , v |I| } represents the supply nodes. We denote by I the set of indices of the supply nodes. The value X vi ≥ 0 represents the volume of supply of node v i for Inria all i ∈ I. The set V = {v 1 , . . . , v |J | } represents the demand nodes. We denote by J the set of indices of the demand nodes. The value Y v j ≥ 0 represents the volume of demand of node v j for all j ∈ J . Let B = (V ∪ V , V × V ) be the complete bipartite graph between the supply and the demand nodes. The real values c vi,v j ≥ 0 represent the linear cost of sending a unit of flow from node v i to node v j for all i ∈ I, j ∈ J . The variable f vi,v j represents the volume of flow sent by node v i to v j for all i ∈ I, j ∈ J . When there is no ambiguity, we abuse the notation writing X i , Y j , c i,j , and f i,j instead of X vi , Y v j , c vi,v j , and f vi,v j , respectively. For all i ∈ I, j ∈ J , the cost of sending a volume of flow f i,j through edge {v i , v j } ∈ E(B) is f i,j c i,j . Given the flows f i,j for all edges of B, the total flow is i∈I j∈J f i,j and the total cost is i∈I j∈J f i,j c i,j .

Minimum-cost flow problem. The classical minimum-cost flow problem (EMD), or transportation problem [START_REF] Dantzig | Application of the simplex method to a transportation problem[END_REF], consists in determining a minimum-cost flow satisfying the demands and respecting the supply constraints. This problem is polynomial since it reduces to solving the following linear program (LP):

               Minimize i∈I j∈J f i,j c i,j subject to i∈I f i,j ≤ Y j ∀j ∈ J , j∈J f i,j ≤ X i ∀i ∈ I, i∈I j∈J f i,j = min( i∈I X i , j∈J Y j ), f i,j ≥ 0 ∀i ∈ I, ∀j ∈ J . (1)
The first line is the objective function, the third line represents the demand constraints, the fourth line describes the supply constraints, the fifth line states that the total amount of flow equals the minimum between the total volume of supplies and the total volume of demands, and the last line guarantees that flows are positive. Note that if i∈I X i ≥ j∈J Y j , then the fifth line can be removed and the inequality constraints of the third line become equality constraints. The number of edges that support flow is at most the total number of nodes (of G and G ) minus one (e.g. see [START_REF] Hillier | Introduction to mathematical programming[END_REF]). Say otherwise, i∈I j∈J |fi,j >0 1 ≤ |I| + |J | -1. Given an optimal solution f for EMD, we define the total number of edges M EMD , the total flow F EMD , the total cost C EMD , and the ratio d EMD (a.k.a. the earth mover distance [START_REF] Rubner | The earth mover's distance as a metric for image retrieval[END_REF]):

M EMD = i∈I,j∈J 1 fi,j >0 , F EMD = i∈I,j∈J f i,j , C EMD = i∈I,j∈J f i,j c i,j , d EMD = C EMD /F EMD . (2)
Connectivity constraints. When the supply nodes (and likewise the demand nodes) are endowed with a graph structure, the distances between different sub-graphs of G must be correlated with the distances between the different sub-graphs of G that receive flow from the former sub-graphs.

To formalize this idea, we introduce the following notations (Fig. 1). Let d G (v, v ) be the number of edges minus one of a minimum shortest path between

v ∈ V and v ∈ V in G . If v = v , we set d G (v, v ) = 0.
In the following, we denote by cc(G ) the set of maximal connected components of G . Let H = {H 1 , . . . , H t } be any set of t ≥ 1 disjoint sub-graphs of G . Note that H i is not necessarily a connected sub-graph. We define

d G (H i , H j ) = min v∈V (H i ),v ∈V (H j ) d G (v, v ). Note that if v and v are not connected in G, then d G (v, v ) = ∞. We define d G (H ) = max i,j∈{1,...,t} d G (H i , H j ). If |H | = 1, then d G (H ) = 0. Furthermore, we define d G (cc(H i )) as follows. Let cc(H i ) = {cc 1 (H i ), . . . , cc s (H i )} be the set of the s ≥ 1 maximal connected Figure 1 Distances between subgraphs. Example of graph G with H = {H 1 , H 2 }: d G (H ) = d G (H 1 , H 2 ) = 1, d G (cc(H 1 )) = 0 because H 1 is connected, and d G (cc(H 2 )) = 3. H'1 H'2 components of H i . We define d G (cc(H i )) = max i,j∈{1,...,s} d G (cc i (H i ), cc j (H j )). As previously said, if |cc(H i )| = 1, that is H i is a connected graph, then d G (cc(H i )) = 0.
Distances within G and G shall be processed by warping functions: • g is non-decreasing: the connectivity constraints are not stronger when the distance d increases.

• g(x) ≥ x for all x ≥ 0: if any two nodes are at distance d in G, then we cannot constrain the distance between the two sub-graphs that receive flow from the former two nodes, to be less than d.

• g(∞) = ∞: if two nodes u, v ∈ V are not in a same maximal connected component in G, then any node that receives flow from u and any node that receives flow from v, can be in two different maximal connected components in G .

We finally arrive at connectivity constraints associated with a distance warping function:

Definition. 2. (Connectivity constraints) Consider a distance warping function g. The connectivity constraints are satisfied for a solution f if and only if for every i, i ∈ I, then

d G (H i , H i ) ≤ g(d G (v i , v i )) and d G (cc(H i )) ≤ g(0), (3)
where H i (H i , respectively) is the sub-graph induced by the nodes that receive flow from v i (v i , respectively).

This definition calls for two comments. First, the two terms of Eq. ( 3) respectively define constraints associated with pairs of vertices and single vertices of the supply graph. In Lemma 1, we shall give an alternative definition of these constraints. Second, the conditions imposed with function g are analogous to Lipschitz conditions on distances. Indeed,

d G (H i , H i ) ≤ g(d G (v i , v i ))
means that the distance between nodes that receive flow from v i and nodes that receive flow from v i must be upper-bounded by a function of the distance between v i and v i in G (constraints for pairs of nodes of G). Furthermore, d G (cc(H i )) ≤ g(0) means that the maximal connected components induced by the set of nodes that receive flow from v i must be at most at distance g(0) each other (constraints for single nodes of G).

As a particular case, we may ask that the two vertices defining an edge from G, and a vertex from G export flow to connected subgraphs of G . Using the two constraints from Eq. (3), we define: Inria Definition. 3. (Strict connectivity constraints) Strict connectivity constraints are defined by a distance warping function g such that g(0) = 0 and g(x) = ∞ for all x ≥ 1.

Transport size constraint. Recall that for EMD, the number of edges supporting flow is at most the total number of nodes minus one. Connectivity constraints may be such that a superlinear number of edges is needed. But since we do not know this number a priori, we impose an upper bound M on it. Formally, given any integer M such that 0 ≤ M ≤ |E(B)|, the transport size constraint is satisfied for f if and only if i∈I j∈J |fi,j >0 1 ≤ M . Practically, our implemented algorithms do not take this upper bound as input, but we carefully analyze the number of edges carrying flow and compare it against the number of nodes.

Problems EMD-FCC and EMD-CCC

Minimum-cost under flow, transport size, and connectivity constraints problem: problem EMD-FCC. This problem (EMD-FCC) consists in computing the smallest total cost when a given amount of total flow must be supported and respecting the transport size, and the connectivity constraints. Formally, given F and M , 0 ≤ F ≤ min( i∈I X i , j∈J Y j ), 0 ≤ M ≤ |E(B)|, EMD-FCC can be written as follows:

                           Minimize i∈I j∈J f i,j c i,j subject to i∈I f i,j ≤ Y j ∀j ∈ J , j∈J f i,j ≤ X i ∀i ∈ I, f i,j ≥ 0 ∀i ∈ I, ∀j ∈ J , i∈I j∈J f i,j ≥ F, i∈I j∈J |fi,j >0 1 ≤ M, d G (H i , H i ) ≤ g(d G (v i , v i )) ∀i, i ∈ I, d G (cc(H i )) ≤ g(0) ∀i ∈ I. (4)
Line 6 represents the flow constraint, Line 7 is the transport size constraint, Lines 8 and 9 describe the connectivity constraints. Fig. 2 describes the solutions of a simple instance of EMD-FCC. Given an optimal solution f for EMD-FCC, we introduce the total number of edges M EMD-FCC , the total flow F EMD-FCC , the total cost C EMD-FCC , and the ratio d EMD-FCC :

M EMD-FCC = i∈I,j∈J 1 fi,j >0 , F EMD-FCC = i∈I,j∈J f i,j , C EMD-FCC = i∈I,j∈J f i,j c i,j , d EMD-FCC = C EMD-FCC F EMD-FCC . (5)
Maximum-flow under cost, transport size, and connectivity constraints problem: problem EMD-CCC. This problem (EMD-CCC) aims at computing the largest volume of flow that can be supported respecting the connectivity constraints, the transport size constraint and such that the total cost is less than a given bound C. We define the following upper bound C max = j∈J Y j max i∈I,j∈J c i,j for the maximum total cost of any admissible flow. Given C

                           Maximize i∈I j∈J f i,j subject to i∈I f i,j ≤ Y j ∀j ∈ J , j∈J f i,j ≤ X i ∀i ∈ I, f i,j ≥ 0 ∀i ∈ I, ∀j ∈ J , i∈I j∈J f i,j c i,j ≤ C, i∈I j∈J |fi,j >0 1 ≤ M, d G (H i , H i ) ≤ g(d G (v i , v i )) ∀i, i ∈ I, d G (cc(H i )) ≤ g(0) ∀i ∈ I. (6)
Given an optimal solution f for EMD-CCC, we introduce the total number of edges M EMD-CCC , the total flow F EMD-CCC , the total cost C EMD-CCC and the ratio d EMD-CCC :

M EMD-CCC = i∈I j∈J |fi,j >0 1 ≤ M, F EMD-CCC = i∈I j∈J f i,j , C EMD-CCC = i∈I j∈J f i,j c i,j , d EMD-CCC = C EMD-CCC F EMD-CCC . (7)
3 The problems are very difficult to solve

We first prove that the connectivity constraints studied in this article are equivalent to constraints with general size and general number of sub-graphs of the source graph (Section 3.1). We then show the NP-completeness of EMD-CCC and EMD-FCC (Section 3.2). Furthermore, we prove a stronger hardness result with strict connectivity constraints (Section 3.3). We finally show that EMD-FCC is hard to approximate within any given constant (Section 3.4).

Inria

The connectivity constraints are equivalent to constraints with general size and general number of sub-graphs

We prove in Lemma 1 that the connectivity constraints studied in this article are equivalent to constraints with general size (taking sub-graphs of G of any size) and with general number of sub-graphs (taking a set {H 1 , H 2 , . . .} of sub-graphs of G of any size).

Lemma. 1. Given any flow f , the connectivity constraints are satisfied for f if and only

if d G (H ) ≤ g(d G (H)) and d G (cc(H i )) ≤ g(d G (cc(H i ))) for all i, 1 ≤ i ≤ t, where H = {H 1 , . . . , H t } is any set of t ≥ 1 disjoint sub-graphs of G, H i is the sub-graph of G induced by
the set of nodes that receive flow from at least one node of H i , and

H = {H 1 , . . . , H t }.
Note in particular that under strict connectivity constraints, the previous lemma implies that the set of nodes of any connected sub-graph of G sends flow to a set of nodes of G that induces a connected sub-graph of G .

NP-completeness of EMD-FCC and EMD-CCC

We prove that EMD-FCC and EMD-CCC are NP-complete even for simple classes of instances. In our following reductions, we use the strongly NP-complete problem 3-Partition [START_REF] Garey | Computers and Intractability: A Guide to the Theory of NP-Completeness[END_REF]. Let m ≥ 1 be any integer. Given a set S = {n 1 , n 2 , . . . , n 3m } of 3m positive integers, 3-Partition problem consists in deciding if S can be partitioned into m subsets such that the sum of the numbers in each subset is equal.

Lemma. 2. For any distance warping function g, the decision version of EMD-FCC is NPcomplete even if:

• the demand graph G is a complete graph;

• and all the volumes of demands are equal (Y j = Y j for all j, j ∈ J );

• and all the unitary costs are equal to one (c i,j = 1 for all i ∈ I, j ∈ J );

• and the volumes of demands and supplies are equal ( i∈I X i = j∈J Y j ).

We deduce in Corollary 1 that the decision version of EMD-CCC is NP-complete. Indeed, given a maximum cost C and a maximum number of edges M , the problem of deciding if there exists a flow f larger than a given F and satisfying all the constraints, is equivalent to the problem of deciding if there is an admissible solution for EMD-FCC with C as input.

Corollary. 1. The decision version of EMD-CCC is NP-complete.

The connectivity constraints are not directly considered in our reduction since G is a complete graph. The choice of M = 3 4 (|I| + |J |) is therefore the main key of the proof of Lemma 2 (and Corollary 1); this number departs from the linear number of edges (at most |I| + |J | -1) for problem EMD [START_REF] Hillier | Introduction to mathematical programming[END_REF]. Allowing a quadratic number of edges is covered in the next section.

Stronger NP-completeness result with strict connectivity constraints

Under strict connectivity constraints, we prove that the decision version of EMD-FCC and EMD-CCC is also NP-complete even if the upper-bound M on the number of edges that can support flow is quadratic in the number of nodes. We first establish the following result:

Lemma. 3. Consider the strict connectivity constraints. There exists an instance of the decision version of EMD-FCC such that there is an admissible solution if and only if

M ≥ |I|(|I| + 1) with |J | = 2|I|.
From which one deduces:

Lemma. 4. Under strict connectivity constraints, the decision version of EMD-FCC and EMD-CCC is NP-complete even if M ∈ Θ(|V ∪ V | 2 ).

Hardness of approximation of EMD-FCC

We now prove in Lemma 5 that EMD-FCC is hard to approximate. More precisely, we show that for any constant k ≥ 1, there is no polynomial time algorithm for EMD-FCC, unless P = NP.

Lemma. 5. EMD-FCC is not in APX even if:

• all the volumes of demands are equal (Y j = Y j for all j, j ∈ J );

• and there are only two different unitary costs for edges of the bipartite graph B (c i,j ∈ {1, K} for all i ∈ I, j ∈ J , where K > 1).

Efficient polynomial time algorithms

In this section, we first describe two classes of instances that admit exact polynomial time algorithms (Section 4.1). We then prove the existence of a Polynomial Time Approximation Scheme for some classes of instances (Section 4.2). We also develop efficient polynomial time greedy algorithms (Section 4.3), and optimize its complexity for strict connectivity constraints (Section 4.4).

Exact polynomial time algorithms for some classes of instances

We prove in Lemma 6 that if the demand graph is complete or if the demand graph is connected and the distance between any two nodes is upper-bounded by the function g (for any value), then there is an exact polynomial time algorithms for EMD-FCC and EMD-CCC.

Lemma. 6. There is a polynomial time algorithm for EMD-FCC and EMD-CCC if

• M = |E(B)|
and G be a complete graph;

• or M = |E(B)|, G is any connected graph, and for every u , v ∈ V , d G (u , v ) ≤ g(x) for all x ≥ 0.

Polynomial Time Approximation Scheme

We prove in Lemma 7 a Polynomial Time Approximation Scheme (PTAS) for EMD-FCC when G is a connected graph and M = |E(B)|. To do that, we first add a volume of flow ε for all edges of the complete bipartite graph B, and then obtain an auxiliary instance (in which we update the supply and the demand volumes for all nodes). By construction, for any distance warping function g, the connectivity constraints are satisfied, and so EMD-FCC is equivalent to EMD for this auxiliary instance, which gets solved by the linear program described in Eq. ( 1). Thus, we get a PTAS for EMD-FCC choosing ε function of the desired approximation ratio. An interesting problem is to determine the minimum number of edges to add in an optimal solution for EMD in order to get an admissible solution for EMD-FCC.

Inria Lemma. 7. Consider any distance warping function g. Let M = |E(B)| and G be any connected graph. Then, for any ε > 0, there is a polynomial time (1 + ε)-approximation algorithm for EMD-FCC.

In Corollary 2, we deduce the same result for EMD-CCC.

Corollary. 2. Consider any distance warping function g. Let M = |E(B)| and G be any connected graph. Then, for any ε > 0, there is a polynomial time (1+ε)-approximation algorithm for EMD-CCC.

Efficient polynomial time algorithms for EMD-CCC

Greedy algorithm. We describe here a polynomial time greedy algorithm for EMD-CCC (Algorithm 1). The inputs are the graphs G = (V, E) and G = (V , E ), a cost upper bound C, a maximum number M of edges that can support flow, and a distance warping function g for the connectivity constraints. Algorithm 1 greedily selects edges of the bipartite graph that can support flow such that the total cost is upper bounded by C, and without violating the transport size and the connectivity constraints. After such a selection, the set of candidate edges for the next step of selection is updated with respect to the connectivity constraints. Algorithm 1 returns the solution f , the total flow F f , the total cost C f , and the total number of edges M f .

Algorithm 1 Greedy algorithm for EMD-CCC. and∃(i, j) such that b i,j = 1, x i > 0, and y j > 0 do 5:

Require: G = (V, E), G = (V , E ), C, M , g. 1: F f := 0; C f := 0; M f := 0; 2: for all i ∈ I, j ∈ J do 3: f i,j := 0; x i := X i ; y j := Y j ; b i,j := 1 4: while C f < C, M F ≤ M -1,
(i t , j t ) = arg min i∈I,xi>0,j∈J ,yj >0 c i,j b i,j ; z := min( C-C f ci t ,j t
, min(x it , y jt ))

6:

F f := F f + z; C f := C f + z.c it,jt ; M F := M F + 1; 7:
f it,jt := f it,jt + z; x it := x it -z; y jt := y jt -z;

8:

Update of b i,j for all i ∈ I, j ∈ J

9: return f , F f , C f , M f .
Variables of Algorithm 1. For all i ∈ I, x i ≥ 0 represents the current volume of supply of node v i ∈ V , and so X i -x i is the current amount of flow sent by v i . For all j ∈ J , y j ≥ 0 represents the current volume of demand of node v j ∈ V , and so Y j -y j is the current amount of flow received by v j . For all v i ∈ V, v j ∈ V , the variable b vi,v j is used to encode if the edge {v i , v j } ∈ E(B) can support flow in respect with the constraints. When there is no ambiguity, we abuse the notation writing b i,j instead of b vi,v j . In other words, b i,j = 1 if the edge {v i , v j } is an edge candidate (b i,j = 0 otherwise). The variable f represents the current solution. In other words, for all i ∈ I, j ∈ J , f i,j represents the current flow sent from v i ∈ V to v j ∈ V . Furthermore, F f is the total volume of the current flow, C f represents the total cost of the current flow, M f is the current number of edges of B that support flow, and F f is the total volume of the current flow. Initially, F f = 0, C f = 0, M f = 0, f i,j = 0, x i = X i , y j = Y j , and b i,j = 1 for all i ∈ I, j ∈ J .

Core of Algorithm 1. While the current cost C f is less than the given cost upper bound C, while the current number M f of edges of B that support flow is strictly less than the given upper bound M , and while there exists an edge candidate {v i , v j } such that x i , y j > 0 (that is such that a positive flow can be supported by {v i , v j } ∈ E(B)), then an edge {v it , v jt } ∈ E(B) is selected (Line 5). Then, the maximum amount of flow z that can be supported by the edge {v it , v jt } is computed. Line 6 updates the values of C f , M f , and F f . Line 7 updates the values of f it,jt , x it , and y jt . Line 8 updates the boolean function b for all i ∈ I, j ∈ J . Algorithm 1 finally returns f , F f , C f , and M f (Line 9). Time complexity of Algorithm 1. The time complexity of our greedy algorithm is a function of the time complexity of updating b i,j , and so function of the connectivity constraints. We formalize that in Lemma 8. Require: G = (V, E), G = (V , E ), (i t , j t ), b, f , x. Ensure: Binary values b i,j for all i ∈ I, j ∈ J .

1: for all i ∈ I do

2:
for all j ∈ J do for all j ∈ J do for all j ∈ J do The following directly follows from the definition of the connectivity constraints: Lemma. 9. Algorithm 2 updates b i,j for all i ∈ I, j ∈ J for general connectivity constraints.

6: if d G (cc(G [V (H i ) ∪ {v j }])) > g(0)
11: if d G (G [V (H i1 ) ∪ {v j }], H i2 ) > g(d G (v i1 , v i2 )) then 12: b i1,j := 0 13: if d G (H i1 , G [V (H i2 ) ∪ {v j }]) > g(d G (v i1 , v i2 )) then
Lemma. 10. The time complexity of Algorithm 2 is O(max(|V | 3 , |V | 2 |V | 3 log(|V |))).
From Lemma 8 and Lemma 10, we deduce in Corollary 3 the time complexity of Algorithm 1 when using Algorithm 2 for the update of b. Remark 1. Selecting a single edge per step in Algorithm 1 calls for two comments. First, the selection of a best set of k > 1 edges at each, does not necessarily guarantee a better solution at Inria the end of the algorithm. Indeed, we can easily construct classes of instances for which such an algorithm (that selects several edges per step) is less efficient than the algorithm presented before (in terms of volume of flow and cost). Second, such a multiple edges selection would increase the time complexity of the update of the boolean function b and so of the greedy algorithm. For these two reasons, we decide to focus on algorithms that select single edge per step.

Iterative algorithm. The maximum cost is an input of Algorithm 1. Since we do not know, a priori, the cost of interesting flow solutions, we must call Algorithm 1 with different input costs, each call generating different statistics as defined by Eq. [START_REF] Garey | Computers and Intractability: A Guide to the Theory of NP-Completeness[END_REF]. To this end, given a cost range [0, C max ], we now describe the iterative algorithm Algorithm 3 which computes different flow solutions for different total costs, by iteratively calling Algorithm 1. Note that practically, C max can be set to the maximum distance between a supply and a demand node, times the total supply. The following lemma proves that Algorithm 4 updates the set of candidate edges that can support flow under strict connectivity constraints: Lemma. 11. Algorithm 4 updates b i,j for all i ∈ I, j ∈ J for the strict connectivity constraints.

Algorithm 3 Iterative algorithm for EMD-CCC. Require: G = (V, E), G = (V , E ), C inf , C sup , g.
The corresponding complexity satisfies:

Lemma. 12. The time complexity of Algorithm 4 is O(|V | + ∆(G) 2 ∆(G )).
From Lemma 8 and Lemma 12, we deduce in Corollary 4 the time complexity of Algorithm 1 when using Algorithm 4 for the update of b with strict connectivity constraints.

(max(|V | 2 |V | 2 , ∆(G) 2 ∆(G )|V ||V |)).
This result shows the decrease of the time complexity when we consider the strict connectivity constraints. Indeed, for general connectivity constraints, recall that the worst case time-complexity of Require: G = (V, E), G = (V , E ), (i t , j t ), b, f , x. Ensure: Binary values b i,j for all i ∈ I, j ∈ J .

1: if X it -x it -f it,jt = 0 then 2:
for all j such that v j ∈ N G (v jt ) do b it,j := 1

3:

for all j such that v j / ∈ N G (v jt ) do b it,j := 0 4:

for all i such that v i ∈ N G (v it ) do 5: if X i -x i = 0 then 6:
for all j such that

v j / ∈ N G [v jt ] do b i,j := 0 7: else 8:
for all j such that v j ∈ N G (v jt ) do b it,j := 1 9:

for all i such that v i ∈ N G (v it ) do 10: if X i -x i = 0 then 11: for all j such that v j ∈ N G [v jt ] do 12: if b i,j = 0 then 13: b i,j := 1 14: for all k such that v k ∈ N G (v i ) do 15: if X k -x k > 0 and b k,j = 0 do b i,j := 0 16: return b.
Note that the maximum degrees of G and G may be linear in the number of nodes. In that case, the time complexity is O(|V | 3 |V | 2 ) but, in general, the time complexity is better. For instance, Corollary 5 proves a better complexity if G and G have bounded degrees. More precisely, the complexity only depends on Algorithm 1.

Corollary. 5. The worst case time-complexity of

Experiments

We benchmark our implementations Alg-EMD-CCC-G (Algorithm 1 solving EMD-CCC) and Alg-EMD-LP (solving EMD) to compare clusterings and analyze molecular data.

Implementations

The Structural Bioinformatics Library is a generic C++/python cross-platform software library targeting complex problems in structural bioinformatics [START_REF] Cazals | The Structural Bioinformatics Library: modeling in biomolecular science and beyond[END_REF]. It combines low-level generic (C++ template based) implementations of various algorithms-in a spirit analogous to CGAL-defining SBL-CORE, and their instantiations to solve specific biophysical problems-defining SBL-APPLICATIONS. Low level generic C++ implementations of our algorithms are available in the following package from SBL-CORE (http://sbl.inria.fr/doc/group__Earth__mover__distance-package. html); instantiations of these methods for PEL are available from the following package from SBL-APPLICATIONS http://sbl.inria.fr/doc/Energy_landscape_comparison-user-manual. html.

Inria

Comparing clusterings

Setup. Clustering is key in data analysis, yet, the variety of options (algorithms and their parameters) makes comparing clusterings a key endeavor. As an illustration, we use clusterings computed by mode-seeking methods [START_REF] Chazal | Persistence-based clustering in riemannian manifolds[END_REF]. These methods identify clusters from the catchment basins of (persistent) local maxima of an estimated density. They inherently embed the notion of persistence for clusters, providing insights on the number of clusters-as opposed to requiring a pre-defined number of clusters. Finally, edges of the clustering graph are naturally defined by the saddles associated to the local maxima defining the clusters [START_REF] Chazal | Persistence-based clustering in riemannian manifolds[END_REF]. Summarizing, we consider clustering graphs (CG) defined as follows:

• Vertices. There is one node per cluster. The cluster is endowed with a representative point of that cluster (e.g. its centroid), used to compute the unit costs c ij . The weight of a node i.e. the supply or demand is the number of input points in that cluster.

• Edges. Selected pairs of clusters are connected, based on a predicate (clustering algorithm dependent).

To compare two clusterings, we enforce strict connectivity constraints, mapping connected regions associated with the first CG to connected regions associated with the second CG.

Data. We compare clusterings of data points drawn from a mixture of 5 anisotropic gaussians. More precisely, each gaussian is defined by its center and the rotation angle α (in degrees) of its principal direction, yielding a triple (c x , c y , α). Results. For a given pair, Alg-EMD-CCC-G is run twice since it is not symmetric. As for Alg-EMD-LP, the demand is always satisfied; but we monitor the fraction of vertices (r c.c. V ) and edges (r c.c. E ) satisfying the strict connectivity constraints. While d EMD-CCC < d EMD always holds, Alg-EMD-CCC-G falls short from satisfying the demand (% flow in the range 13-32%). On the other hand, Alg-EMD-LP fails from satisfying connectivity constraints, in particular for edges (r c.c. E ∈ 0.24 -0.45%). See Table 1. Table 1 Comparisons for the 3 clusterings of Fig. 3: performances. The 3 rows correspond to the three pairwise comparisons. Top table: Alg-EMD-CCC-G (two calculations since the algorithm is not symmetric) Bottom table: Alg-EMD-LP

Analysis of molecular energy landscapes

Setup. The potential energy landscape (PEL) of a molecular system is the graph of the function associating a potential energy to each conformation [START_REF] Wales | Energy Landscapes[END_REF]. PEL codes all thermodynamic and kinetic properties, and of particular importance are local minima and the (index one) saddles connecting them. Because PEL usually exhibits a number of minima exponential in the number of degree of freedom (3n of them for n atoms, e.g. n = 5, 000 for a medium sized protein), a pruning of local minima is generally in order. Topological persistence proved instrumental in this respect [START_REF] Cazals | Conformational ensembles and sampled energy landscapes: Analysis and comparison[END_REF][START_REF] Carr | Energy landscapes and persistent minima[END_REF]. Summarizing, we consider so-called transition graphs (TG):

• Vertices. One vertex for each persistent local minimum. The unit cost c ij between two vertices is taken as the root mean square deviation between the two conformations associated with the local minima. The weight (supply or demand) is theoretically defined by the integral of Boltzmann's factor over is catchment basin. Practically, using thermodynamic ensembles, it is estimated from the number of points in the basin.

• Edges. Two vertices sharing an index one saddle on the PEL.

We use our algorithms to compare two TG generated by two Monte Carlo runs. The comparisons aim at assessing the coherence between two explorations of the PEL, enforcing strict connectivity constraints to map connected regions of the two TG with one another. Practically, we use BLN69, a model protein with three types of pseudo amino-acids (a.a.) namely hydrophobic (B), hydrophylic (L) and neutral (N) [START_REF] Oakley | Energy landscape and global optimization for a frustrated model protein[END_REF]. This system has 69 a.a., whence a total of 207 Cartesian coordinates. The corresponding PEL has been thoroughly studied, with of the order of one million of local minima reported by state-of-the-art Monte Carlo methods [START_REF] Oakley | Energy landscape and global optimization for a frustrated model protein[END_REF][START_REF] Roth | Hybridizing rapidly growing random trees and basin hopping yields an improved exploration of energy landscapes[END_REF], and 10 low lying local minima identified.

We define a sampling as N = 10 4 conformations generated by T-RRT [START_REF] Jaillet | Randomized tree construction algorithm to explore energy landscapes[END_REF], a randomized incremental algorithm favoring the exploration of not-visited-yet regions. We further quench each sample p i to its local minimum q(p i ), by performing a gradient descent on the PEL. To convert the sampling into a TG, we use the Tomato algorithm [START_REF] Chazal | Persistence-based clustering in riemannian manifolds[END_REF]. This algorithm, which is reminiscent of clustering by mode seeking, identifies saddles at which the basins of nearby local minima merge.

Inria

The resulting graph is called the transition graph (TG). Furthermore, we simplify the TG using topological persistence, to retain the 50 most significant local minima. For this simplified TG, the mass of each basin is the fraction of samples located in its basin. This operation was carried out once for each local minimum in the top ten. We refer to this data set as TRRT-top10, and to a particular graph as TRRT-top10-i, i = 1, . . . , 10.

Our ten TG yield 45 pairs, whence 45 instances for Alg-EMD-LP (which is symmetric), and 90 instances for Alg-EMD-CCC-G.

Results: algorithm Alg-EMD-LP and constraint satisfaction. Since Alg-EMD-LP is oblivious to critical point connectivity, we compute the fraction r c.c. V and r c.c.

E of vertices and edges of the input graph inducing through the flow a connected subgraph of the demand graph. Out of the 45 instances of the dataset TRRT-top10, the min, median and max values for vertices and edges are (0.1, 0.62, 1.), and (0.03, 0.89, 1.), respectively. That is, transport plans obtained from solutions of the linear program do disrupt connectivity constraints.

Results: algorithm Alg-EMD-CCC-G and demand satisfaction. The connectivity constraints may prevent Alg-EMD-CCC-G to fully satisfy the demand. For each instance, we therefore monitor the total flow F Alg-EMD-CCC-G provided by the transport plan, the ideal value being one. On these 90 instances, a worst-case of 0.99 is observed. Further inspection shows that such performances owe to the distribution of weights in the basins. Indeed, for each transition graph TRRT-top10-i, it turns out that the local minimum from which the exploration was started takes most of the mass. Therefore, in comparing two such graphs, a transportation plan essentially reduces to moving the mass in-between the two prominent local minima. For this particular application, the fact that connectivity constraints are lenient shows that all runs discovered the same local minima and transitions. This stability is informative and gives confidence for physical analysis carried out downstream. Results: transport costs. To assess transport costs, we compute the linear correlation between three sets of 45 values, namely the transport costs of Alg-EMD-LP of the 45 instances, and those of Alg-EMD-CCC-G on the 45 × 2 pairs (recall that Alg-EMD-CCC-G is not symmetric). The three coefficients obtained are equal to 0.999, a property again owing to the structure of the basins, as just discussed.

Conclusion

This paper introduces optimal transportation problems which depart from classical ones as they embed connectivity constraints, and shows that these problems are in general hard to solve. A greedy polynomial time algorithm is also proposed for one of them, which is of particular interest to compare graphs representing clusterings and molecular energy landscapes. Our experiments show that optimizing the transport plan and respecting connectivity constraint can be competing objectives. On the other hand, comparable transport plans observed without and with connectivity constraints show that the supply and demand graphs are comparable, both in terms of embedding of vertices and edges. This information is especially interesting for selected applications, in particular (molecular) simulation, as it shows the stability of the data generation processes.

On the theoretical side, future work will aim at understanding whether problems with cost and connectivity constraints are hard to approximate, to possibly develop efficient algorithms. On the applied side, the importance of distance warping functions will be studied in the context of our two applications (comparing clusterings, analyzing molecular data).

Appendix 7.1 A simple example with strict connectivity constraints

In the following, we detail the example of Fig. 2: • Fig. 2 (a) describes a simple instance. Let V = {v 1 , v 2 , v 3 } be a set of three supply nodes such that X

1 = 8, X 2 = 5, X 3 = 4. Let V = {v 1 , v 2 , v 3 } be a set of three demand nodes such that Y 1 = 4, Y 2 = 3, Y 3 = 6.
Integers on nodes represent these supply and demand values. The graph G = (V, E) is a path, where

E = {{v 1 , v 2 }, {v 2 , v 3 }}. The graph G = (V , E ) is also a path, where E = {{v 1 , v 2 }, {v 2 , v 3 }}.
Integers on edges of the complete bipartite graph B represent unitary costs c i,j for all i, j ∈ {1, 2, 3}. The unitary costs are

c 1,1 = 1, c 1,2 = 7, c 1,3 = 1, c 2,1 = 6, c 2,2 = 1, c 2,3 = 9, c 3,1 = 9, c 3,2 = 5, and c 3,3 = 1.
• Fig. 2 (b) represents an optimal solution for EMD:

f 1,1 = 4, f 1,3 = 2, f 2,2 = 3, f 3,3 = 4, and f 1,2 = f 2,1 = f 2,3 = f 3,1 = f 3,2 = 0.
Only links of cost 1 are used and so the cost of the solution is C EMD = j∈{1,2,3} Y j = 13. This solution is not admissible for EMD-FCC. Indeed node v 1 ∈ V sends flow only to demand nodes v 1 ∈ V and v 3 ∈ V (that is f 1,1 > 0, f 1,2 = 0, and f 1,3 > 0), and the nodes v 1 and v 3 do not induce a connected sub-graph because {v 1 , v 3 } / ∈ E . One can observe that there does not exist an admissible solution of cost 13 for EMD-FCC even when M = |E(B)| = 9. In the following, we consider EMD-FCC with F = j∈J Y j = 13.

• Fig. 2 (c) represents an admissible solution for EMD-FCC for any M ∈ [START_REF] Cheng | Mean shift, mode seeking, and clustering[END_REF][START_REF] Hillier | Introduction to mathematical programming[END_REF], and for any real number x ∈ ]0, 3]:

f 1,1 = 4, f 1,2 = x, f 1,3 = 2, f 2,2 = 3 -x, f 3,3 = 4, and f 2,1 = f 2,3 = f 3,1 = f 3,2 = 0. The total cost is C EMD-FCC = i∈I j∈J f i,j c i,j = 6x + 13.
Thus, lim x→0 i∈I j∈J f i,j c i,j = 13 but we cannot obtain an admissible solution of cost 13 because x > 0.

• Fig. 2 (d) shows an optimal solution for EMD-FCC for M = 4:

f 1,1 = 4, f 2,2 = 3, f 2,3 = 2, f 3,3 = 4, and f 1,2 = f 1,3 = f 2,1 = f 3,1 = f 3,2 = 0. The total cost is C EMD-FCC = i∈I j∈J f i,j c i,j = 29.
• Fig. 2 (e) describes an optimal solution for EMD-FCC for M = 3:

f 1,3 = 6, f 2,2 = 3, f 3,1 = 4, and f 1,1 = f 1,2 = f 2,1 = f 2,3 = f 3,2 = f 3,3 = 0. The total cost is C EMD-FCC = i∈I j∈J f i,j c i,j = 45.
One can observe that there does not exist an admissible solution for EMD-FCC when 0 ≤ M ≤ 2.

Proofs for hardness -section 3

Proof of Lemma 1. We prove that the connectivity constraints are satisfied (for every H) if and only if these constraints are satisfied for

H = {H 1 , H 2 } such that |H 1 | = |H 2 | = 1 and H 1 = H 2 .
In the following, we set

H 1 = {v a } and H 2 = {v b } such that v a = v b .
(⇒) If the connectivity constraints are satisfied, then, clearly, the connectivity constraints are satisfied for every

H = {H 1 , H 2 } such that H 1 = {v a } and H 2 = {v b }, v a = v b .
(⇐) Suppose that the connectivity constraints are not satisfied for some H = {H 1 , . . . , H t }, t ≥ 1. Let H = {H 1 , . . . , H t }, where H i is the set of nodes of G that receive flow from at least one node of H i , for all i, 1 ≤ i ≤ t. We prove that there exists two nodes v 1 ∈ V (G), and

v 2 ∈ V (G), v 1 = v 2 such that either d G (H v1 , H v2 ) > g(d G (v 1 , v 2 )), or d G (cc(H v1 )) > g(0), or d G (cc(H v2 )) > g(0)
, where H vi is the sub-graph induced by the set of nodes that receive flow from v i , i ∈ {1, 2}. Without loss of generality, we have denoted such two nodes v 1 and v 2 . There are two cases:

• First, suppose that the connectivity constraints are not satisfied because d G (H ) > g(d G (H)).

It means that there exist i, j,

1 ≤ i < j ≤ t, such that d G (H i , H j ) > g(d G (H)). Observe that g(d G (H)) ≥ g(d G (H i , H j ))
. Furthermore, for every two nodes v 1 ∈ V (H i ), and

v 2 ∈ V (H j ), v 1 = v 2 , then d G (v 1 , v 2 ) ≤ d G (H i , H j ).
(d G (H)) ≥ g(d G (H i , H j )) and d G (v 1 , v 2 ) ≤ d G (H i , H j ), then g(d G (H)) ≥ g(d G (v 1 , v 2 )) because g is non-decreasing. Finally, we get that d G (H v1 , H v2 ) > g(d G (v 1 , v 2 )).
• Second, suppose that the connectivity constraints are not satisfied because there exists i,

1 ≤ i ≤ t, such that d G (cc(H i )) > g(d G (cc(H i ))). There exist j, j , 1 ≤ j < j ≤ s , such that d G (cc j (H i ), cc j (H i )) > g(d G (cc(H i ))
). Suppose first that there exist two nodes

u, v ∈ V (G) such that u ∈ cc k (H i ), v ∈ cc k (H i ) with 1 ≤ k < k ≤ s = |cc(H i )|
, and such that both u and v send flow to at least one node of cc j (H i ) and at least one node of cc j (H i ). Then,

d G (cc l (H u ), cc l (H v )) > g(d G (u, v)), where cc l (H u ) (resp, cc l (H v ))
is a maximal connected component of the sub-graph induced by the set of nodes that receive a flow from at least u or v that is contained in

cc j (H i ) (resp. cc j (H i )). In- deed, g(d G (u, v)) ≥ g(d G (cc(H i )) because g is non-decreasing and d G (cc l (H u ), cc l (H v )) > d G (cc j (H i ), cc j (H i ))
. Thus, it is done for this case (the connectivity constraints are not satisfied when choosing these two nodes). Then, suppose that such two nodes do not exist. Thus, it necessarily means that there exists one node v ∈ V (H i ) such that v sends flow to a set

H v ⊆ H i such that d G (cc(H v )) > g(d G ({v}))) = g(0)
. Thus the connectivity constraints are not satisfied for one node.

We have proved that the connectivity constraints are satisfied (for every H) if and only if these constraints are satisfied for every

H = {H 1 , H 2 } such that |H 1 | = |H 2 | = 1 and H 1 = H 2 .
Proof of Lemma 2. Consider an instance of 3-Partition problem. Let m ≥ 1 be any integer and let S = {n 1 , n 2 , . . . , n 3m } be a set of 3m positive integers. We construct the instance of EMD-FCC as follows. Set |I| = 3m and |J | = m. Set X i = n i for all i ∈ I. Let Z = i∈I X i . Without loss of generality, let Y j = Y with Z = mY . Set c i,j = 1 for all i ∈ I, j ∈ J . Let G = (V, E) be any connected graph and let G = (V , V × V ). Let F = Z, M = 3m, and D = m. Since G is a complete graph, the connectivity constraints are always satisfied (for any function g). We prove that there is an admissible solution for EMD-FCC if and only if there is a solution for the instance of 3-Partition problem.

(⇐) Assume there is a solution for the instance of 3-Partition problem, that is S can be partitioned into m subsets S 1 , S 2 , . . . , S m such that the sum of the numbers in each subset is equal. We construct our solution for EMD-FCC as follows. For all i ∈ I, j ∈ J , if n i ∈ S j , then set f i,j := n i , otherwise set f i,j := 0. By construction, we have 0 ≤ f i ≤ X i for all i ∈ I, j ∈ J . Since S 1 , S 2 , . . . , S m is a solution for the instance of 3-Partition problem, then i∈I f i,j = Y j for all j ∈ J . Finally we prove that the number of edges of B that support flow is (at most) M = 3m. By construction, f i,j1 f i,j2 = 0 for all i ∈ I, j 1 , j 2 ∈ J . Thus, for all i ∈ I, there is at most one edge adjacent to v i that supports flow. Thus, the solution is admissible because |I| = 3m.

(⇒) Assume there is an admissible solution for EMD-FCC. Since i∈I X i = j∈J Y j , then j∈J f i,j > 0 for all i ∈ I. In other words, there is at least one edge adjacent to v i that supports

Inria flow for all i ∈ I. Furthermore there is at most one edge adjacent to v i that supports flow for all i ∈ I because M = 3m = |I|. Thus, for all i ∈ I, there is exactly one edge adjacent to v i that supports flow. We construct a solution for the instance of 3-Partition problem as follows. For all i ∈ I, j ∈ J , if f i,j > 0, then n i ∈ S j , otherwise set n i / ∈ S j . By hypothesis (existence of an admissible flow), the sum of the numbers of S j is Y because i∈I X i = j∈J Y . Thus, there is a solution for the instance of 3-Partition problem.

In conclusion, the decision version of EMD-FCC is strongly NP-complete because it is in NP and 3-Partition problem is strongly NP-complete [START_REF] Garey | Computers and Intractability: A Guide to the Theory of NP-Completeness[END_REF].

Proof of Lemma 3. Let t ≥ 1 be any integer. Set |I| = t + 1 and set |J | = 2(t + 1). Let a ≥ t be any real number and let ε < 1/(p.t). Set X i := a + ε for all i ∈ I. Set Y j := a/2 + (j -1)(ε/t) for all j ∈ {1, . . . , t + 1}. Set

Y j := a/2 + (t -(j -(t + 2)))(ε/t) = a/2 + (2t -j2)(ε/t) for all j ∈ {t+2, . . . , 2(t+1)}. Set F = (a+ε)(t+1). Let G = (V, E) be any graph and let G = (V , E ) be such that E = {{v j , v j+1 }|j = 1, . . . , |J | -1} ∪ {{v |J | , v 1 }}.
Observe that G is a cycle. Let p ≥ t be any integer. We now define the unitary cost as follows. First, we denote by h(j) the index of the unique node v h(j) that is at distance t + 1 of node v j in G , for all j, 1 ≤ j ≤ |J |/2. Set c i,i = c i,h(i) = 1 for all i ∈ I and set c i,j = p for all i ∈ I and for all j ∈ J \ {i, h(i)}. Set C = (t + 1)(a + εp).

We now prove that there exists an admissible solution if and only if

M ≥ |I|(|I| + 1) = (t + 1)(t + 2) = (|J |/2)(|J |/2 + 1).
(⇐) Assume that M = (t + 1)(t + 2). We prove that there exists an admissible solution. Consider the following flow. For all i ∈ I, f i,i := a/2, f i,j := ε/t for every j ∈ {i + 1, . . . , t + 1}, f i,j := ε/t for every j ∈ {t + 2, . . . , h(i) -1}, f i,h(i) := a/2, and f i,j := 0 for every j ∈ {h(i) + 1, . . . , |J |}. Every node u ∈ V sends flow to a sub-path of size t + 2 of G . Thus, the connectivity constraints are satisfied and the number of edges that support flow is exactly (t+2)(t+1) = M because |I| = (t+1). Furthermore, by construction, each node v j ∈ V , j ∈ J , receives a flow of Y j . Thus, i∈I j∈J f i,j = F = j∈J Y j . Finally, i∈I j∈J f i,j c i,j = i∈I a + εp = (t + 1)(a + εp) by construction of f . (⇒) Assume that there exists an admissible solution. We prove that M ≥ (t + 1)(t + 2). By contradiction. Suppose that M < (t + 1)(t + 2). It means that there exists i * ∈ I such that j∈J |f i * ,j >0 1 < t + 2 because |I| = t + 1. For such a node v i * ∈ V , it is so impossible to send flow to v i * ∈ V and flow to v h(i * ) ∈ V . Indeed, since any shortest path between v i * and v h(i * ) is composed of t nodes, if node v i * sends flow to both v i * and v h(i * ) , then v i * must send flow to at least t other nodes because otherwise the connectivity constraints would not be satisfied. Thus, since j∈J |f i * ,j >0 1 < t + 2, then either v i * or v h(i * ) does not receive flow from v i * . Without loss of generality, suppose that such a node is v i * . By construction of the instance, it means that i∈I\{i * } f i,i * c i,i * ≥ p.a 2 . We get that i∈I j∈J f i,j c i,j ≥ p.a 2 + t.a because j∈J \{i * } f i,j c i,j ≥ t.a (the other nodes of G receive at least a volume of flow of a/2 and the cost unitary cost is at least 1). Since the solution is admissible, then we get that p.a 2 + t.a ≤ C = (t + 1)(a + εp). Set p := t, a := t, and ε := 1 p.t . We get that t 2 /2 ≤ t + 1 + 1/t. It is clearly false for sufficient large value of t. Thus, the solution is not admissible and we get a contradiction. Thus, it means that M ≥ (t + 1)(t + 2).

Proof of Lemma 4. Consider the instance I 1 of the proof of Lemma 3 and the instance I 2 of the proof of Lemma 2. We construct the instance I that is the union of I 1 and I 2 such that the unitary cost between any supply (demand, respectively) node of I 1 and any demand (supply, respectively) node of I 2 is infinite or sufficiently large. We set M = M I1 +M I2 = (t+1)(t+2)+3m. We also set F = F I1 + F I2 and C = C I1 + C I2 . (The index I x means that the parameter concerns the instance x ∈ {1, 2}.) By Lemma 3 and Lemma 2 (and Corollary 1), we get that the problem G , and B and the cost c i,j for all i ∈ I, j ∈ J are those of the original instance. Let ε > 0 be a real value such that X i -|J |ε > 0 for all i ∈ I and such that Y j -|I|ε > 0 for all j ∈ I. We denote by X i the volume of supply for all i ∈ I in the auxiliary instance. Set X i = X i -|J |ε for all i ∈ I. We denote by Y i the volume of demand for all j ∈ J in the auxiliary instance. Set Y j = Y j -|I|ε for all j ∈ I. Let f be an optimal solution for this auxiliary instance for EMD. Recall that this can be done in polynomial time since it reduces to solve a linear program. The cost of f is i∈I j∈J f i,j c i,j .

We now construct an admissible solution f for the original instance for EMD-FCC as follows. Set f i,j = f i,j + ε for all i ∈ I, j ∈ J . All the connectivity constraints are satisfied because for every v ∈ V , we have H v = G , where H v is the sub-graph induced by the of nodes of G that receive flow from v. Recall that M = |E(B)|. Thus, the solution is admissible. The cost of f is i∈I j∈J f i,j c i,j = |I||J |ε + i∈I j∈J f i,j c i,j . Let f * be an optimal solution for the original instance of EMD-FCC. Observe that Proof of lemma 11. Lines 1-6 of Algorithm 4 update the boolean function b if X i -(x it +f it,jt ) = 0, that is if the supply node v it sends flow for the first time. In that case, all the neighbors of v jt in G can receive flow from v it (Line 2) and all the other nodes of G cannot receive flow from v it (Line 3). Furthermore, all the neighbors of v it in G that do not have sent flow, cannot send flow to the nodes of G that are not neighbors of v jt in G (Lines 4-6). Lines 7-15 update the boolean function b if X i -(x it + f it,jt ) = 0, that is if the supply node v it has already sent flow before the current step. All the neighbors of v jt in G can receive flow from v it (Line 8). Every neighbor v i of v it in G that does not have sent flow, can send flow to every neighbor v j of v jt in G if b i,j = 0, that is if the edge {v i , v j } does not support flow (Lines 9-13). Furthermore, for every neighbor v i of v it in G that does not have sent flow, for every neighbor v j of v jt in G such that b i,j = 0, and for every neighbor v k of v i in G that has already sent flow and such that v k cannot send flow to v j , then we set that v i cannot send flow to v j , that is b i,j = 0 (Lines 9-15).

Algorithm 4 finally returns the variables b i,j for all i ∈ I, j ∈ J (Line 16).

Definition. 1 .

 1 (Distance warping function or DWF) A distance warping function is a function g : [0, |V |] ∪ ∞ → R + , mapping distances in G to distances in G , such that:

Figure 2

 2 Figure 2 Optimal transportation on graphs: a simple example. (a) The supply and demand graphs are paths; supply and demand are indicated in the nodes, while the edges are decorated with the unitary costs. (b) An optimal solution for EMD. (c) An admissible solution for EMD-FCC for any M ∈ [5, 9] and for any real number x ∈ ]0, 3]. (d) An optimal solution for EMD-FCC for M = 4. (e) An optimal solution for EMD-FCC for M = 3. See full detail in the supplemental section 7.1.

  0 ≤ C ≤ C max , 0 ≤ M ≤ |E(B)|, EMD-CCC is defined as follows:

Lemma. 8 .Algorithm 2

 82 The time complexity of Algorithm 1 is O(max(|V | 2 |V | 2 , |V ||V |C(g))), where C(g) is the time complexity of updating b with g representing the connectivity constraints. We propose in Algorithm 2 a polynomial time algorithm for updating the boolean function b used in Algorithm 1. Update of the boolean function b used in Algorithm 1 for general connectivity constraints.

4 :

 4 for all i ∈ I do 5:

  then

8 :

 8 for all i 1 ∈ I do 9:for all i 2 ∈ I do 10:

14 :

 14 b i2,j := 0 15: for all i ∈ I do 16:for all j ∈ J do 17: b i,j := b i,j 18: return b.

Corollary. 3 .

 3 The worst case time-complexity of Algorithm 1 with Algorithm 2 for the update of b, is O(max(|V | 4 |V |, |V | 3 |V | 4 log(|V |))).

1 : 4 : 5 :

 145 F inf := total flow returned by Algorithm 1 with G, G , C := C inf , M := |E(B)|, g 2: F sup := total flow returned by Algorithm 1 with G, G , C := C sup , M := |E(B)|, g 3: if F inf < F sup then Algorithm 3 with G, G , C inf , C sup := C, g Algorithm 3 with G, G , C inf := C, C sup , g4.4 Results for strict connectivity constraintsIn the sequel, we describe Algorithm 4, a polynomial time algorithm improving the generic update of the boolean function b (algorithm 2) for the particular case of strict connectivity constraints (Def. 3).We use the following notations. For any subset S ⊆ V , the open neighborhood N G (S) of S is the set of nodes in V \ S having a neighbor in S and the closed neighborhood of S, denoted by N G [S], is defined as N (S) ∪ S. If S = {v}, we use N G (v) and N G [v] instead of N G ({v}) and N G [{v}], respectively. We denote by ∆(G) the maximum degree of G.

Corollary. 4 .

 4 The worst case time-complexity of Algorithm 1 with Algorithm 4 for the update of b, is O

  Algorithm 1 is O(max(|V | 4 |V |, |V | 3 |V | 4 log(|V |))) (Corollary 3).

Algorithm 4

 4 Update of the boolean function b used in Algorithm 1 for the strict connectivity constraints.

  Algorithm 1 with Algorithm 4 for the update of b, is O(|V | 2 |V | 2 ) if ∆(G) = O(1) and ∆(G ) = O(1).

  The five gaussian used have parameters (-d, d, 45), (-d, -d, -45), (0, d/2., 0), (2 * d, 2 * d, 135), (2 * d, 2 * d, 45). Using d = {0, 40, 50} yields three mixtures. Clustering N = 5000 points drawn at random yields three clusterings, whence three pairwise comparisons.

Figure 3

 3 Comparing clustering of 5000 points drawn according to a mixture of 5 gaussians. From left to right: d = 50; middle: d = 40; right d = 0 (see text for details). Algorithm Tomato [4] yielded 44, 38 and 41 clusters, respectively, using a persistence threshold of 7.5 * 10 -3 .

dFigure 4

 4 Figure 4 Conformation of the BLN69 model.

  i∈I j∈J f i,j c i,j ≤ i∈I j∈J f * i,j c i,j and that i∈I j∈Jf * i,j c i,j ≤ |I||J |ε + i∈I j∈J f i,j c i,j = i∈I j∈J f i,j c i,j .We finally choose ε > 0 such that|I||J |ε + i∈I j∈J f i,j c i,j ≤ (1 + ε) i∈I j∈J f i,j c i,j .Thus, we geti∈I j∈J f i,j c i,j = |I||J |ε + i∈I j∈J f i,j c i,j and |I||J |ε + i∈I j∈J f i,j c i,j ≤ (1 + ε) i∈I j∈J f * i,j c i,j .We get a polynomial time (1 + ε)-approximation algorithm for EMD-FCC because f is obtained by solving a linear program and f is directly deduced from f . Proof of lemma 8. The number of steps (number of iterations of the while loop of Algorithm 1) is at most |V ||V |. At the beginning of each step, we select an edge. The time-complexity of Line 5 of Algorithm 1 to perform such edge selection is O(|V ||V |). We get the first term O(|V | 2 |V | 2 ). At the end of each step, we update b i,j, for all i ∈ I, j ∈ J . The time complexity of such a computation is C(g). We get the second term O(|V ||V |C(g)). To conclude, note that other computations are negligible. Proof of lemma 10. Before calling the main algorithm, we can initially compute the shortest paths between any two nodes of G by using the classical Floyd Warshall algorithm. We also do such a computation for any two nodes of G. The time complexity of such first step is O(|V | 3 + |V | 3 ). With such an initial computation, the time complexity of Lines 11-14 of Algorithm 2 is O(|V | 2 log(|V |) + log(|V |)). We get that the time complexity of Lines 8-14 is O(|V | 2 |V |(|V | 2 log(|V |) + log(|V |))) and that computation dominates the complexity of Lines 4-7. Then, the global time complexity is O(max(|V | 3 , |V | 3 , |V | 2 |V |(|V | 2 log(|V |) + log(|V |)))) because of the three for loops (Lines 8-10). After simplification, we get that the time complexity is O(max(|V | 3 , |V | 2 |V | 3 log(|V |))).

Proof of lemma 12 .

 12 The time complexity of the first part (Lines 1-6) is O(|V | + ∆(G)∆(G )), and the time complexity of the second part (Lines 7-15) is O(∆(G) 2 ∆(G )). Inria RESEARCH CENTRE SOPHIA ANTIPOLIS -MÉDITERRANÉE 2004 route des Lucioles -BP 93 06902 Sophia Antipolis Cedex Publisher Inria Domaine de Voluceau -Rocquencourt BP 105 -78153 Le Chesnay Cedex inria.fr

  Observe also that H v1 ⊆ H i andH v2 ⊆ H j because v 1 ∈ H 1 and v 2 ∈ H 2 . Indeed a node of G that receive flow from v 1 (resp. v 2 ), receive flow from at least one node of H i (resp. H j ). Thus, d G (H v1 , H v2 ) ≥ d G (H i , H j ). Recall that it is assumed that d G (H i , H j ) > g(d G (H)). Then, it follows that d G (H v1 , H v2 ) > g(d G (H)). Furthermore, since we have proved before that g

RR n°8991

of deciding if there exists an admissible solution is NP-complete. Finally, the number of nodes of I is |V ∪ V | = 3(t + 1) + 4m. By setting t = m -1, we get that M ∈ Θ(|V ∪ V | 2 ).

Proof of lemma 5. Suppose there is a constant k > 1 such that there is a polynomial time kapproximation algorithm for EMD-FCC.

Consider an instance of 3-Partition problem. Let m ≥ 1 be any integer and let S = {n 1 , n 2 , . . . , n 3m } be a set of 3m positive integers. Set |I| = 3m + 1 and

There exists a solution for EMD-FCC such that j∈J -f 3m+1,j + i∈I -f i,m+1 = 0 if and only if there is a solution for the instance of 3-Partition problem (Lemma 2). The cost of this solution is Y + i∈I,j∈J f i,j c i,j = Y + Z. We prove that if there does not exist a solution for the instance of 3-Partition problem, then the cost of any admissible solution for EMD-FCC is at least Z -1 + K. Suppose that there does not exist a solution for the instance of 3-Partition problem. Thus, we have j∈J f 3m+1,j > 0. There are two cases.

As we have supposed that there exists a polynomial time k-approximation algorithm for EMD-FCC, then if there is a solution of cost Y + Z, the kapproximation algorithm returns such a solution (otherwise the approximation ratio would be wrong); otherwise (solution of cost at least Z -1 + K), the k-approximation ratio would return a solution with cost at least Z -1 + K. Thus, the polynomial time (k-approximation) algorithm solves 3-Partition problem which is a strongly NP-complete problem [START_REF] Garey | Computers and Intractability: A Guide to the Theory of NP-Completeness[END_REF]. A contradiction, unless P=NP.

Proofs for algorithms -section 4

Proof of Lemma 6. First case. The connectivity constraints are always satisfied because G is a complete graph. Indeed, for every

, where H i is the sub-graph induced by the nodes that receive flow from the unique node that belongs to H i , i ∈ {1, 2}. Furthermore, the transport size constraint is also always satisfied because M = |E(B)|. Thus, EMD-FCC and EMD-CCC can be solved by linear programs that are variants of the linear program described in Eq. 1 for EMD.

Second case. The connectivity constraints are always satisfied because G has diameter at most g(x) for all x ≥ 0. Indeed for every

, for any i ∈ {1, 2}. Furthermore, the transport size constraint is also always satisfied because M = |E(B)|. Thus, EMD-FCC and EMD-CCC can be solved by linear programs that are variants of the linear program described in Eq. 1 for EMD.

Proof of lemma 7. In this proof, F is always chosen as the sum of the volumes of demands. Consider an instance of EMD-FCC. We construct an auxiliary instance as follows. The graphs G, Inria


Saber Salah
email: saber.salah@inria.fr

Reza Akbarinia
email: reza.akbarinia@inria.fr

Florent Masseglia
email: florent.masseglia@inria.fr

Mining Maximally Informative k-Itemsets in Massively Distributed Environments

The discovery of informative itemsets is a fundamental building block in data analytics and information retrieval. While the problem has been widely studied, only few solutions scale. This is particularly the case when i) the data set is massive, calling for large-scale distribution, and/or ii) the length k of the informative to be discovered is high. In this paper, we address the problem of parallel mining of maximally informative k-itemsets (miki) based on joint entropy. We propose PHIKS (Parallel Highly Informative K-ItemSet) a highly scalable, parallel miki mining algorithm. PHIKS renders the mining process of large scale databases (up to terabytes of data) succinct and effective. Its mining process is made up of only two efficient parallel jobs. With PHIKS, we provide a set of significant optimizations for calculating the joint entropies of miki having different sizes, which drastically reduces the execution time of the mining process. PHIKS has been extensively evaluated using massive real-world data sets. Our experimental results confirm the effectiveness of our proposal by the significant scale-up obtained with high itemsets length and over very large databases.

La découverte d'itemsets informatifs est un élément fondamental dans l'analyse de donnés et la recherche d'information. Bien que le problème a été largement étudié, il y a peu de solutions qui passent à l'échelle. Ceci est particulièrement le cas lorsque i) les données sont de très grane taille, ce qui demande une distribution à grande échelle, et / ou ii) la longueur k itemsets informatifs à découvrir est élevée. Dans cet article, nous abordons le problème de la fouille des k iems les plus informatifs (appelé miki) qui est calculé en considérant l'entropie conjointe des items. Nous proposons PHIKS (Parallel Highly Informative K-itemset), un algorithme parallèle d'extraction de miki. PHIKS rend le processus d'extraction de grandes bases de données à grande échelle (jusqu'à plusieurs téraoctets de données) rapide et efficace. Son processus d'extraction est constitué de seulement deux jobs parallèles. Avec * Saber Salah -This work was partially supported by the Inria Project Lab Hemera.

PHIKS, nous proposons un ensemble d'optimisations significatives

pour le calcul des entropies conjointes de miki ayant des tailles différentes, ce qui réduit considérablement le temps d'exécution du processus. PHIKS a été largement évaluée en utilisant de grands jeux de données réels. Nos résultats expérimentaux confirment l'efficacité de notre proposition qui passe à l'échelle pour de grandes valeurs de k et sur de très grandes bases de données.

THE PHIKS APPROACH

Featureset, or itemset, mining [START_REF] Han | Data mining : concepts and techniques[END_REF] is one of the fundamental building bricks for exploring informative patterns in databases. Features might be, for instance, the words occurring in a document, the score given by a user to a movie on a social network, or the characteristics of plants (growth, genotype, humidity, biomass, etc.) in a scientific study in agronomics. A large number of contributions in the literature has been proposed for itemset mining, exploring various measures according to the chosen relevance criteria. The most studied measure is probably the number of co-occurrences of a set of features, also known as frequent itemsets [START_REF] Agrawal | Fast algorithms for mining association rules in large databases[END_REF]. However, frequency does not give relevant results for a various range of applications, including information retrieval [START_REF] Greengrass | Information retrieval: A survey[END_REF], since it does not give a complete overview of the hidden correlations between the itemsets in the database. This is particularly the case when the database is sparse [START_REF] Heikinheimo | Finding low-entropy sets and trees from binary data[END_REF]. Using other criteria to assess the informativeness of an itemset could result in discovering interesting new patterns that were not previously known. To this end, information theory [START_REF] Cover | Elements of information theory[END_REF] gives us strong supports for measuring the informativeness of itemsets. One of the most popular measures is the joint entropy of an itemset. An itemset X that has higher joint entropy brings up more information about the objects in the database.

We study the problem of Maximally Informative k-Itemsets (miki for short) discovery in massive data sets, where informativeness is expressed by means of joint entropy and k is the size of the itemset [START_REF] Gray | Entropy and information theory[END_REF][START_REF] Knobbe | Maximally informative k-itemsets and their efficient discovery[END_REF][START_REF] Zhang | Discovering highly informative feature sets from data streams[END_REF]. Miki are itemsets of interest that better explain the correlations and relationships in the data. Example 1 gives an illustration of miki and its potential for real world applications such as information retrieval. EXAMPLE 1. In this application, we would like to retrieve documents from Table 1, in which the columns d1, d10 are documents, and the attributes A, B, C, D, E are some features (items, keywords) in the documents. The value "1" means that the feature occurs in the document, and "0" not. It is easy to observe that the itemset (D, E) is frequent, because features D and E occur together in almost every document. However, it provides little help for document retrieval. In other words, given a document dx in our data set, one might look for the occurrence of the itemset (D, E) and, depending on whether it occurs or not, she will not be able Features Documents d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 Miki mining is a key problem in data analytics with high potential impact on various tasks such as supervised learning [9], unsupervised learning [10] or information retrieval [START_REF] Greengrass | Information retrieval: A survey[END_REF], to cite a few. A typical application is the discovery of discriminative sets of features, based on joint entropy, which allows distinguishing between different categories of objects. Unfortunately, it is very difficult to maintain good results, in terms of both response time and quality, when the number of objects becomes very large. Indeed, with massive amounts of data, computing the joint entropies of all itemsets in parallel is a very challenging task for many reasons. First, the data is no longer located in one computer, instead, it is distributed over several machines. Second, the number of iterations of parallel jobs would be linear to k (i.e., the number of features in the itemset to be extracted [START_REF] Knobbe | Maximally informative k-itemsets and their efficient discovery[END_REF]), which needs multiple database scans and in turn violates the parallel execution of the mining process. We believe that an efficient miki mining solution should scale up with the in the size of the itemsets, calling for cutting edge parallel algorithms and high performance evaluation of an itemset's joint entropy in massively distributed environments.

A 1 1 1 1 1 0 0 0 0 0 B 0 1 0 0 1 1 0 1 0 1 C 1 0 0 1 0 1 1 0 1 0 D 1 0 1 1 1 1 1 1 1 1 E 1 1 1 1 1 1 1 1 1 1
We propose a deep combination of both information theory and massive distribution by taking advantage of parallel programming frameworks such as MapReduce [START_REF] Dean | Mapreduce: simplified data processing on large clusters[END_REF] or Spark [START_REF] Zaharia | Spark: Cluster computing with working sets[END_REF]. To the best of our knowledge, there has been no prior work on parallel informative itemsets discovery based on joint entropy. We designed and developed an efficient parallel algorithm, namely Parallel Highly Informative K-itemSet (PHIKS in short), that renders the discovery of miki a very large database (up to Terabytes of data) simple and effective. It performs the mining of miki in two parallel jobs. PHIKS cleverly exploits available data at each mapper to efficiently calculate the joint entropies of miki candidates. For more efficiency, we provide PHIKS with optimizations that allow for very significant improvements of the whole process of miki mining. The first technique estimates the upper bound of a given set of candidates and allows for a dramatic reduction of data communications, by filtering unpromising itemsets without having to perform any additional scan over the data. The second technique reduces significantly the number of scans over the input database of each mapper, i.e., only one scan per step, by incrementally computing the joint entropy of candidate features. This reduces drastically the work that should be done by the mappers, and thereby the total execution time.

PHIKS has been extensively evaluated using massive real-world data sets. Our experimental results show that PHIKS significantly outperforms alternative approaches, and confirm the effectiveness of our proposal over large databases containing for example one Terabyte of data.

  (c) 2016, Copyright is with the authors. Published in the Proceedings of the BDA 2016 Conference (15-18 November, 2016, Poitiers, France). Distribution of this paper is permitted under the terms of the Creative Commons license CC-by-nc-nd 4.0. (c) 2016, Droits restant aux auteurs. Publié dans les actes de la conférence BDA 2016 (15 au 18 Novembre 2016, Poitiers, France). Redistribution de cet article autorisée selon les termes de la licence Creative Commons CCby-nc-nd 4.0. .

Table 1 :

 1 Features in the documents

	to decide which document it is. By contrast, the itemset (A, B, C)
	is infrequent, as its member features rarely or never appear to-
	gether in the data. And it is troublesome to summarize the value
	patterns of the itemset (A, B, C). Providing it with the values
	(1, 0, 0) we could find the corresponding document d3; similarly,
	given the values (0, 1, 1) we will have the corresponding document
	d6. Although (A, B, C) is infrequent, it contains lots of useful in-
	formation which is hard to summarize. By looking at the values of
	each feature in the itemset (A, B, C), it is much easier to decide
	exactly which document they belong to. (A, B, C) is a maximally
	informative itemset of size k = 3.

ADDITIONAL AUTHORS