Julien Ah-Pine
email: julien.ah-pine@eric.univ-lyon2.fr

Graph Clustering by Maximizing Statistical Association Measures

Keywords: Graph clustering, Community detection, Statistical association measures, Modularity

We are interested in objective functions for clustering undirected and unweighted graphs. Our goal is to define alternatives to the popular modularity measure. To this end, we propose to adapt statistical association coefficients, which traditionally measure the proximity between partitions, for graph clustering. Our approach relies on the representation of statistical association measures in a relational formulation which uses the adjacency matrices of the equivalence relations underlying the partitions. We show that graph clustering can then be solved by fitting the graph with an equivalence relation via the maximization of a statistical association coefficient. We underline the connections between the proposed framework and the modularity model. Our theoretical work comes with an empirical study on computer-generated graphs. Our results show that the proposed methods can recover the community structure of a graph similarly or better than the modularity.

Introduction

Many real-world problems can be designed using graphs where entities of the studied system are represented as nodes and their relationships as edges between nodes. In many domains such as biology, ecology, social network analysis . . . , graph theory tools are employed as means for representing complex systems. In this context, graph clustering consists in partitioning nodes into groups such that vertices belonging to the same group are better interconnected to each other than to vertices outside of the group. Discovering such clusters can lead to new and important insights. In biology for example, clustering a protein-protein interaction network helps to find proteins with the same biological function. Another example is in social network analysis, where graph clustering leads to the detection of community structures [START_REF] Papadopoulos | Community detection in social media[END_REF]. Such knowledge can help to better understand the social system and its related phenomenons.

There exist many graphs clustering techniques. We particularly focus on methods that optimize an objective function. The benefit criterion aims at reflecting the quality of a clustering. In this context, density-based objective functions are well-known approaches. In these cases, clusters are defined as subgraphs with high densities of edges. The modularity measure proposed by Newman and Girvan in [START_REF] Newman | Finding and evaluating community structure in networks[END_REF] is a popular density-based objective function. It assumes that two nodes belong to the same community if the number of edges between them is greater than the expected number of edges under a null random model.

We address the graph clustering task from a viewpoint different from the one underlying the modularity. We suppose that an undirected and unweighted graph can be seen as a perturbed equivalence relation and finding groups of nodes can be interpreted as fitting the graph with a partition. To this end, we need to quantify the proximity between two partitions. In the statistical literature there are numerous coefficients addressing this exact problem. These criteria are known as statistical association measures (SAM) between categorical variables or partitions. Our proposal is thus to fit a given graph with a partition by maximizing a SAM. However, using such measures in this context is not straightforward. Indeed, these coefficients are typically defined by using contingency tables over the set of categories of the two partitions. Yet, the contingency table between a given graph (which is not an equivalence relation) and a partition does not exist. To overcome this drawback, we review the research works of Marcotorchino who showed in [START_REF] Marcotorchino | Utilisation des comparaisons par paires en statistique des contingences partie I[END_REF][START_REF] Marcotorchino | Utilisation des comparaisons par paires en statistique des contingences partie II[END_REF], that many SAM can be equivalently expressed through the adjacency matrices of the equivalence relations underlying the categorical variables. Based on this approach, we show how we can convert SAM to define new density-based quality functions for graph clustering.

In section 2, we recall some density-based objective functions for graph clustering. We particularly emphasize the modularity concept. Then, in section 3, we introduce our framework. We recall SAM both in their contingency and their relational formulations. Then we show how these measures can lead to graph clustering methods. Moreover, we study the relationships between the modularity and SAM. Next, in section 4, we empirically examine the behaviors of the proposed objective functions on artificial graphs and we compare their results with the ones provided by the modularity. We finally conclude and sketch some future works in section 5.

Related work : modularity optimization

There are several types of density-based benefit functions for graph clustering [START_REF] Schaeffer | Graph clustering[END_REF][START_REF] Fortunato | Community detection in graphs[END_REF]. One first family is based on graph cuts measures which iteratively split the set of nodes of a graph into two, providing that the density of edges between the two clusters is low. To apply such methods, one can generally use any max-flow/min-cut algorithm such as the Ford-Fulkerson one. Another method is spectral clustering which computes the Fiedler eigenvector of the Laplacian of the graph. Edges cuts criteria and the aforementioned algorithms are particularly used to tackle graph partitioning problems. These tasks are slightly distinct from graph clustering problems. In graph partitioning, the number of clusters and their sizes are known and one has to recover the correct partition given these pieces of information. In contrast, in graph clustering, we do not assume any information about the number nor the shape of the communities.

In order to better deal with the graph clustering task, Newman and Girvan proposed the modularity concept [START_REF] Newman | Finding and evaluating community structure in networks[END_REF]. Their approach has the advantage to better formalize the concept of community and to avoid setting the number of clusters manually. This measure is denoted by Q and it can be expressed as follows [START_REF] Newman | Finding community structure in networks using the eigenvectors of matrices[END_REF] : Q = Number of edges within communities -Expected number of such edges.

More formally, let us assume that we are given a graph with n vertices and m edges. Its adjacency matrix is denoted by A. Let us denote by P the pairwise matrix of general term P ij which is the expected number of edges between nodes i and j. Since we are concerned with undirected and unweighted graph, P ij can be interpreted as the probability to have an edge between i and j. The modularity can be formulated by the equation below [START_REF] Newman | Finding community structure in networks using the eigenvectors of matrices[END_REF] :

Q(A, δ) = 1 2m n i=1 n j=1 (A ij -P ij ) δ(g i , g j ) . ( 1
)
where g i is the cluster of i and δ(g i , g j ) = 1 if g i = g j and 0 otherwise. From this general formulation, Newman adopted different assumptions which led to the definition of a specific null random model. His hypothesis are the following ones [START_REF] Newman | Finding community structure in networks using the eigenvectors of matrices[END_REF] : (i) since the graph is undirected then P should satisfy the relation, ∀i, j : P ij = P ji ; (ii) Q should be null when all vertices are in a single group and thus 1 i,j A ij = i,j P ij = 2m; (iii) the degrees distribution of the random model should be approximately the same as the one of the given graph which leads to the following constraint, ∀i : j P ij = k i where k i = j A ij is the observed degree of node i; (iv) edges should be placed at random meaning that the probability of observing an edge between i and j should be independent from the probability of observing an edge involving i and the probability of observing and edge involving j. Under these assumptions, the simplest null random model is when ∀i, j : P ij = k i k j /2m [START_REF] Newman | Finding community structure in networks using the eigenvectors of matrices[END_REF]. Accordingly, the following modularity formulation is the one which is commonly used in the literature :

Q(A, δ) = 1 2m n i=1 n j=1 A ij - k i k j 2m δ(g i , g j ) . (2)
It is worthwhile to mention that apart from (2), other coefficients relying on the modularity concept could be designed from [START_REF] Papadopoulos | Community detection in social media[END_REF]. In that perspective, Newman suggested that the assumptions (i) and (ii) are fundamental and should be considered as axioms of the modularity framework unlike (iii) and (iv) [START_REF] Newman | Finding community structure in networks using the eigenvectors of matrices[END_REF].

Adopting [START_REF] Newman | Finding and evaluating community structure in networks[END_REF], one can optimally solve the graph clustering problem via modularity maximization with the following integer linear program (see for e.g. [START_REF] Agarwal | Modularity-maximizing graph communities via mathematical programming[END_REF]) :

max Y 1 2m i,j A ij - k i k j 2m (1 -Y ij ) subject to : Y ik ≤ Y ij + Y jk ∀i, j, k Y ij ∈ {0, 1} ∀i, j . (3)
1 In order to lighten the notations we write i,j as a shortcut for n i=1 n j=1 .

where Y ij = 1 if i and j are not in the same cluster and 0 otherwise. However, optimizing the modularity (and any other objective functions) over the set of possible partitions is an NP-hard problem. As a result, many research works have been devoted to approximately maximize the modularity with different strategies and heuristics [START_REF] Schaeffer | Graph clustering[END_REF][START_REF] Fortunato | Community detection in graphs[END_REF].

The application of modularity to the graph clustering task has demonstrated good performances both on artificial and real-world networks. However, some recent works have highlighted certain limits of this method [START_REF] Lancichinetti | Limits of modularity maximization in community detection[END_REF]. In particular, optimizing the modularity tends to split large groups while small communities below a certain threshold are not correctly detected.

In the next section, we introduce new quality functions for graph clustering which provide alternatives to the modularity criterion given in (2).

The proposed approach : statistical association measures (SAM) optimization

Density-based techniques typically rely on the definition of a community and use heuristics to discover sub-graphs satisfying this definition. From our viewpoint, since clustering a given graph consists in detecting a hidden community structure, we can interpret the graph as an equivalence relation perturbed by noise. Thereby, we argue that graph clustering can be thought of as recovering the real community structure and this can be achieved by fitting the graph with a partition. This approach assumes there is a way to assess the proximity between the graph and a partition. In what follows, we introduce some statistical association measures which aim at measuring the similarity between two partitions by means of contingency tables. Then, we recall the relational formulation of these coefficients due to Marcotorchino. Using the latter expressions of SAM, we show how to use these coefficients as benefit functions for the graph clustering task. In that perspective, we underline some theoretical links between SAM and the modularity concept in order to bring into light some conceptual similarities between these two frameworks in the context of graph clustering.

SAM and their relational representation

Let us assume two categorical variables V k and V l with respectively p k and p l categories. Note that a categorical variable infers a set of disjoint groups of items which in turn can be interpreted as a partition or a clustering or an equivalence relation 2 . In categorical data analysis, in order to analyse the relationship between two categorical variables, we use the contingency table of dimensions (p k × p l ) denoted by N whose general term is defined by : N uv = Number of items belonging to both category u of V k and category v of V l . Then, a core concept in categorical data analysis is the deviation from the statistical independence situation which occurs when for all pairs of categories 2 Therefore, we will use these different terms interchangeably.

(u, v), the probability of jointly observing u and v equals the product between the probability of observing u and the probability of observing v. Using N, this principle translates into the following formula3 : ∀(u, v) : N uv /n = (N u. N .v )/n 2 where N u. = v N uv . In this context, the greater the difference between N uv and (N u. N .v )/n (for all pairs (u, v)), the stronger the relationship between the categorical variables.

Accordingly, we propose to study the following coefficients :

B(V k , V l ) = p k u=1 p l v=1 N uv - N u. N .v n 2 . (4)
E(V k , V l ) = p k u=1 p l v=1 N 2 uv - N 2 u. N 2 .v n 2 .
(5)

J(V k , V l ) = 1 n p k u=1 p l v=1 N uv N uv - N u. N .v n . ( 6
)
LM (V k , V l ) = p k u=1 p l v=1 N 2 uv N u. - 1 n p l v=1 N 2 .v . (7)
The SAM B, E, J and LM are respectively the Belson [10], Marcotorchino's square independence deviation [START_REF] Marcotorchino | Utilisation des comparaisons par paires en statistique des contingences partie I[END_REF], the Jordan 4 [11] and the Light-Margolin [START_REF] Light | An analysis of variance for categorical data[END_REF] criteria. They are all null in case of statistical independence. B and LM can only have positive values while E and J can be either positive or negative [START_REF] Marcotorchino | Utilisation des comparaisons par paires en statistique des contingences partie I[END_REF]. Given V k , these coefficients achieve their maxima when V l is exactly the same partition as V k [START_REF] Marcotorchino | Utilisation des comparaisons par paires en statistique des contingences partie I[END_REF].

The contingency representation is the usual way to introduce SAM. However, there exists an equivalent representation of these coefficients which emphasizes the relational nature of categorical variables. Indeed, as we mentioned beforehand, categorical variables are equivalence relations and such algebraic structures can be represented by graphs. This point of view was adopted by Marcotorchino and enabled him to formulate SAM with adjacency matrices5 [START_REF] Marcotorchino | Utilisation des comparaisons par paires en statistique des contingences partie I[END_REF][START_REF] Marcotorchino | Utilisation des comparaisons par paires en statistique des contingences partie II[END_REF]. Let us denote by C k the adjacency matrix6 associated to V k . Its general term is defined by C k ij = 1 if i and j belong to the same category and 0 otherwise. Marcotorchino provided correspondence formulas between the contingency table N on the one hand and the relational representations C k and C l on the other hand [START_REF] Marcotorchino | Utilisation des comparaisons par paires en statistique des contingences partie I[END_REF][START_REF] Marcotorchino | Utilisation des comparaisons par paires en statistique des contingences partie II[END_REF]. Here are some of these transformation formulas : (i)

p k u=1 p l v=1 N 2 uv = n i=1 n j=1 C k ij C l ij ; (ii) u N 2 u. = i,j C k ij ; (iii) u,v N uv N u. N .v = i,j ((C k i. + C k .j )/2)C l ij ; (iv) u,v (N 2 uv /N u. ) = i,j (2C k ij /(C k i. + C k .j ))C l ij where C k i. = j C k ij and C k i. = C k .i since C k is symmetric.
Applying these correspondence formulas enables the following expressions of SAM in terms of C k and C l :

B(C k , C l ) = n i=1 n j=1 C k ij - C k i. n - C k .j n + C k .. n 2 C l ij . (8)
E(C k , C l ) = n i=1 n j=1 C k ij - C k .. n 2 C l ij . (9)
J(C k , C l ) = 1 n n i=1 n j=1 C k ij - 1 2
C k i. n + C k .j n C l ij . ( 10
)
LM (C k , C l ) = n i=1 n j=1 2C k ij C k i. + C k .j - 1 n C l ij . (11)
It is noteworthy that the different formulations of the statistical independence deviation with contingency tables in (4), ( 5), ( 6) and ( 7), translate into different types of deviation concepts in the relational representation ( 8), ( 9), ( 10) and [START_REF] Jordan | Les coefficients d'intensité relative de Korosy[END_REF]. Such properties were examined in [START_REF] Ah-Pine | Cluster analysis based on the central tendency deviation principle[END_REF] and led to the formalization of the central tendency deviation principle in cluster analysis. Indeed, one can observe the following central tendencies : in (9) C k .. /n 2 is the mean average over all the terms of C k ; in [START_REF] Belson | Matching and prediction on the principle of biological classification[END_REF]

(C k i. + C k .j )/2n is the arithmetic mean of C k i.
/n and C k .j /n and in (11) 1/n is the mean average over all terms of the matrix of general term

2C k ij /(C k i. + C k .j ) (which is equivalent to C k ij /C k i.
). Regarding [START_REF] Agarwal | Modularity-maximizing graph communities via mathematical programming[END_REF], the central tendency concept is of geometrical nature. Since C k is a dot product matrix (or Gram matrix) the transformation of

C k ij into C k ij -C k i. /n-C k .j /n+C k .
. /n 2 is known as the double centering (or Torgerson) transformation. This operation results in dots products between vectors centered with respect to the mean vector. Now that we have provided the expression of SAM using the graph relations underlying partitions, we show in the next paragraph how to employ such criteria for clustering graphs.

Graph clustering by maximizing SAM

We interpret a given undirected and unweighted graph as a perturbed equivalence relation and clustering the graph can be seen as attempting to recover the real partition. To solve the graph clustering task, we thus propose to fit the graph encoded by its adjacency matrix A with an equivalence relation by maximizing one of the SAM introduced previously. In other words, we want to find the partition that is the most similar to A according to a given SAM. More formally, we introduce the following benefit functions for clustering graphs :

B(A, X) = n i=1 n j=1 A ij - A i. n + A .j n - A .. n 2 X ij . (12) E(A, X) = n i=1 n j=1 A ij - A .. n 2 X ij . (13)
J(A, X) = 1 n n i=1 n j=1 A ij - 1 2
A i. n + A .j n X ij . (14)
LM (A, X) = n i=1 n j=1 2A ij A i. + A .j - 1 n X ij . ( 15
)
where X is the adjacency matrix of the partition we want to recover and whose general term is X ij = 1 if nodes i and j are in the same cluster and 0 otherwise. X represents an equivalence relation which, from an algebraic standpoint, is a binary relation with the following properties : (i) reflexivity (X ii = 1, ∀i); (ii) symmetry (X ij = 1 ⇔ X ji = 1, ∀i, j) and (iii) transitivity (X ij = 1 ∧ X jk = 1 ⇒ X ik = 1, ∀i, j, k). Marcotorchino and Michaud showed that these relational properties can be formulated as linear constraints through X [START_REF] Marcotorchino | Heuristic approach of the similarity aggregation problem[END_REF] and this finding allowed them to model the clustering problem as an integer linear program : max X ∆(A, X) subject to :

       X ii = 1 ∀i X ij -X ji = 0 ∀i, j X ij + X jk -X ik ≤ 1 ∀i, j, k X ij ∈ {0, 1} ∀i, j . (16)
where, in our case, ∆(A, X) is either ( 12) or ( 13) or ( 14) or [START_REF] Marcotorchino | Heuristic approach of the similarity aggregation problem[END_REF]. It is important to mention that this integer linear program allowed Marcotorchino to design the maximal association model for clustering data described by categorical variables. In Marcotorchino's works, A was considered as an equivalence relation [START_REF] Marcotorchino | Heuristic approach of the similarity aggregation problem[END_REF] or the sum over several equivalence relations [START_REF] Marcotorchino | Maximal association theory as a tool of research[END_REF]. Our proposal can thus be understood as the extension of the maximal association framework to graph clustering by considering A to be a general adjacency matrix without any particular property (except being undirected).

Before moving to the section dedicated to the experiments, we establish some interesting relationships between the modularity framework and our proposal based on SAM.

Some relationships between modularity and SAM

Firstly, using the notations introduced previously, the standard modularity defined in (2) can be reformulated as below :

Q(A, X) = 1 A .. n i=1 n j=1 A ij - A i. A .j A .. X ij . ( 17
)
In addition to the correspondence formulas given in paragraph 3.1, let us introduce the following identity : [START_REF] Marcotorchino | Utilisation des comparaisons par paires en statistique des contingences partie II[END_REF]. From this equation, if we identify C k and C l to A and X respectively and if we assume A and X to be partitions associated to two categorical variables V k and V l , we can easily show that the modularity given in ( 17) can be expressed by means of a contingency table as follows :

(v) v ( u N u. N uv ) 2 = i,j C k i. C k .j C l ij
Q(V k , V l ) = 1 p k u=1 N 2 u.   p k u=1 p l v=1 N 2 uv - 1 p k u=1 N 2 u.   p l v=1 p k u=1 N u. N uv 2     .
(18) Then, one can easily check that Q(V k , V l ) is null in case of statistical independence between V k and V l . This outcome shows the potential application of the modularity measure in categorical data analysis.

Let us now place the SAM in the context of the modularity concept developed by Newman. Let us formally introduce the following central tendencies :

µ Q ij = A i. A .j /A .. ; µ B ij = A i. /n+A .j /n-A .. /n 2 ; µ E ij = A .. /n 2 ; µ J ij = A i.
/(2n)+A .j /(2n) and µ LM ij = 1/n. In that case, ( 17), ( 12), ( 13), ( 14) and ( 15) can all be reformulated as :

α n i=1 n j=1 (A ij -µ Z ij ) with Z ∈ {Q, B, E, J, LM }, α = 1 when Z ∈ {B, E, LM }, α = 1/A .. when Z = Q, α = 1/n when Z = J and by substi- tuting A ij with A ij = 2A ij /(A i. + A .j
) when Z = LM . This expression of SAM better underlines the connections between the general modularity framework given in (1) and cluster analysis methods based on the central tendency deviation principle [START_REF] Ah-Pine | Cluster analysis based on the central tendency deviation principle[END_REF]. Furthermore, one can easily check that Newman's axioms we recalled in section 2 are both satisfied by all SAM under study except the LM method : (i) ∀i, j :

µ Z ij = µ Z ji for Z ∈ {B, E, J, LM }; (ii) i,j A ij = i,j µ Z ij for Z ∈ {B, E, J}.
As a result, B, E and J fit in the modularity model. However, hypothesis (iii) and (iv) are not satisfied by any SAM under study except B for which we have (iii) ∀i : j µ B ij = k i = A i. . In such a context, it is also interesting to notice that the SAM E given in [START_REF] Ah-Pine | Overview of the relational analysis approach in data-mining and multi-criteria decision making[END_REF] corresponds to another suggested modularity model which assumes a Bernoulli distribution for P ij in [START_REF] Papadopoulos | Community detection in social media[END_REF] and which boils down to the following constant 7 , P ij = A .. /n 2 , ∀i, j [START_REF] Newman | Finding community structure in networks using the eigenvectors of matrices[END_REF].

After having introduced the proposed objective functions and some properties about the relationships between the modularity concept and SAM, we examine in the next section if our proposals lead to interesting graph clustering methods from an empirical standpoint. Another goal of these experiments is to enable us to initiate a comparison between the modularity framework and SAM based optimization with regard to the hypothesis underlying each method.

Experiments

Our experiments are based on computer-generated graphs of different sizes. We give below the details about the algorithm we used to maximize the different density-based objective functions presented previously. We explain the tool we employed and the parameters we set to generate the artificial graphs. We then analyse the quality of the graph clustering results obtained with the different techniques.

Greedy optimization by agglomerative hierarchical clustering

The optimization problems given in ( 16) and ( 3) are NP-hard8 , and thus, numerous heuristics attempting to provide sub-optimal solutions have been proposed (see for e.g. the surveys [START_REF] Schaeffer | Graph clustering[END_REF][START_REF] Fortunato | Community detection in graphs[END_REF][START_REF] Papadopoulos | Community detection in social media[END_REF]). In our experiments, we used the greedy optimization algorithm proposed by Newman in [START_REF] Newman | Fast algorithm for detecting community structure in networks[END_REF] in order to maximize the modularity criterion given in [START_REF] Newman | Finding and evaluating community structure in networks[END_REF].

This heuristic is based on a simple agglomerative hierarchical clustering strategy. It starts with n distinct clusters and at each iteration it merges the two clusters that allow the best improvement of the modularity value. The merging process goes on until there is no pair of clusters whose fusion enables increasing the modularity value. This heuristic has the advantage to avoid fixing the number of clusters as a parameter.

This algorithm can be applied to other kinds of quality measures and in order to provide a fair comparison between the different objective functions, we thus used this technique to maximize ( 12), ( 13), ( 14) and (15) as well.

LFR benchmark graphs

The computer-generated graphs we analyzed in our experiments rely on the LFR benchmarks proposed by Lancichinetti, Fortunato and Radicchi in [START_REF] Lancichinetti | Benchmark graphs for testing community detection algorithms[END_REF][START_REF] Fortunato | Community detection algorithms: a comparative analysis: invited presentation[END_REF]. These benchmarks aim at providing the research community with graphs whose properties reflect real-world cases. Indeed, observed complex networks are characterized by heterogeneous distributions both for node degrees and cluster sizes. As a consequence, the authors developed a model that generates graphs which satisfy these features. They also implemented a freely available tool 9 that we used to conduct our empirical work.

Their approach is based on the planted l-partition model in which node degrees follow a power law distribution with exponent τ 1 and clusters size a power law distribution with exponent τ 2 . Overall the parameters of their model are : (i) n the number of nodes; (ii) the average degree of nodes; (iii) the maximum degree of nodes; (iv) τ 1 ; (v) τ 2 and (vi) µ ∈ [0, 1] the mixing parameter. The latter parameter µ is the one that allows gradually monitoring the presence or the absence of a community structure in the graph. It represents the percentage of edges that a node shares with vertices that do not belong to its group. Therefore, as µ grows, the community structure progressively degrades and the limit case µ = 1 corresponds to the situation where edges are totally placed randomly. Typically, we assume that there is a strong community structure within the graph as long as µ < 0.5.

Experiments settings and results

We studied graphs of different sizes : 500, 1000 and 2000 nodes. We used the same parameters values as in [START_REF] Fortunato | Community detection algorithms: a comparative analysis: invited presentation[END_REF] : the average degree was set to 20; the maximum degree to 50; τ 1 = 2; τ 2 = 1 and we vary µ from 0 to 0.7 with a 0.1 step.

The LFR benchmark graphs also provide a built-in community structure which allowed us to compare the clustering results with the real partition. To assess the proximity between the found partition and the ground-truth, the malized mutual information (N M I) measure was used. Let us denote by V k the clustering output of our algorithm with p k clusters and by V l the real clustering with p l groups. Let P (V k = u, V l = v) = P (u, v) = N uv /n be the probability of jointly observing u and v and let P (u) = N u. /n and P (v) = N .v /n be the probability of observing u and v respectively. The mutual information measure between V k and V l denoted by I(V k , V l ) is defined as follows : I(V k , V l ) = p k u=1 p l v=1 P (u, v) log(P (u, v)/(P (u)P (v))). Its normalized version denoted by

N M I(V k , V l ) is then given by N M I(V k , V l ) = 2I(V k , V l )/(H(V k ) + H(V l ))
where H(V k ) =u P (u) log(P (u)) is the entropy of V k . The N M I measure ranges from 0 and 1. It equals 1 when V k = V l and it is null when V k and V l are statistically independent. Note that we used the N M I coefficient to assess our clustering models because this measure is often used in the graph clustering literature (see for e.g. [START_REF] Fortunato | Community detection algorithms: a comparative analysis: invited presentation[END_REF]). In that way, we can also position our contributions with respect to other papers and graph clustering techniques.

To have a better estimation of the performances, we generated 5 different graphs for each distinct parameter setting and we took the median value. The experimental results obtained for N M I measures are shown in the first row of Fig. 1. We also computed the number of clusters found by each method in order to examine if the different techniques are able to recover the right number of clusters. In this case, we also took the median over the 5 trials. These results are shown in the second row of Fig. 1.

The first row of Figure 1 allows us to compare the quality of the clustering outputs found by the different methods. We claim the following outcomes : (i) as expected, all the methods have their quality diminishing as µ grows; (ii) the LM coefficient dominates the modularity Q and other methods and this superiority seems to grow with the size of the graphs; (ii) B and J perform similarly than Q whatever the size of the graphs; (iii) E is the less good approach and it particularly performs the worst when µ ∈ [0.3, 0.7].

Concerning the number of clusters found, we can make the following observations from the second row of Fig. 1 : (i) the number of real communities provided by the LFR benchmarks is stable and varies around 20, 40 and 80 for graphs of size 500, 1000, 2000 respectively; (ii) except for E with graphs of size 500, all the techniques produce less clusters than the correct number of groups; (iii) as µ grows the number of clusters decreases for all the methods and beyond a certain point (approximately 0.6 for Q, B, LM ; 0.5 for J and 0.3 for E), it tends to either grow or stabilize; (iv) the LM criterion tends to produce more clusters than Q, B and J.

Overall, B and J are comparable to Q while LM is clearly a better objective function than the other ones. One reason that could explain the superiority of LM is the fact that it implicitly transforms the binary matrix A into a non negative one, A whose general term is A ij = 2A ij /(A i. + A .j ). Then its related central tendency scheme, µ LM ij = 1/n, gives the same value for all pairs of nodes (i, j). Such an approach is indeed different from the other quality functions we examined, since they all keep the binary matrix but what changes from one function to the other is the underlying central tendency scheme µ Z , Z ∈ {Q, B, J, E} as we underlined in paragraph 3.3.

Moreover, our experimental results invite us to further analyze the hypothesis underlying the different objective functions. Regarding Newman's assumptions for the modularity given in [START_REF] Newman | Finding and evaluating community structure in networks[END_REF], it is interesting to notice that B and J perform similarly than Q despite the fact they do not satisfy all the hypothesis assumed by the latter criterion. More importantly, LM violates most of Newman's hypothesis but outperforms all other methods including Q.

Conclusion

We have proposed new objective functions for clustering undirected and unweighted graphs. Our method consists in maximizing SAM represented in their relational representation, in order to fit the given graph with a partition. Our empirical study on artificial graphs has shown encouraging results. Most of the proposed SAM perform equivalently or better than the modularity criterion. In particular, the Light-Margolin coefficient dominates the latter approach. As for future work, we plan to develop the analysis provided in paragraph 3.3 and leverage the empirical results presented previously by further comparing the modularity concept and SAM from a theoretical viewpoint. We also plan to extend our experiments on larger graphs both for computer-generated cases and for real-world networks in order to further validate our findings.

Fig. 1 .

 1 Fig. 1. First row : N M I values (vertical axis) versus mixing parameters µ (horizontal axis). Second row : Number of clusters found (vertical axis) versus mixing parameters µ (horizontal axis). We show the curves for each objective function. From left to right, the plots correspond to graphs with 500, 1000 and 2000 nodes respectively.

In order to lighten the notations we write u as a shortcut for p k u=1 .

It is actually an interpretation of Jordan's measure given by Marcotorchino in[START_REF] Marcotorchino | Utilisation des comparaisons par paires en statistique des contingences partie I[END_REF].

The study of association and aggregation of binary relations using graph theory and mathematical programming led to the Relational Analysis method developed by Marcotorchino and which has many applications in statistics, data-mining and multiple-criteria decision making (see for e.g.[START_REF] Ah-Pine | Overview of the relational analysis approach in data-mining and multi-criteria decision making[END_REF] and references therein).

Also called relational matrix in the Relational Analysis framework.

Note that in the case of E, we assume that the graph is reflexive unlike Q. In the latter case, the constant is A../(n(n -1)).

Note that the constraints in[START_REF] Marcotorchino | Maximal association theory as a tool of research[END_REF] and in (3) are equivalent : the former models the properties of an equivalence relation while the latter models the properties of the complementary of an equivalence relation which is a distance relation satisfying the triangle inequality constraint.

http://santo.fortunato.googlepages.com/inthepress2


Julien Ah-Pine
email: julien.ah-pine@univ-lyon2.fr

Brice Mayag
email: brice.mayag@dauphine.fr

Antoine Rolland
email: antoine.rolland@univ-lyon2.fr

Identification of a 2-additive bi-capacity by using mathematical programming

Keywords: 2-additive bi-capacity identification, Bipolar Choquet integral, Preference elicitation

In some multi-criteria decision making problems, it is more convenient to express the decision maker preferences in bipolar scales. In such cases, the bipolar Choquet integral with respect to bi-capacities was introduced. In this paper, we address the problem of eliciting a bipolar Choquet integral with respect to a 2-additive bi-capacity. We assume that we are given a set of examples with (i) their scores distribution in regard to several criteria and (ii) their overall scores. We propose two types of optimization problems that allow identifying the parameters of a 2-additive bi-capacity such that the inferred bipolar Choquet integral is consistent with the given examples as much as possible. Furthermore, since the elicitation process we study has many relationships with problems in statistical machine learning, we also present the links between our models and concepts developed in the latter field.

Introduction

Multi-criteria decision making (MCDM) aims at representing the preferences of a decision maker (DM) over a set of options (or alternatives) and in regard to several criteria. It then seeks to formalize the DM's decision process through mathematical tools in order to help him make decisions over the set of alternatives. The DM's decision process is assumed to be guided by the importance and the relationships he wants to take into account regarding the criteria. Concerning the preferences representations, one possible model is the Multi-Attribute Utility Theory (MAUT) which assumes that each attribute (or criterion) provides a utility value (or score) over the set of alternatives. Then, an aggregation function is used to combine, for each option, its scores distribution (or profile) in an overall score. The latter global utility values are then employed to make decisions. There are many types of aggregation functions to model a decision process. The Choquet integral has been proved to be a versatile tool to construct overall scores (see for example [START_REF] Choquet | Theory of capacities[END_REF][START_REF] Grabisch | k-order additive discrete fuzzy measures and their representation[END_REF][START_REF] Grabisch | Multiple Criteria Decision Analysis: State of the Art Surveys[END_REF][START_REF] Meyer | On the expressiveness of the additive value function and the choquet integral models[END_REF]). This aggregation function is intimately based on the concept of a capacity (or fuzzy measure). In particular, it assumes that partial utilities belong to non-negative or unipolar scales. Unipolar scales are not always appropriate to represent the DM's preferences (see the motivating example in [START_REF] Grabisch | A decade of application of the Choquet and Sugeno integrals in multi-criteria decision aid[END_REF]). In some problems, bipolar scales are more convenient. This type of scales is typically composed of a negative, a positive and a neutral part which respectively allow representing a negative, a positive and a neutral affect towards an option. To apply the Choquet integral in the case of bipolar scales, the bipolar Choquet integral (BCI) was introduced in [START_REF] Grabisch | Bi-capacities for decision making on bipolar scales[END_REF] and [START_REF] Greco | Bipolar Sugeno and Choquet integrals[END_REF]. In this paper, we particularly focus on BCI which use the concept of a bi-capacity (BC) introduced in [START_REF] Grabisch | Bi-capacities for decision making on bipolar scales[END_REF] and which was further studied in [START_REF] Grabisch | Bi-capacities-I: definition, Möbius transform and interaction[END_REF][START_REF] Grabisch | Bi-capacities-II: the Choquet integral[END_REF] and in [START_REF] Fujimoto | New characterizations of k-additivity and k-monotonicity of bicapacities[END_REF][START_REF] Fujimoto | Some characterizations of k-monotonicity through the bipolar möbius transform in bi-capacities[END_REF]. The BCI typically requires the DM to set 3 n -1 values where n is the number of attributes. When n exceeds some units, it is impossible for the DM to set all parameters of his decision model. In order to better cope with this combinatorial burden, the BCI with respect to (w.r.t.) a 2-additive bi-capacity (2A-BC) was introduced in the following papers [START_REF] Grabisch | The choquet integral for 2-additive bi-capacities[END_REF][START_REF] Fujimoto | New characterizations of k-additivity and k-monotonicity of bicapacities[END_REF]. The 2-additivity property implies that only the interactions between at most two criteria are taken into account in a BC, and it enables reducing the number of parameters from 3 n -1 to 2n 2 + 1. This has facilitated the use of this aggregation function in practice.

Even though there have been many papers studying 2A-BC, most of them have focused on theoretical aspects. In this contribution, we study the practical problem of identifying the parameters of a 2A-BC on the basis of information provided by the DM. This problem is also known as preference elicitation. There are different contexts in which we can proceed to the elicitation of the preference model of a DM. In our case, we assume that the DM provides the bipolar scores for a subset of (real or fictitious) options w.r.t. all criteria of the decision problem. In addition, he provides the overall bipolar scores of the same set of alternatives. These evaluated examples constitute the only data we have at our disposal. Then, the elicitation model consists in inferring the parameters of a 2A-BC such that the associated BCI is consistent with the preferences given by the DM on these examples. We propose optimization models that address this kind of preference elicitation problems.

Eliciting preference models is a research topic that has been studied by many researchers (see for example [START_REF] Jacquet-Lagreze | Assessing a set of additive utility functions for multicriteria decision-making, the UTA method[END_REF][START_REF] Grabisch | A review of methods for capacity identification in Choquet integral based multi-attribute utility theory: Applications of the kappalab R package[END_REF]). However, the BCI w.r.t. 2A-BC has not been studied very much so far. To our knowledge, the only paper that addressed this exact problem is [START_REF] Mayag | Elicitation of a 2-additive bi-capacity through cardinal information on trinary actions[END_REF]. Yet, in the latter paper, the authors assumed an elicitation process in which the DM was asked to provide cardinal information on trinary actions. This setting is different from the one considered in this paper. Our approach is also in line with the work detailed in [START_REF] Grabisch | A review of methods for capacity identification in Choquet integral based multi-attribute utility theory: Applications of the kappalab R package[END_REF] 1 about the preference elicitation using unipolar Choquet integrals with mathematical programming. In this latter work, the authors also assume that the only information provided are examples evaluated by the DM.

The elicitation process we deal with has many relationships with the problems addressed in statistical machine learning. The interconnections between preference elicitation on the one hand and machine learning on the other hand were highlighted in [START_REF] Waegeman | Kernel-based learning methods for preference aggregation[END_REF]. There has been a growing interest for the last years about cross-fertilizing these two domains by studying how the concepts devel-1 And with the papers cited therein.

oped in one field can be applied in the other one. In line with this research topic, we also discuss the links between the concepts of these two domains that our models involve.

The rest of the paper is organized as follows. We recall in section 2 some basic definitions about BC and the properties of 2A-BC by using the bipolar Möbius transform defined by [START_REF] Fujimoto | New characterizations of k-additivity and k-monotonicity of bicapacities[END_REF]. We then propose in section 3 two identification methods of a 2A-BC using linear programming and quadratic programming. In order to illustrate our proposals, we apply the different methods to a numerical example. Next, in section 4, we underline the relationships of our approaches with concepts developed in the field of machine learning. We finally conclude this paper and sketch some future works in section 5.

Bi-capacities and bipolar Choquet integrals

Let us denote by N = {1, . . . , n} a finite set of n criteria and X = X 1 × • • • × X n the set of possible alternatives, where X 1 , . . . , X n represent the attributes. For all i ∈ N , the function u i : X i → R is called a utility function. Given an element x = (x 1 , . . . , x n ), we denote by U (x) = (u 1 (x 1 ), . . . , u n (x n )), the element's profile or its scores distribution. We will often write ij, ijk instead of {i, j} and {i, j, k} respectively.

2-additive bi-capacities

Let us denote by 2 N := {S ⊆ N } the set of subsets of N and 3 N := {(A, B) ∈ 2 N ×2 N : A∩B = ∅} the set of couples of subsets of N with an empty intersection. We define on 3 N the following relation ,

∀(A 1 , A 2 ), (B 1 , B 2 ) ∈ 3 N : (A 1 , A 2 ) (B 1 , B 2 ) ⇔ [A 1 ⊆ B 1 and B 2 ⊆ A 2 ]
Definition 1 (Bi-capacity (BC) [START_REF] Grabisch | Bi-capacities-II: the Choquet integral[END_REF], [START_REF] Grabisch | A decade of application of the Choquet and Sugeno integrals in multi-criteria decision aid[END_REF]) A function ν : 3 N → R is a BC on 3 N if it satisfies the following two conditions :

ν(∅, ∅) = 0 (1) ∀(A 1 , A 2 ), (B 1 , B 2 ) ∈ 3 N : [(A 1 , A 2 ) (B 1 , B 2 ) ⇒ ν(A 1 , A 2 ) ≤ ν(B 1 , B 2 )] (2)
Note that (2) is called the monotonicity condition.

In addition, a BC is said to be normalized if it satisfies :

ν(N, ∅) = 1 and ν(∅, N ) = -1 (3)
A BC is also said to be additive if the following relation holds :

∀(A 1 , A 2 ) ∈ 3 N : ν(A 1 , A 2 ) = i∈A1 ν(i, ∅) + j∈A2 ν(∅, j) (4)
An additive BC assumes that the attributes are independent from each other and this kind of BC boils down to linear decision models.

In order to better formalize some of the properties of BC, the following definition of a (bipolar) Möbius transform2 of a BC was proposed.

Definition 2 (Bipolar Möbius transform of a bi-capacity [START_REF] Fujimoto | New characterizations of k-additivity and k-monotonicity of bicapacities[END_REF][START_REF] Fujimoto | k-additivity and c-decomposability of bi-capacities and its integral[END_REF]) Let ν be a BC on 3 N . The bipolar Möbius transform of ν is a set function b : 3 N → R defined for any (A 1 , A 2 ) ∈ 3 N by :

b(A 1 , A 2 ) := B 1 ⊆ A 1 B 2 ⊆ A 2 (-1) |A1\B1|+|A2\B2| ν(B 1 , B 2 ) (5) = (∅,A2) (B1,B2) (A1,∅) (-1) |A1\B1|+|A2\B2| ν(B 1 , B 2 )
Conversely, for any (A 1 , A 2 ) ∈ 3 N , it holds that :

ν(A 1 , A 2 ) := B 1 ⊆ A 1 B 2 ⊆ A 2 b(B 1 , B 2 ). (6)
Note that using b, (1) is equivalent to :

b(∅, ∅) = 0 (7)
BC on 3 N generally require 3 n -1 parameters. In order to reduce this number, [START_REF] Grabisch | Bi-capacities-I: definition, Möbius transform and interaction[END_REF][START_REF] Grabisch | Bi-capacities-II: the Choquet integral[END_REF] and [START_REF] Grabisch | A decade of application of the Choquet and Sugeno integrals in multi-criteria decision aid[END_REF] proposed the concept of k-additivity of a BC. This concept translates as follows in terms of the bipolar Möbius transform.

Proposition 1 ( [START_REF] Fujimoto | k-additivity and c-decomposability of bi-capacities and its integral[END_REF]) Given a positive integer k < n, a BC ν is k-additive if and only if the two following conditions are satisfied :

∀(A 1 , A 2 ) ∈ 3 N : |A 1 ∪ A 2 | > k ⇒ b(A 1 , A 2 ) = 0 (8) ∃(A 1 , A 2 ) ∈ 3 N : |A 1 ∪ A 2 | = k ∧ b(A 1 , A 2 ) = 0 (9)
To avoid a heavy notation, we use the following shorthands for all i, j ∈ N , i = j :

• ν i| := ν(i, ∅), ν |j := ν(∅, j), ν i|j := ν(i, j), ν ij| := ν(ij, ∅), ν |ij := ν(∅, ij), • b i| := b(i, ∅), b |j := b(∅, j), b i|j := b(i, j), b ij| := b(ij, ∅), b |ij := b(∅, ij).
Whenever we use i and j together, it always means that they are different.

Using the above definitions, we propose the following properties of a 2A-BC ν and its bipolar Möbius transform b : Proposition 2 1. Let ν be a 2A-BC and b its bipolar Möbius transform. For any (A 1 , A 2 ) ∈ 3 N we have :

ν(A 1 , A 2 ) = i∈A1 b i| + j∈A2 b |j + i ∈ A 1 j ∈ A 2 b i|j + {i,j}⊆A1 b ij| + {i,j}⊆A2 b |ij (10) 2. If the coefficients b i| , b |j , b i|j , b ij| , b |ij are
given for all i, j ∈ N , then the necessary and sufficient conditions to get a 2A-BC generated by [START_REF] Fujimoto | New characterizations of k-additivity and k-monotonicity of bicapacities[END_REF] are :

∀(A, B) ∈ 3 N , ∀k ∈ A : b k| + j∈B b k|j + i∈A\k b ik| ≥ 0 (11) ∀(A, B) ∈ 3 N , ∀k ∈ A : b |k + j∈B b j|k + i∈A\k b |ik ≤ 0 ( 12
)
3. The inequalities ( 11) and ( 12) can be respectively reformulated in terms of the BC ν as follows :

∀(A, B) ∈ 3 N , ∀k ∈ A : j∈B ν k|j + i∈A\k ν ik| ≥ (|B|+|A|-2)ν k| + j∈B ν |j + i∈A\k ν i| ∀(A, B) ∈ 3 N , ∀k ∈ A : j∈B ν j|k + i∈A\k ν |ik ≤ (|B|+|A|-2)ν |k + j∈B ν j| + i∈A\k ν |i Proof. (Sketch of)
1. Because ν is 2-additive, the proof of ( 10) is given by using the relation ( 6) between ν and b. 2. The proof of the second point is based on the expression of ν(A 1 , A 2 ) given in [START_REF] Fujimoto | New characterizations of k-additivity and k-monotonicity of bicapacities[END_REF] and on these equivalent monotonicity properties (which are easy to check) :

∀(A, B) ∈ 3 N and ∀A ⊆ A , (a) ν(A, B) ≤ ν(A , B) ⇔ {∀k ∈ A : ν(A \ k, B) ≤ ν(A, B)}; (b) ν(B, A ) ≤ ν(B, A) ⇔ {∀k ∈ A : ν(B, A) ≤ ν(B, A \ k)}.
3. These inequalities are obtained departing from [START_REF] Fujimoto | Some characterizations of k-monotonicity through the bipolar möbius transform in bi-capacities[END_REF] and [START_REF] Grabisch | The choquet integral for 2-additive bi-capacities[END_REF] and by using the relation ( 6) between ν and b.

Hence, according to proposition 2 and (10), the computation of a 2A-BC ν only requires the values of b on the elements (i, ∅), (∅, i), (i, j), (ij, ∅), (∅, ij), ∀i, j ∈ N . However, in order to satisfy the monotonicity condition given in (2) a 2A-BC should also satisfy the inequalities [START_REF] Fujimoto | Some characterizations of k-monotonicity through the bipolar möbius transform in bi-capacities[END_REF] and [START_REF] Grabisch | The choquet integral for 2-additive bi-capacities[END_REF]. Moreover, we have the following conditions in order to obtain a normalized 2A-BC :

ν N | = i∈N b i| + {i,j}⊆N b ij| = 1 and ν |N = i∈N b |i + {i,j}⊆N b |ij = -1 (13)
C ν (x) := n i=1 |x σ(i) | ν(N σ(i) ∩N + , N σ(i) ∩N -)-ν(N σ(i+1) ∩N + , N σ(i+1) ∩N -) (14) where N + = {i ∈ N |x i ≥ 0}, N -= N \ N + , N σ(i) := {σ(i), . . . , σ(n)} and σ is a permutation on N such that |x σ(i) | ≤ |x σ(i+1) | ≤ . . . ≤ |x σ(n) |.
We also have the following equivalent expression of the BCI w.r.t. b, given by [START_REF] Fujimoto | Some characterizations of k-monotonicity through the bipolar möbius transform in bi-capacities[END_REF] :

C b (x) = (A1,A2)∈3 N b(A 1 , A 2 ) i∈A1 x + i ∧ j∈A2 x - j (15)
where

x + i = x i if x i > 0 x + i = 0 if x i ≤ 0 and x - i = -x i if x i < 0 x - i = 0 if x i ≥ 0 . Note that C ν (x) = C b (x)
and the subscript is meant to clarify whether it is ν or b which is used in the calculation. Besides, the BCI of x w.r.t. a 2A-BC represented by b reduces to :

C b (x) = n i=1 b i| x + i + n i=1 b |i x - i + n i,j=1 b i|j (x + i ∧ x - j ) (16)
+ {i,j}⊆N b ij| (x + i ∧ x + j ) + {i,j}⊆N b |ij (x - i ∧ x - j )
We have introduced the basic tools related to 2A-BC. In the next section, we focus on the problem of identifying a 2A-BC.

Identifying a 2-additive bi-capacity

We establish mathematical programming problems that enable the identification of a 2A-BC. First, we detail the type of elicitation process we are concerned with. Next, we state all the constraints that allow the representation of a 2A-BC. Then, we define objective functions that reflect the quality of the identified 2A-BC in regard to the information provided by the DM. We end this section by illustrating the results obtained with the proposed models on a numerical example.

Elicitation process

In MCDM, there are two types of paradigms for elicitation processes : direct and indirect methods. In the former case, the DM is able to provide the parameters of his decision model directly. However, when using a BCI w.r.t. a BC, the direct method is infeasible if the number of attributes n exceeds some units (typically 4), since a BC requires 3 n -1 values to be set. We argued in the introduction that in order to reduce this complexity, 2A-BC were introduced. Nevertheless, this latter case cannot be applied in practice neither since, even if the number of parameters reduces to 2n 2 + 1, the complexity remains very high. Moreover, a BC is a too complex aggregation operator to ensure that a DM will understand the influence of each parameter on the final result. Even with more simple aggregation rules, it has been shown that there is no clear link between the parameters values provided by the DM and the way these values are used in the decision model [START_REF] Bouyssou | Evaluation and Decision models with multiple criteria: stepping stones for the analyst[END_REF].

Therefore, we follow the indirect paradigm. In that case, the DM does not give information about his decision model, instead, he provides information on the outputs of his decision model. In our setting, we suppose that the DM gives for some examples x ∈ X ⊆ X, their partial utilities for all criteria (U (x)) and also their overall scores (S(x)). We then assume, that there is no further interaction with the DM. Given the judged examples, we have to infer a decision model based on the BCI w.r.t. a 2A-BC. The estimated BCI should predict overall scores, C b (x), that are consistent with the preference relations provided by the DM. In other words, if S(x) ≥ S(x ), which means that x is preferred or equivalent to x , then the inferred decision model should also satisfy

C b (x) ≥ C b (x ).
However, it might happen that this condition is not fulfilled for some pairs (x, x ). There are two main reasons for such situations : either the judgements provided by the DM himself are not consistent or the restriction of the decision model to 2A-BC does not allow fitting the DM preferences correctly. In MCDM, inconsistencies are usually treated in an interaction loop with the DM. It is assumed that the DM preferences can change in order to fix these incoherences when they are encountered. In our setting, the interaction loop is not permitted. Consequently, in order to cope with incoherences, we propose two versions of our models : the first one does not deal with inconsistencies and thus will return that the problem is infeasible if any incoherence is encountered whereas the second one allows inconsistencies and attempts to infer a model that minimizes errors due to such situations as much as possible.

Mathematical programming problems

We propose two types of optimization problems to identify a 2A-BC in the context we have described in the previous paragraph. We base our work on some of the elicitation methods detailed in [START_REF] Grabisch | A review of methods for capacity identification in Choquet integral based multi-attribute utility theory: Applications of the kappalab R package[END_REF] in the case of unipolar Choquet integral. Before introducing the objective functions of our optimization problems, we start by enumerating the different sets of constraints that need to be satisfied.

We represent the unknown 2A-BC ν, via its associated bipolar Möbius transform b. There are two reasons for this. Firstly, equations ( 5) and ( 6) state that it is equivalent to work with either ν or b. Secondly, since we restrict the BC to be 2-additive and since this property, given in ( 8) and ( 9), is defined in terms of b, it is thus necessary to use the latter representation in our optimization problems. However, we do not take (9) into account in our set of constraints. This equation ensures that ν is exactly 2-additive and by discarding it, we explicitly allow b to be either 2-additive or simply additive. To summarize this first set of constraints, we need to integrate the following relations in our optimization problems in order to have a normalized 2A-BC in terms of b : ( 7), ( 8) with k = 2, (11), ( 12) and [START_REF] Jacquet-Lagreze | Assessing a set of additive utility functions for multicriteria decision-making, the UTA method[END_REF].

Next, we have to take into account the preference relations provided by the DM on the subset of examples X . If S(x) ≥ S(x ) then the BCI should be in concordance with this inequality. Accordingly, we have the following second set of constraints :

∀x, x ∈ X , x = x : S(x) -S(x ) ≥ 0 ⇒ C b (x) -C b (x ) ≥ δ c ( 17
)
where δ c is a non-negative indifference threshold which is a parameter of the model. Note that this set of constraints does not allow incoherences. Indeed, the inferred 2A-BC b could not be flexible enough to satisfy C b (x) -C b (x ) ≥ δ c for some pairs (x, x ). In that case the optimization problem is infeasible.

As discussed previously, in order to overcome this drawback, we transform the previous constraints as follows :

∀x, x ∈ X , x = x : S(x) -S(x ) ≥ 0 ⇒ C b (x) -C b (x ) ≥ δ c -ξ xx (18)
where ξ xx are non-negative slack variables which allow inconsistencies. However, we want ξ xx to be has low as possible and thus there should be a term in the objective function seeking to minimize x,x :S(x)≥S(x ) ξ xx . Note that when the latter term is null, it means that the inferred model does not produce any incoherence. On the contrary, if for some pairs (x, x ), ξ xx > δ c then the optimal solution has not been able to satisfy the preference relations on these pairs. The third set of constraints is related to the computation of the BCI. Indeed, in [START_REF] Fujimoto | k-additivity and c-decomposability of bi-capacities and its integral[END_REF] or [START_REF] Bouyssou | Evaluation and Decision models with multiple criteria: stepping stones for the analyst[END_REF], we need to calculate C b (x) for each x ∈ X . As a consequence, we need to add the constraints provided by [START_REF] Waegeman | Kernel-based learning methods for preference aggregation[END_REF] in our models. Note that despite the fact that the latter equations involve the minimum function, we can pre-compute the terms (x ± i ∧ x ± j ) since they are parameters of the models. Consequently, the constraints ( 16) are linear equations.

The fourth set of constraints is optional. It simply consists in adding upper and lower bounds for the BCI values :

∀x ∈ X : lb ≤ C b (x) ≤ ub (19)
where lb and ub are two real parameters.

After having introduced the constraints, we now focus on the different objective functions and the resulting optimization problems.

We propose two kinds of optimization models. In the first approach, we extend the maximum split method introduced in [START_REF] Marichal | Determination of weights of interacting criteria from a reference set[END_REF]. This model assumes the following constraints in place of (17) :

∀x, x ∈ X , x = x : S(x) -S(x ) ≥ 0 ⇒ C b (x) -C b (x ) ≥ δ c + ε ( 20
)
where ε is a variable of the problem unlike δ c which is a parameter. The objective function consists in maximizing ε. In other words, we want to maximize the difference (split) C ν (x) -C ν (x ) for any x = x ∈ X such that S(x) ≥ S(x ). We refer to the following optimization problem as the split method : max ε subject to ( 7), ( 8) with k = 2, ( 11), ( 12), ( 13), ( 16), ( 19) and [START_REF] Cortes | Support-vector networks[END_REF].

However, the split model does not address incoherences. Hence, as explained previously, we propose the split f lex approach which uses the following third set of constraints instead of (20) :

∀x, x ∈ X , x = x : S(x) -S(x ) ≥ 0 ⇒ C b (x) -C b (x ) ≥ δ c + ε -ξ xx (21)
More formally, the split f lex model is defined by : max εx,x :S(x)≥S(x ) ξ xx subject to ( 7), ( 8) with k = 2, ( 11), ( 12), ( 13), ( 16), ( 19) and [START_REF] Siskos | UTA methods. In: Multiple Criteria Decision Analysis: State of the Art Surveys[END_REF]. Note that the split and split f lex optimization problems have linear objective functions and linear constraints. Therefore, there are linear programs.

We now present the second type of model for 2A-BC identification. This approach is a regression-like method and yields to quadratic programs. We propose to minimize the sum of square errors between S and C b which results in the following objective function : min x,x ∈X (S(x) -C b (x)) 2 . Accordingly, we named rss (for Residual Sum of Square) the following problem : min x,x ∈X (S(x) -C b (x)) 2 subject to ( 7), ( 8) with k = 2, ( 11), ( 12), ( 13), ( 16) and [START_REF] Fujimoto | k-additivity and c-decomposability of bi-capacities and its integral[END_REF]. Similarly to split, the rss method does not permit incoherences.

As a consequence, we introduce a flexible version of rss that we call rss f lex and which is given by : min x,x ∈X (S(x) -C b (x)) 2 + x,x :S(x)≥S(x ) ξ xx subject to ( 7), ( 8) with k = 2, ( 11), ( 12), ( 13), ( 16) and ( 18).

An illustrative example

We applied the four different mathematical programming problems defined previously on a numerical example taken from [START_REF] Grabisch | A review of methods for capacity identification in Choquet integral based multi-attribute utility theory: Applications of the kappalab R package[END_REF]. It concerns the grades (utilities) obtained by 7 students (alternatives) for n = 5 subjects (attributes) : statistics (S), probability (P), economics (E), management (M), and English (En). The grades globally belong to [0, 20] but in this example, the scores only vary in [START_REF] Fujimoto | Some characterizations of k-monotonicity through the bipolar möbius transform in bi-capacities[END_REF][START_REF] Bouyssou | Evaluation and Decision models with multiple criteria: stepping stones for the analyst[END_REF]. In our perspective, we transformed them in order to have a bipolar scale by simply applying a translation of -14 to the original scores. Therefore, in this bipolar scale, the scores belong to [-3, 4]. Suppose that a student is delivered his diploma with honors providing that his overall grade is greater or equal to 14. Hence, the translated scores in the bipolar scale allow us to deal with the decision problem of delivering honors as follows : the student is attributed the honors if and only if his overall grade in the bipolar scale is non-negative.

The performance table is given in Table 1 (a). In Table 1 (b), the first column S corresponds to the (translated) overall grades as given in [START_REF] Grabisch | A review of methods for capacity identification in Choquet integral based multi-attribute utility theory: Applications of the kappalab R package[END_REF]. Then, in the subsequent columns of Table 1 (b), we show the different estimated scores. Note that for all models we set δ c = 0.5, lb = -3 and ub = 4. For split and split f lex, even if the inferred overall scores are not the same, the two solutions are actually equivalent since they give the same objective function value (the problem is not strictly concave). However, in terms of decisions, the sign of the overall grade for c is not the same for the two models. This example exhibits some limits of this type of model. Besides, we precise that in split f lex outputs, all slack variables ξ xx are null which means that there is no inconsistency. Regarding rss and rss f lex, we obtain a null objective function value and in the latter case, the slack variables are also all null as expected.

Student S P E M En a 4 -3 -3 -3 4 b 4 -3 4 -3 -3 c -3 -3 4 -3 4 d 4 4 -3 -3 -3 e -3 -3 4 4 -3 f -3 -3 4 -3 -3 g -3 -3 -3 -3 4 S 1 0.5 0 -0.5 -1 -1.5 -2 split 1.
To illustrate the case with incoherences, we modified the overall score S into S . We simply change the global grade of g from -2 to 0.5. This new score is an example of inconsistent preferences provided by the DM since his decision model is not monotonic in that case. Indeed, if we compare the profiles of c and g in regard to their overall score S (c) and S (g), we can observe that g is preferred to c while the scores distribution of the former student is Pareto dominated by the latter one. This situation is not consistent with a rational decision. Table 1 (c) presents the estimated BCI. As expected, split and rss returned an infeasible problem. On the contrary, split f lex and rss f lex provide interesting results. In both cases, we precise that ξ gc = 0.5 while the other slack variables are null.

We present in Table 2 the elicited bipolar Möbius transform b for each optimization problem when S is the targeted overall grades vector. In Table 3, it is the estimated b for the split f lex and rss f lex methods when S is put in place of S, which are shown. Notice that an empty cell in these two tables means that the solver returned a null value for the corresponding elements.

We show the bipolar Möbius transform b and not the 2A-BC ν because the latter set function requires 3 5 -1 = 242 non-null values which represents a too large table. In contrast, because ν is 2-additive, we have at most 2 × 5 2 + 1 = 51 values for its related bipolar Möbius transform b.

The utility of presenting Tables 2 and3, is that they allow one to check that the constraints [START_REF] Greco | Bipolar Sugeno and Choquet integrals[END_REF], [START_REF] Grabisch | Bi-capacities-I: definition, Möbius transform and interaction[END_REF] with k = 2, (11), ( 12), ( 13) are indeed satisfied. We thus obtain a 2A-BC. However, it is difficult to interpret the bipolar Möbius transform b with regard to the underlying elicited preference model. Moreover, from this illustrative example, it is not straightforward to understand the impact of taking into account inconsistencies either when we compare the regular and the f lex versions of split and rss models in Table 2, or when we look at the obtained b for S in Table 2 and the one estimated for S in Table 3. Accordingly, further experiments should be undertaken in order to have a better understanding of the behaviors of the elicited ν and b for each proposed optimization problem but such an experimental study is out of the scope of this paper.

split split f lex rss rss f lex A1 A2 b A 1 |A 2 b A 2 |A 1 b A 1 |A 2 b A 2 |A 1 b A 1 |A 2 b A 2 |A 1 b A 1 |A 2 b A 2 |A 1 ∅ ∅ ∅ S 0.

Relationships with machine learning

The preference elicitation problem have many common points with supervised learning problems (SL) in statistical machine learning (ML). In the latter field, we are given a training set which consists of items described in a feature space and each element also comes with a value in regard to a target variable. In ML, the goal of SL is to infer from the training set a mapping from the feature space to the target variable. This description if similar to the preference elicitation setting we have described previously : the set of examples with their profiles U (x) and their overall scores S(x) are the equivalent of the training set in ML and identifying the parameters of the decision model in order to reproduce the preference relations given by the DM is the same as inferring a mapping from a feature space (X) and a target variable (S). Despite these straightforward similarities, the roots of these two domains are distinct, and it is only recently, that there has been a growing interest in applying concepts or/and techniques developed in MCDM to ML problems and vice-versa [START_REF] Waegeman | Kernel-based learning methods for preference aggregation[END_REF].

split f lex rss f lex A1 A2 b A 1 |A 2 b A 2 |A 1 b A 1 |A 2 b A 2 |A 1 ∅ S -
To contribute in that direction, we discuss the relationships between our work and some concepts defined in ML.

Firstly, in MCDM, incoherences are not generally allowed from the DM viewpoint who has to fix them during the course of the elicitation procedure. In ML, on the contrary, such situations are typically observed in real-world applications. Thus, the models split and rss are typical of MCDM whereas split f lex and rss f lex are closer to SL problems.

Secondly, it is noteworthy that the restriction to 2A-BC can be interpreted as an explicit regularization of the decision model based on the BCI. Indeed, in SL, regularization is a concept that aims at dealing with the bias-variance trade-off of predictive models. Typically, in order to avoid an over-fitting effect, we allow the estimated predictive model to make more errors on the training set (bias) but in return we want it to be less variable (variance). The goal of regularization is to enhance the ability of the inferred model to predict correctly the overall score of observed examples (the training set) but also and in particular the global score of unseen examples (the test set). In our case, 2A-BC are less flexible than unconstrained BC and this could lead to more inconsistencies as explained beforehand. However, the 2-additivity property makes the BC less complex and we intuitively expect 2A-BC to be less variable than unconstrained BC.

Thirdly, we stated that the constraint (9) ensuring that ν is exactly 2-additive was not part of our set of constraints. In that case ν could thus be either 2additive or simply additive. We can make the correspondence between this approach and the Occam's Razor principle often used in ML. This concept states that one should prefer simpler models than more complex ones because they allow a better understanding of the phenomenon under study. In our case, we can transpose this statement as follows : if an additive BC fits better the DM judgements than a 2A-BC then we should go for the former one. By discarding (9) from our constraints we make the latter statement possible. However, if it happens that a 2A-BC and an additive BC yield to the same optimal objective function value then it is not guaranteed that the optimization solver will provide the simplest solution.

Finally, in order to make the models split and rss more flexible regarding inconsistencies, we have proposed to integrate slack variables in the constraints and in the objective function as well. Our approach is inspired from the Support Vector Machine (SVM) method developed in SL in order to deal with non linearly separable cases in binary classification [START_REF] Cortes | Support-vector networks[END_REF]. However, it is noteworthy that the UTA (UTilités Additives) framework is a MCDM methodology that also addresses inconsistencies by integrating overestimation and underestimation error variables in the elicitation model (see for example [START_REF] Siskos | UTA methods. In: Multiple Criteria Decision Analysis: State of the Art Surveys[END_REF]). This approach is similar to adding slack variables.

Conclusion

We have proposed optimization problems that allow the identification of the parameters of a 2A-BC. The decision model is inferred from examples that the DM evaluated both regarding their partial utilities and their global scores. We have considered the traditional preference elicitation setting where no inconsistency is allowed. But, we have also extended the models to the more flexible case where we have to cope with such incoherences. In this context, we have emphasized the relationships between our approaches and concepts in ML.

In our future work, we intend to integrate other kinds of information provided by the DM such as the importance of criteria and the interactions between them.

Moreover, as pointed out previously, additional experiments should be conducted in order to better characterize the elicited bipolar Möbius transform b and its associated 2A-BC ν, we obtain for each type of mathematical program.

As regard to our ongoing work, we are investigating optimization problems for preferences learning that integrate other concepts developed in the ML literature. We are currently working on objective functions that involve a penalty term that favors sparse BC. Such an approach adds a trade-off between the accuracy of the model and its simplicity (Occam's Razor principle). The expected advantage is to elicit preference models that are easier to interpret and such a feature is of great importance in MCDM.

2. 2

 2 Bipolar Choquet integral w.r.t. a 2-additive bi-capacity Definition 3 (Bipolar Choquet integral (BCI) (w.r.t. a BC) [9]) Let ν be a BC on 3 N and x = (x 1 , . . . , x n ) ∈ R n . The expression of the BCI of x w.r.t. ν is given by

Table 1 .

 1 (a) Performance table; (b) Results obtained with the original (translated) overall score S; (c) Results obtained with the modified overall score S that presents inconsistencies.

		split		rss			split		rss
		f lex	rss	f lex	S	split	f lex	rss	f lex
	68	1.02	1	1	1	.	0.22	.	1.12
	1.04	0.38	0.5	0.5	0.5	.	-0.28	.	0.62
	0.41	-0.25	0	0	0	.	-0.78	.	0.12
	-0.23	-0.89	-0.5	-0.5	-0.5	.	-1.28	.	-0.5
	-0.86	-1.53	-1	-1	-1	.	-1.78	.	-1
	-1.5	-2.16	-1.5	-1.5	-1.5	.	-2.28	.	-1.5
	-2.14	-2.8	-2	-2	0.5	.	-0.78	.	0.12
	(a)	(b)					(c)

Table 2 .

 2 Values of the elicited bipolar Möbius transform b when S is the overall score.

		26 -0.145	-0.039 0.077	0.056
	∅ P		-0.077 0.038 -0.056
	∅ E		0.145 -0.038 0.077	0.056
	∅ M		-0.08 0.038
	∅ En	-0.066	-0.039 0.008 -0.056
	∅ SP	-0.56	-0.33	-0.24
	∅ SE -0.076		0.11 -0.004 0.11	0.19
	∅ SM -0.21		-0.008
	∅ SEn -0.29 0.48	0.43 0.039 0.47	0.53
	∅ PE		-0.29	-0.43
	∅ PM		-0.03
	∅ PEn		-0.008
	∅ EM -0.42	-0.23	0.004
	∅ EEn	0.26	0.32 -0.14 0.23 -0.056 0.17
	∅ MEn			-0.17
	S P
	S E		-0.07 -0.039	-0.056
	S M		0.145 0.004	-0.02
	S En		0.039 0.039 0.056
	P E
	P M		0.039
	P En		0.038	0.056
	E M		0.004
	E En	-0.145	-0.039 0.034
	M En	0.066	0.073

Table 3 .

 3 Values of the elicited bipolar Möbius transform b when S is the overall score with incoherences.

		0.077 0.47
	∅ P -0.19	-0.53 0.077
	∅ E	-0.077 0.023
	∅ M	0.031 -0.14 0.077
	∅ En -0.031		0.19
	∅ SP -0.26	0.077 0.09
	∅ SE	0.19	0.263
	∅ SM	0.73
	∅ SEn -0.16 0.055	-0.018
	∅ PE
	∅ PM	0.14
	∅ PEn 0.031		-0.077
	∅ EM
	∅ EEn -0.24	-0.39 -0.023
	∅ MEn -0.16		-0.077
	S P	0.06
	S E
	S M	-0.031 -0.33
	S En
	P E
	P M		0.14
	P En	0.19 -0.077 0.31
	E M
	E En		0.077
	M En	-0.077

Note that[START_REF] Grabisch | The choquet integral for 2-additive bi-capacities[END_REF] was the first paper to define the Möbius transform of a BC. Their definition is different from the one given in[START_REF] Fujimoto | New characterizations of k-additivity and k-monotonicity of bicapacities[END_REF]. However, there is a one-to-one correspondence between the two Möbius transform definitions. This equivalence was established in[START_REF] Fujimoto | Some characterizations of k-monotonicity through the bipolar möbius transform in bi-capacities[END_REF].