Pairwise Decomposition for Combinatorial Optimization in Graphical Models

1 Introduction

Probabilistic Graphical Models (PGM) provides a general approach for automated reasoning under uncertainty [Koller and Friedman, 2009]. PGMs cover Bayesian networks, Markov networks, but also deterministic frameworks such as constraint networks. This paper focuses on discrete PGMs and the MPE optimization problem which consists in maximizing a product of nonnegative functions over a set of discrete variables. In Section 2, we show how, using log-space transformation, one can reformulate the MPE problem as a WCSP, a general constraint optimization framework which consists in minimizing a sum of nonnegative cost functions over a set of discrete variables [Meseguer et al., 2006].

Exact algorithms for MPE (or WCSP) are based either on Depth-First Branch and Bound (DFBB) search algorithm or on inference methods, including variable elimination, jointree, and compilation. Inference methods generally run out of memory on large and complex networks (with large treewidth). The current state-of-the-art MPE solvers combine DFBB and memory-bounded inference methods (see for example [S´anchez et al., 2008; Marinescu and Dechter, 2009]). DFBB is a complete tree search method using linear memory space. During search, for minimization problems, the algorithm maintains an upper bound UB on the minimal cost solution, such as the cost of the best solution found so far. Moreover, each node in the search tree is associated with a lower bound LB, an underestimation of a minimum cost solution to the subproblem below the node. If LB ≥ UB, DFBB prunes the search space below the current node. Among the different techniques used to provide strong lower bound, soft local consistencies have been introduced in WCSPs. Since their worstcase time complexities are exponential in the maximum arity (number of variables) of the cost functions [Cooper et al., 2010], they are usually enforced on small arity (2 or 3) cost functions only. It is therefore desirable to decompose cost functions into smaller arity cost functions. Note that variable elimination may also benefit from this decomposition (as the number of neighboring variables may decrease).

This work has been partially funded by the french "Agence nationale de la Recherche", reference ANR-10-BLA-0214.

It is well-known that conditional independences (CI) allows to factorize probability distributions. However neither Bayesian networks nor Markov networks can explicitly represent all CI perfectly [Koller and Friedman, 2009] 1 . There may also be finer-grain context-specific independences, such as (X ⊥⊥ Y | Z = z), which is not explicit in the network and which can be exploited after entering evidence Z = z or performing variable elimination or value pruning by soft local consistency. In Markov networks, several factorizations may exist in the case of nonpositive distributions. In Section 3, we show how to exploit efficiently pairwise independence in WCSPs in the context of DFBB enhanced with soft local consistency and variable elimination.

2 Preliminaries

A Probabilistic Graphical Model (PGM) (or Gibbs distribution) is defined by a product of nonnegative functions F, over a set of discrete variables X , conveying probabilistic or deterministic information [Koller and Friedman, 2009]. Subsets of X will be denoted using bold letters such as S.

1 A famous example is the segregation network [Fishelson et al., 2005], where P(Gi,jp | Ga,jp,Ga,jm, Si,jp) = 2P(Gi,jp | Ga,jp, Si,jp)P(Gi,jp | Ga,jm, Si,jp) because (Ga,jp ⊥⊥ Ga,jm | Gi,jp, Si,jp), but this cannot be represented in a BN formulation.

Definition 1 (PGM). A PGM is a triplet (X , D, F) with X = {X1, . . . ,Xn}, a set of variables, D = {DX1 , . . . ,DXn }, a set of finite domains of maximum size d = maxni=1 dXi (dXi = |DXi |), and F = {f1, . . . , fe}, a set of nonnegative real valued functions, each defined over a subset of variables Si ⊆ X (i.e., the scope). It defines a joint distribution:

Bayesian networks represent a specific case using one conditional probability per variable (e = n) and a unit normalizing constant. The Most Probable Explanation (MPE) problem is to find the most likely assignment to all variables in X maximizing P(X ).

A Constraint Satisfaction Problem (CSP) is a triplet (X,D, C) with X,D defined as for PGM, and C, a set of constraints. Each constraint is defined as a boolean function over a subset of variables Si ⊂ X indicating which tuples in DSi = ∏X∈Si DX are authorized or not. The problem is to find a feasible assignment of X , satisfying all the constraints. A Weighted Constraint Satisfaction Problem (WCSP) is an extension of CSP targeted towards optimization.

Definition 2 (WCSP). A WCSP is a triplet (X,D,W) with X,D defined as for PGM, and W = {w1, . . . ,we}, a set of cost functions. Each cost function is defined over a subset of variables and takes nonnegative integer values in E+ = N ∪ {⊤} (⊤ = +∞ being associated to forbidden assignments 2 ). The goal is to find a feasible assignment of X minimizing ∑ei=1 wi(Si).

A PGM (X,D,F) can be translated into a WCSP (X , D, W) with ∀i ∈ [1, e], wi(Si) = ⌈-M log(fi(Si))+C⌉ (M,C being two positive constants for real-to-integer precision and nonnegativity of wi) which preserves the set of optimal solutions if a sufficiently large M value is used.

In the rest of the paper, we use f to denote cost functions. Given an assignment t ∈ DX , t[S] denotes the subassignment of t to the variables in S. Let f = f1 + f2 denote the sum (join) of two cost functions defined as f (t) = f1(t[S1]) + f2(t[S2]), ∀t ∈ DS1∪S2 . Let f = f1 - f2 denote the subtraction of two cost functions such that S2 ⊆ S1 and f (t) = f1(t) - f2(t[S2]), ∀t ∈ DS1 with ⊤ - ⊤ = ⊤. Let f [S′ ] denote the projection of a cost function over a subset S′ of its variables S such that S′ ⊆ S and ∀t′ ∈ DS′f [S′ ](t′ ) = mint∈DS s.t. t[S′ ]=t′ f (t). A cost function over two (resp. one) variable(s) is called a binary (resp. unary) cost function. An empty cost function has all its costs equal to zero and is removed from the WCSP.

Depth-First Branch and Bound (DFBB) is a complete tree search method to solve WCSPs. Soft local consistency methods produce strong lower bounds for DFBB by applying Equivalence-Preserving Transformations (EPT) [Cooper et al., 2010]. An arc EPT adds the projection f [X] of a cost function f (S) to a unary cost function f1(X),X ∈ S and subtracts f [X] from f (i.e. replaces f by f - f [X]) in order to get an equivalent problem. Enforcing soft arc consistency consists in applying arc EPTs for all cost functions and all directions (∀f (S) ∈ W, ∀X ∈ S) until all the projections are empty. A directional arc EPT adds the projection (f + f2)[X] from binary and unary cost functions f (X, Y ) and f2(Y ) to f1(X) and subtracts the result of this projection from f . Enforcing soft directional arc consistency (DAC) consists in applying directional arc EPTs on all cost functions in one direction only (X < Y ), as indicated by a total order on X , thereby ensuring termination. DFBB can be further improved using on-the-fly variable elimination [Larrosa, 2000]. Let F = {fi ∈ W s.t. X ∈ Si} and S = ⋃fi∈F Si, the elimination of variable X consists in replacing X and F by the projection (∑f∈F f )[S∖{X}] in the current WCSP. The elimination of a variable is exponential in the size of S and is cheap if |S| is small but also if X is assigned or if it is connected to the rest of the problem either through a single cost function or through at least one bijective (equality) constraint. DFBB enhanced with on-the-fly i-bounded elimination (all variables with |S| ≤ i are eliminated) and soft local consistency (EDAC for binary&ternary cost functions [S´anchez et al., 2008]) is denoted as DFBB-VE(i) in the experiments.

2 For simplicity, we use infinite ⊤, but our results are still valid for finite ⊤.

3 Pairwise cost function decomposition

We first define the notion of pairwise decomposition.

Definition 3. A pairwise decomposition of a cost function f (S) with respect to two variables X, Y ∈ S (X ̸= Y ), is a rewriting of f into a sum of two (nonnegative) cost functions f1(S ∖ {Y }) and f2(S ∖ {X}) such that:

An example of decomposable function is given in Figure 1. A pairwise decomposition replaces a cost function of arity r = |S| with smaller-arity (r - 1) cost functions. This decomposition process can be recursively repeated on the resulting cost functions f1, f2 until no pairwise decomposition can be found for each remaining cost function. A cost function f (S) that cannot be pairwise decomposed for any choice of X, Y ∈ S is said to be non decomposable. A WCSP is pairwise decomposed if all its cost functions are non decomposable. This property is enforced by applying the previously described iterative process to all cost functions.

The following theorem shows that testing if a cost function can be pairwise decomposed w.r.t X, Y is equivalent to testing pairwise independence between X and Y in an appropriate probability distribution.

Theorem 1 (Equivalence of pairwise independence and decomposition). Let S = {X, Y } ∪ Z a set of random variables, f (S) a cost function on S with at least one feasible assignment, and a distribution Pf = 1 ∑S exp(-f(S)) exp(-f (S)) (a Gibbs distribution parameterized by f). (X ⊥⊥ Y | Z) denotes that X and Y are pairwise independent given all other variables in Pf . Then,

Instead of testing pairwise independence in Pf , which involves summations and multiplications of real numbers, we propose a simpler test for pairwise decomposability based on equalities of cost differences. For a decomposable cost function f (X, Z, Y ) = f1(X, Z) + f2(Z, Y ) that takes only finite costs, and for any tuple z ∈ DZ, if we consider any fixed pair of values k, l for Y then the difference f (x, z, k) - f (x, z, l) = (f1(x, z) + f2(z, k)) - (f1(x, z) + f2(z, l)) = (f2(z, k) -f2(z, l)) does not depend on x. When f (X, Z, Y ) is not limited to finite costs however, subtracting an infinite cost from another cost is ill-defined. For example, if f2(z, k) and f2(z, l) are both equal to ⊤ then f (x, z, k) = f (x, z, l) = ⊤ for all x. Therefore, f1(x, z) can take any value and the decomposition will still hold. Infinite costs offer additional freedom in decomposition and must be specifically handled.

To design a test that can exploit infinite costs, we introduce an extended cost valuation structure E = Z ∪ {-⊤, ⊤, Ω} including negative costs and a special absorbing element Ω which captures the freedom generated by subtraction of infinite costs. Let a ⊟ b denote the difference of two elements a, b ∈ E: a ⊟ b = a - b, except for ⊤⊟⊤ = -⊤ ⊟ -⊤ = a ⊟ Ω = Ω ⊟ b = Ω, ⊤ ⊟ -⊤ = ⊤⊟ c = c ⊟ -⊤ = ⊤, -⊤ ⊟ ⊤ = -⊤ ⊟ c = c ⊟⊤ = -⊤ with c ∈ Z. The comparison of two elements a, b ∈ E denoted by a ⊜ b is true if and only if (1) a, b ∈ Z ∪ {⊤, -⊤} and a = b, or (2) a = Ω or b = Ω. We then have:

Theorem 2. A cost function f (X, Z, Y ) is pairwise decomposable w.r.t. X, Y iff ∀z ∈ DZ, ∀k, l ∈ DY (k < l):

Example 1. In Figure 1, we have f (W, X, Y, U ) being pairwise decomposable w.r.t. W and U due to the fact that: f (1111) ⊟ f (1112) ⊜ f (2111) ⊟ f (2112) : 5 ⊟⊤ ⊜ 4 ⊟ ⊤ f (1121) ⊟ f (1122) ⊜ f (2121) ⊟ f (2122) : ⊤⊟ 3 ⊜ ⊤⊟ 2 f (1211) ⊟ f (1212) ⊜ f (2211) ⊟ f (2212) : 6 ⊟ 5 ⊜ 2 ⊟ 1 f (1221) ⊟ f (1222) ⊜ f (2221) ⊟ f (2222) : ⊤⊟⊤ ⊜ 3 ⊟ 1

Next, we show how to construct f1, f2 of a pairwise decomposable cost function using projections and subtractions of cost functions (see Figure 1 for its application).

Theorem3. Let f (X, Z, Y ) be pairwise decomposable w.r.t. X, Y . Then f1(X,Z) = f [X,Z] and f2(Z, Y ) = (f - f1)[Z, Y ] is a valid decomposition of f, i.e. f = f1 + f2.

Enforcing pairwise decomposition on a cost function f (S) with r = |S| may require r(r-1) 2 tests, i.e. checking for all the pairs of variables in S. And it may produce two (non empty) (r - 1)-ary cost functions. Repeating pairwise decomposition on the resulting cost functions will produce at most min((r-p r ), 2p ) cost functions of arity (r - p) with a different scope. Each test is linear in the size of the (r - p)-ary cost functions, in O(d(r-p)+1 ). So the overall worst-case complexity of enforcing pairwise decomposition for e original cost functions with maximum arity r is O(e ∑r-1 p=0 2p (r - p)2d(r-p)+1 ) = O(e ∑r-1 p=0 (r - p)2dr+1 ) = O(e r(r+1)(2r+1) 6 dr+1 ) = O(er3dr+1 ) in time and O(e maxr-1 p=0 2pd(r-p) ) = O(edr ) in space.

At each step, the choice of the variables used for decomposition and the way costs are spread across f1 and f2 may influence the decomposability of f1 and f2 for next iteration. Under certain conditions, it is possible to find a sequence of pairwise decompositions (which pair of variables to select in f and how to spread f across f1 and f2) which will lead to an optimal decomposition (with minimal arity).

Property 1 (Unique minimal decomposition for finite cost functions [Hammersley and Clifford, 1971]). A cost function with only finite costs admits a unique minimal decomposition.

Sketch of proof. If a cost function f (S) has only finite costs, then its associated probability distribution Pf (see Theorem 1) is positive. The Hammersley & Clifford theorem [Hammersley and Clifford, 1971] 3 says that a positive distribution factorizes over a unique minimal Markov network H. The network is minimal in the sense that removing an edge will create a new conditional independence not satisfied by Pf . By taking the maximal cliques in H, each clique gives the scope of a factor (a nonnegative function) and Pf is then equal to the product of these factors divided by a normalizing constant. This factorization is equivalent to a sum of cost functions by taking - log(Pf ) (see the proof of Theo 1). From this resulting decomposition, it is easy to construct a reverse sequence of valid (nonnegative) cost function decompositions until the original cost function f is reached.

However, as soon as there are infinite costs, the previous theorem fails to apply. In this case, we propose a heuristic approach whose goal is to accumulate costs on the first variables in the DAC variable ordering. This should provide stronger DAC-based lower bounds. We therefore test pairs of variables in f (S) in a reverse DAC order 4 and project cost functions (Theorem 3) on the scope containing the first variables in the DAC order first. In Figure 1, assuming a lexicographic order, this corresponds to project first on {W, X, Y }. We show that our heuristic approach is able to find an optimal decomposition for some particular functions.

Property 2 (Tree structure identification for hard constraints [Meiri et al., 1990]). A cost function with only infinite costs that admits a decomposition into a tree of binary constraints is identifiable by our decomposition strategy.

3 See also Theorem 4.5 page 121 in [Koller and Friedman, 2009]. 4 We test Y,Z before W,X if max(Y,Z) > max(W,X) ∨ (max(Y, Z) = max(W, X) ∧ min(Y, Z) > min(W, X)).

This result follows from [Meiri et al., 1990] and the fact that a pairwise decomposition w.r.t. X, Y will remove a redundant edge {X, Y } in the minimal network 5 without affecting the set of solutions. For instance, a global equality constraint on S will be decomposed into a tree of binary equality constraints: f (S) = ∑Y ∈S∖{X} fY (X, Y ) with fY (X, Y ) ≡ (X = Y ).

Similarly, a linear cost function f (X1, . . . ,Xr) = ∑ri=1 aiXi will be decomposed into a sum of r unary cost functions 6 .

Related Works

Related works have considered the factorization of specific probability functions such as the noisy-or and its generalizations [Heckerman and Breese, 1996], or more general factorizations dedicated to probabilistic inference [Savicky and Vomlel, 2007]. Compared to our approach, both works add extra variables. Note that another possible approach to decrease the arity of cost functions is to switch to the dual representation [Meseguer et al., 2006]. This has been tested in [de Givry et al., 2003] for Max-3SAT and was apparently ineffective. This should likely worsen for larger domains.

4 Project&Subtract on binary cost functions

When a WCSP is pairwise decomposed, it is still possible to infer valuable information from the resulting nonbinary cost functions by projecting & subtracting them on smaller arity cost functions. We choose to project each nonbinary cost function on all the possible binary cost functions inside its scope 5 by following an order of projections&subtractions compatible with the DAC variable ordering (i.e. project first on the pair of variables with the smallest DAC positions). Each projection fi(X, Y ) = f [X, Y ] is followed by a subtraction f = f - fi in order to preserve the problem equivalence. Note that f may be empty after these subtractions and is then removed from the problem. Note also that the same binary cost function may receive projected costs from several nonbinary cost functions having overlapping scopes, resulting in stronger inference. The worst-case complexity of project&subtract applied on e cost functions with maximum arity r is O(er2dr ) in time and O(er2d2 ) in space.

5 Experimental results

DFBB-VE(i) 7 maintaining EDAC and using a dynamic conflict-based variable ordering [Lecoutre et al., 2009] obtained the best result in the UAI'08 8 (and UAI'10 9 for the 20minutes track) MPE Evaluation contest, except for two hard benchmarks: linkage analysis and grid networks.

The genetic linkage analysis networks. The problem instances have 334 to 1289 variables. The maximal domain size d is between 3 and 7 and the maximal arity is 5. The benchmark family consists of 22 problems.

The grid networks. Each problem is an l×l grid and each 3-ary CPT is generated uniformly randomly. 50, 75 or 90% of the CPTs are deterministic and all variables are Boolean. The instances have l ranging from 12 to 50 (i.e. n = 2500 variables). The benchmark family has 32 problems.

Experimentations 7 were performed on a 2.6 GHz Intel Xeon computer with 4GB running Linux 2.6. Total CPUtimes are in seconds and limited to 1 hour for each instance ("-" means time out). No initial upper bound is given. We tested enforcing pairwise decomposition (DFBBVE(i)+dec), project&subtract of n-ary cost functions (DFBBVE(i)+ps), and the combination of both techniques (DFBBVE(i)+dec+ps). Pairwise decomposition, project&subtract, and variable elimination of degree |S| ≤ i are performed in preprocessing only. On-the-fly variable elimination for degree at most 2 is used during search. All methods use MinFill variable elimination ordering heuristic and the reverse order is used for DAC. The following table summarizes the number of problems solved by the various versions of DFBBVE(i). Pairwise decomposition and project&subtract solved more problems than DFBB-VE(i) alone, and their combination obtained the best result, except for Linkage with i = 2, 3. Remarkably, all grids instances were solved with large ibounded variable elimination, showing the effect of pairwise decomposition, especially when there is sufficient determinism. A similar factorization was observed in [S´anchez et al., 2004] for CSPs.

We then compared the best method with AND/OR Branch and Bound with caching and static mini-bucket j-bound (AOBB-C+SMB(j)) which got the best results for Linkage and Grids at UAI'08 Evaluation. This method is implemented in aolibWCSP 10 . The value j is chosen as in [Marinescu and Dechter, 2009]. Table 1 gives problem sizes (n, d), treewidth (w), and summarizes the total CPU-time used by the various methods with the corresponding values i (choosing i ∈ [2, 7] which provides the best result) and j. The preprocessing "dec+ps" globally improves the results of DFBB-VE(i) and AOBB-C+SMB(j). The preprocessing time was always less than a second in our experiments (n′ , number of variables after preprocessing). DFBB-VE(i)+dec+ps gave the best results for all problems except ped44 and 90-38-1. Note that the range of variation of i for DFBB-VE(i)+dec+ps is smaller than the range of j for AOBB-C+SMB(j) 11 .

6 Conclusion

Pairwise decomposition combined with projections on binary cost functions and variable elimination is a powerful technique inside DFBB for MPE. Further work should be done on approximate pairwise decomposition and experiments performed on other problems, including probabilistic inference.

10 graphmod.ics.uci.edu/group/aolibWCSP

11 A good way of tuning i is to take the value corresponding to the maximum of the degree distribution for the considered graphical model (i.e. i = 5 for Linkage and i = 6 for Grids).

A Proof of Theorem 1

In the following, we use P as a shorthand for Pf and C = ∑S exp(-f (S)) > 0, a normalizing constant. By definition of conditional independence [Lauritzen, 1996], we have: (X ⊥⊥ Y | Z) ⇐⇒P(X, Y,Z) = P(X | Z)P(Y | Z)P(Z) ⇐⇒ P(X, Y, Z) = P(X,Z) P(Z) P(Y,Z), with P(Z) > 0.

=⇒ Assume P(X, Y, Z) = P(X | Z)P(Y | Z)P(Z). We have f (X, Z, Y ) = - log(CP(X, Y, Z)) = - log(CP(X | Z)P(Y | Z)P(Z)). Let define f ′1 (X,Z) = - log(P(X | Z)) and f ′2 (Z, Y ) = - log(P(Y | Z)P(Z)), two nonnegative cost functions. Thus f (X, Z, Y ) = f ′1 (X, Z) + f ′2 (Z, Y ) - log(C ). Because f is nonnegative, we have ∀x ∈ DX, ∀y ∈ DY , ∀z ∈ DZ, f ′1 (x, z) + f ′2 (z, y) - log(C ) ≥ 0. If log(C ) is negative, we just add - log(C ) to f ′1 or f ′2 in order to obtain f1, f2. Otherwise, for each z ∈ DZ, we decompose log(C ) = c1z + c2z into two positive numbers. Let ˆxz = argminx∈DX f ′1 (x, z) and ˆ yz = argminy∈DY f ′2 (z, y). Moreover, we have f ′1 (ˆxz, z) + f ′2 (z, ˆ yz ) - log(C ) = f (ˆxz , z, ˆ yz) ≥ 0. A feasible solution is c1z = min(f ′1 (ˆxz, z), log(C )) and c2z = log(C ) - c1z . Either c1z = log(C ) ≤ f ′1 (ˆxz, z) and c2z = 0, or c1z = f ′1 (ˆxz , z) and c2z = log(C ) - f ′1 (ˆxz, z). In this case, f ′2 (z, ˆ yz ) - c2z = f ′2 (z, ˆ yz )+f ′1 (ˆxz , z)-log(C ) ≥ 0, thus ∀y ∈ DY , f ′2 (z, y)- c2z ≥ 0. We define ∀x, y, z, f1(x, z) = f ′1 (x, z) - c1z and f2(z, y) = f ′2 (z, y) - c2z , which are nonnegative numbers. Finally, we deduce f (X, Z, Y ) = f1(X, Z) + f2(Z, Y ).

⇐= Assume we have three cost functions f (X, Z, Y ), f1(X, Z), f2(Z, Y ) such that f (X, Z, Y ) = f1(X,Z) + f2(Z, Y ). We have exp(-f (X, Z, Y )) = exp(-f1(X,Z)) exp(-f2(Z, Y )). By marginalization, we have P(X,Z) = ∑Y P(X, Y,Z) = 1C ∑Y exp(-f (X, Z, Y )) = 1C exp(-f1(X,Z)(∑Y exp(-f2(Z, Y ))), P(Y, Z) = ∑X P(X, Y, Z) = 1C exp(-f2(Z, Y ))(∑X exp(-f1(X, Z))) and P(Z) = ∑X,Y P(X, Y,Z) = 1C (∑X exp(-f1(X, Z)))(∑Y exp(-f2(Z, Y ))). Finally, we deduce P(X,Z) P(Z) P(Y,Z) = 1C exp(-f1(X, Z)) exp(-f2(Z, Y )) = P(X, Y, Z).

B Proof of Theorem 2

It is always possible to decompose a cost function f into three (nonnegative) cost functions f (X, Z, Y ) = f1(X, Z) + f2(Z, Y ) + Δ(X, Z, Y ) (f1, f2 possibly being equal to zero constant functions). A nonzero cost function is called reducible if f1 or f2 are nonzero cost functions. Otherwise it is called irreducible. We give a specific condition for testing reducibility.

Property 3. A nonzero cost function f (X, Z, Y ) is reducible if ∃x ∈ DX, ∃z ∈ DZ s.t. miny∈DY f (x, z, y) > 0 (i.e. f1 = f [X,Z] ̸= 0) or ∃y ∈ DY , ∃z ∈ DZ s.t. minx∈DX f (x, z, y) > 0 (i.e. f2 = f [Z, Y ] ̸= 0).

Property 4. If f1(X,Z), f2(Z, Y ), Δ(X,Z, Y ), three (nonnegative) cost functions, is a maximal reducing of f (X, Z, Y ) then Δ(X, Z, Y ) is irreducible.

Using Property 3, we prove the following lemma.

Lemma 1. If f (X, Z, Y ) is a nonzero irreducible cost function, then ∃z ∈ DZ, ∃k, l ∈ DY (k ̸= l), ∃u, v ∈ DX(u ̸= v) s.t. f (u, z, k) ⊟ f (u, z, l) ̸⊜ f (v, z, k) ⊟ f (v, z, l).

Proof by contradiction. Assume ∀z ∈ DZ, ∀k, l ∈ DY (k ̸= l), ∀u, v ∈ DX(u ̸= v) s.t. f (u, z, k) ⊟ f (u, z, l) ⊜ f (v, z, k) ⊟ f (v, z, l). Remember that cost functions always have costs in E+ = N ∪ {⊤} and not in E. We have either: First case: ∃u, z, k, l s.t. f (u, z, k) ⊟ f (u, z, l) ̸⊜ 0

So 0 ⩽ f (u, z, l) < f (u, z, k). Moreover ∀v ∈ DX, f (u, z, k) ⊟ f (u, z, l) ⊜ f (v, z, k) ⊟ f (v, z, l). Thus ∀v ∈ DX, f (v, z, k) ⊟ f (v, z, l) > 0, or f (v, z, k) ⊟ f (v, z, l) = ⊤, or f (v, z, k) ⊟ f (v, z, l) = Ω. Consequently 0 ⩽ f (v, z, l) < f (v, z, k) or f (v, z, k) = f (v, z, l) = ⊤. So min v∈DX f (v, z, k) > 0.

f (u, z, k) ⊟ f (u, z, l) < 0 or f (u, z, k) ⊟ f (u, z, l) = -⊤

So 0 ⩽ f (u, z, k) < f (u, z, l). Moreover ∀v ∈ DX, f (u, z, k) ⊟ f (u, z, l) ⊜ f (v, z, k) ⊟ f (v, z, l). Thus ∀v ∈ DX, f (v, z, k) ⊟ f (v, z, l) < 0, or f (v, z, k) ⊟ f (v, z, l) = -⊤, or f (v, z, k) ⊟ f (v, z, l) = Ω. Consequently 0 ⩽ f (v, z, k) < f (v, z, l) or f (v, z, k) = f (v, z, l) = ⊤. So min v∈DX f (v, z, l) > 0.

Second case: ∀u, z, k, l s.t. f (u, z, k) ⊟ f (u, z, l) ⊜ 0

f (u, z, k) ⊟ f (u, z, l) = Ω

So f (u, z, k) = f (u, z, l) = ⊤. Moreover ∀p ∈ DY , f (u, z, p) ⊟ f (u, z, l) = 0 or f (u, z, p) ⊟ f (u, z, l) = f (u, z, p) ⊟ ⊤ = Ω. Thus f (u, z, p) = ⊤. So f (u, z, l) = f (u, z, k) = f (u, z, p) = ⊤, thus min p∈DY f (u, z, p) = ⊤ > 0.

f

(u, z, k) =

f

(u, z, l)

>

0

So f (u, z, k) ⊟ f (u, z, l) = 0 and by assumption, ∀p ∈ DY , f (u, z, p) ⊟ f (u, z, l) = 0. Thus f (u, z, k) ⊟ f (u, z, l) = f (u, z, p) ⊟ f (u, z, l) and 0 < f (u, z, l) = f (u, z, k) = f (u, z, p). So min p∈DY f (u, z, p) > 0.

f (u, z, k) = f (u, z, l) = 0

By assumption, ∀p ∈ DY , f (u, z, p) ⊟ f (u, z, l) ⊜ 0. Then f (u, z, p) ⊟ 0 ⊜ 0. Thus f (u, z, p) ⊟ f (u, z, l) ̸= Ω and so f (u, z, k) ⊟ f (u, z, l) = f (u, z, p) ⊟ f (u, z, l). Thus ∀p ∈ DY , f (u, z, p) = 0. From all these cases, we conclude using Property 3 that f (X, Z, Y ) is reducible or equal to the zero function (if always being in the last sub-case).

Let a ⊞ b denote the addition of two elements a, b ∈ E : a ⊞ b = a + b except for ⊤⊞⊤ = ⊤⊞ c = c ⊞⊤ = ⊤, -⊤ ⊞ -⊤ = -⊤ ⊞ c = c ⊞ -⊤ = -⊤, ⊤ ⊞ -⊤ = -⊤ ⊞ ⊤ = a ⊞ Ω = Ω ⊞ b = Ω with c ∈ Z.

Now, we are able to prove Theorem 2.

=⇒ Assume f (X, Z, Y ) = f1(X,Z) + f2(Z, Y ) and f1(X, Z), f2(Z, Y ) nonnegative functions (0 ⩽ f1(X, Z) ⩽ f (X, Z, Y )). Then the following systems Sz of linear equations (using ⊞ and ⊜ operators) must have a solution.

∀z ∈ DZ, ∀k, l ∈ DY , ∀u ∈ DX, from Sz, we deduce f(uzk) ⊟ f(uzl) ⊜ (f1(uz) ⊞ f2(zk)) ⊟ (f1(uz) ⊞ f2(zl)) ⊜ (f1(uz) ⊟ f1(uz)) ⊞ (f2(zk) ⊟ f2(zl))

1st case: f1(uz) ∈ N

f (uzk) ⊟f (uzl) ⊜ (f1(uz) ⊟f1(uz)) ⊞ (f2(zk) ⊟f2(zl)) ⊜ 0 ⊞ (f2(zk) ⊟ f2(zl)) ⊜ f2(zk) ⊟ f2(zl)

2nd case: f1(uz) = ⊤

f (uzk) ⊟f (uzl) ⊜ (f1(uz) ⊟f1(uz)) ⊞ (f2(zk) ⊟f2(zl)) ⊜ Ω ⊞ (f2(zk) ⊟ f2(zl)) ⊜ Ω

We conclude that ∀z ∈ DZ, ∀k, l ∈ DY (k ̸= l) ⊜u∈DX f (uzk) ⊟ f (uzl)

⇐= Assume ∀z ∈ DZ, ∀k, l ∈ DY (k < l) : ⊜u∈DX f (uzk) ⊟ f (uzl).

Moreover, we have f (X, Z, Y ) = f1(X, Z) + f2(Z, Y ) + Δ(X, Z, Y ) f1(X,Z), f2(Z, Y ) and Δ(X,Z, Y ) nonnegative functions which defines the following linear systems:

f1(uz) ⊞ f2(zk) ⊞ Δ(uzk) ⊜ f(uzk) ∀u, k

f1(dXz) ⊞ f2(zdY ) ⊞ Δ(dXzdY ) ⊜ f(dXzdY )

(S ′z ) { f1(1z) ⊞ f2(z1) ⊞ Δ(1z1) ⊜ f(1z1)

f (uzk) ⊟ f (uzl)

⊜ (f1(uz) ⊞ f2(zk) ⊞ Δ(uzk)) ⊟ (f1(uz) ⊞ f2(zl) ⊞ Δ(uzl))

⊜ f1(uz) ⊞ f2(zk) ⊞ Δ(uzk) ⊞ (-f1(uz)) ⊞ (-f2(zl))

⊞(-Δ(uzl))

⊜ (f1(uz) ⊟ f1(uz)) ⊞ (f2(zk) ⊟ f2(zl)) ⊞ (Δ(uzk) ⊟ Δ(uzl))

Assume f (X, Z, Y ) does not decompose in {f1(X, Z), f2(Z, Y )} and with the property 4 we can find Δ(X,Z, Y ) such that Δ(X,Z, Y ) is a nonzero irreducible cost function. With the theorem 1 we have ∃z ∈ DZ, k, l ∈ DY , k < l and u, v ∈ DX, u < v s.t. Δ(uzk) ⊟ Δ(uzl) ̸⊜ Δ(vzk) ⊟ Δ(vzl).

We deduce that Δ(uzk) ⊟ Δ(uzl) ̸= Ω and Δ(vzk) ⊟ Δ(vzl) ̸= Ω with the definition of ⊜. Also, we have f (uzk) ⊟ f (uzl) ̸= Ω ̸= f (vzk) ⊟ f (vzl). Thus f (uzk) ⊟ f (uzl) ⊜ f (vzk) ⊟ f (vzl) ⊜ ϕ s.t. ϕ ∈ Z ∪ {⊤,-⊤} so f1(uz) ⊟ f1(uz) = 0 = f1(vz) ⊟ f1(vz) et f2(zk) ⊟ f2(zl) ̸= Ω

We use these properties below:

Δ(uzk) ⊟ Δ(uzl) ̸= Δ(vzk) ⊟ Δ(vzl)

f2(zk) ⊟ f2(zl) ⊞ Δ(uzk) ⊟ Δ(uzl) ̸= f2(zk) ⊟ f2(zl) ⊞ Δ(vzk) ⊟ Δ(vzl) f1(uz) ⊟ f1(uz) ⊞ f2(zk) ⊟ f2(zl) ⊞ Δ(uzk) ⊟ Δ(uzl) ̸= f1(vz) ⊟ f1(vz) ⊞ f2(zk) ⊟ f2(zl) ⊞ Δ(vzk) ⊟ Δ(vzl) f1(uz) ⊞ f2(zk) ⊞ Δ(uzk) ⊟ (f1(uz) ⊞ f2(zl) ⊞ Δ(uzl)) ̸= f1(vz) ⊞ f2(zk) ⊞ Δ(vzk) ⊟ (f1(vz) ⊞ f2(zl) ⊞ Δ(vzl)) f(uzk) ⊟ f(uzl) ̸= f(vzk) ⊟ f(vzl)

This is in contradiction with the first hypothesis : ∀ k, l ∈ DY , ∀ z ∈ DZ f (uzk) ⊟ f (uzl) ⊜ f (vzk) ⊟ f (vzl).

C Proof of Theorem 3

Lemma 2. Let a cost function f (X, Z, Y ) be a pairwise decomposition w.r.t. X, Y. If f1(x, z) = miny∈DY f (x, z, y) then ∀z ∈ DZ ∃k ∈ DY such that ∀x ∈ DX, f1(x, z) = f (x, z, k) .

Proof. Let z ∈ DZ.

∃x ∈ Dx, ∃y, f (x, z, y) ̸= ⊤ Let DK = {k1, k2, . . . , kdY } such that ∀i < j, f (x, z, ki) ⩽ f (x, z, kj). We deduce ∀i ∈ [1, dY ], f (x, z, k1) ⩽ f (x, z, ki) and f1(x, z) = miny∈DY f (x, z, y) = f (x, z, k1).

Using Theorem 2, we find ∀kp ∈ DK ∀u ∈ DX, 0 ⩽ f (x, z, kp) ⊟ f (x, z, k1) ⊜ f (u, z, kp) ⊟ f (u, z, k1). So f (u, z, k1) ⩽ f (u, z, kp) or f (u, z, k1) = f (u, z, kp) = ⊤, thus ∀u ∈ DX, f1(x, z) = miny∈DY f (u, z, y) = f (u, z, k1).

Now, we are able to prove Theorem 3. We have f1(x, z) = miny∈DY f (x, z, y) and f2(z, y) = minx∈DX [f (x, z, y) - f1(x, z)] so f2(z, y) ⩽ f (x, z, y) - f1(x, z)

If f2(z, y) = f (x, z, y) - f1(x, z), f (x, z, y) = f1(x, z) + f2(z, y).

If f2(z, y)

< f

(x, z, y)

- f1(x, z) ,

f (x, z, y) = f1(x, z) = ⊤

We have f1(x, z) + minx∈DX [f (x, z, y) - f1(x, z)] = ⊤ = f (x, z, y) so f (x, z, y) = f1(x, z) + f2(z, y).

f (x, z, y) ̸= ⊤ or f1(x, z) ̸= ⊤

u = argminu∈DX [f (u, z, y)-f1(u, z)]. Using Lemma 2, we have l = argmink∈DY f (x, z, k) = argmink∈DY f (u, z, k), thus f (u, z, l) - f (u, z, l) < f (x, z, y) - f (x, z, l), but also f (u, z, y) ⊟ f (u, z, l) ⊜ f (x, z, y) ⊟ f (x, z, l) because f is pairwise decomposable w.r.t. X, Y . Moreover f (x, z, y) ̸= ⊤ or f (x, z, l) ̸= ⊤, thus f (u, z, y)-f (u, z, l) = f (x, z, y)- f (x, z, l). Finally f2(z, y) = f (x, z, y) - f1(x, z), thus f (x, z, y) = f1(x, z) + f2(z, y).