Dino Ienco
email: ienco@di.unito.it

Céline Robardet
email: celine.robardet@insa-lyon.fr

Ruggero G Pensa

Rosa Meo
email: meo@di.unito.it

Parameter-Less Co-Clustering for Star-Structured Heterogeneous Data

The availability of data represented with multiple features coming from heterogeneous domains is getting more and more common in real world applications. Such data represent objects of a certain type, connected to other types of data, the features, so that the overall data schema forms a star structure of inter-relationships. Co-clustering these data involves the specification of many parameters, such as the number of clusters for the object dimension and for all the features domains. In this paper we present a novel co-clustering algorithm for heterogeneous star-structured data that is parameter-less. This means that it does not require either the number of row clusters or the number of column clusters for the given feature spaces. Our approach optimizes the Goodman-Kruskal's τ , a measure for cross-association in contingency tables that evaluates the strength of the relationship between two categorical variables. We extend τ to evaluate co-clustering solutions and in particular we apply it in a higher dimensional setting. We propose the algorithm CoStar which optimizes τ by a local search approach. We assess the performance of CoStar on publicly available datasets from the textual and image domains using objective external criteria. The results show that our approach outperforms state-of-the-art methods for the co-clustering of heterogeneous data, while it remains computationally efficient.

Introduction

Clustering is by far one of the most popular techniques among researchers and practitioners interested in the analysis of data from different sources. It is often employed to obtain a first schematic view of what the data looks like, since it does not require much background knowledge, usually incorporated into the class labels, nor involves pre-or post-treatment of the data. Applications of clustering techniques range over a wide variety of domains, from biology to physics, from e-commerce to social network analysis. As such, there exists a wide literature which investigates new algorithms, measures, and application domains for data clustering.

Usually, data come in the form of a set of objects described by a set of attributes that can be heterogeneous for type and underlying properties. For instance, the most common ones are numeric, categorical, boolean, textual, sequential, or networked. In the most classic settings, objects are represented by features obtained by a unique source of information or view. For instance, census data are produced by a national statistical institute, textual data are extracted from a document repository, tagged items come from a social media website, and so on.

However, in some cases, computing clusters over a single view is not sufficient to capture all the intrinsic similarities or differences between objects. Consider, for instance, an e-commerce website specialized in high-tech products, where each product is described by a general description, a technical note, and reviews are posted by customers. It is quite common to look for certain products and be unable to decide because of some incomplete technical details, uninformative descriptions or misleading reviews. The typical customer behavior is then to merge all the different sources of information, by considering them for what they are: technical notes give details on components and features (usually taken from the producer's website); general descriptions provide synthetic views of the main characteristics and advantages of the products; customers' reviews possibly complete missing details and add further useful information.

As another example consider a scientific digital library where each paper is described by its textual content (e.g., title, abstract and/or keywords), by the list of co-authors, and by a list of cited papers. If the goal is, for instance, to find the scientific communities interested in the same research field, one can try to apply a clustering algorithm over one of the feature spaces. However, intersections between communities are very common: for instance, people interested in knowledge discovery and data mining may work and publish together with people coming from the information retrieval community. Abstracts or titles on their own are not sufficiently informative: they might focus on the applicative aspects of the work, or be too generic or specific. Citations can be misleading as well: in some cases papers refer to other works dealing with similar topics but in other fields of research; sometimes papers are cited for another purpose than reporting closely related works. Therefore, clustering papers using only one feature space could lead to imprecise results.

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z Instead, taking into account all the feature space may help in identifying more interesting and relevant clusters.

A straightforward solution consists in merging all the feature spaces into a single one, and applying a clustering algorithm on the overall table. However this solution would not always work as expected, and this for several reasons: first the feature space may contain different types of attributes; second, even though the attribute type is homogeneous, each feature may fit a specific model, show properties or be affected by certain problems such as noise, which might affect in a different form the other feature space. Third, the form of data may vary significantly from one space to another: document-keyword data are typically sparse, while categorical data are usually dense. Networked data, such as social networks, contains highly dense portions of nodes, and sparse ones. Fourth, even though the spaces do not fall into the above-mentioned cases, their cardinality may vary drastically having the consequence that larger feature spaces have more influence on the clustering process than smaller ones.

To cope with the limitation of traditional clustering algorithms on multiple views data, a new clustering formulation was introduced in 2004, the so-called multi-view clustering approach [START_REF] Bickel | Multi-view clustering[END_REF]. Instead of considering a standard clustering approach over each single-view, the authors proposed an algorithm based on EM that clusters each view on the basis of the partitions of the other views. A strongly related approach is the so-called star-structured high-order heterogeneous co-clustering [START_REF] Gao | Star-structured high-order heterogeneous data coclustering based on consistent information theory[END_REF]. Though it has never been explicitly linked to multi-view clustering, this family of algorithms addresses the same issues of a co-clustering approach that simultaneously clusters objects and features from the different spaces. Co-clustering methods have two main advantages: they are more effective in dealing with the well-known problem of the curse of dimensionality, and they provide an intensional description of the object clusters by associating a cluster of features to each cluster of objects. When applied on multi-view data, the several views, or feature spaces, as well as the object set, are simultaneously partitioned. A common problem of existing co-clustering approaches (as well as of other dimensionality reduction techniques such as non-negative matrix factorization [START_REF] Lee | Algorithms for non-negative matrix factorization[END_REF]), is that the number of clusters/components is a parameter that has to be specified before the execution of the algorithm. Choosing a correct number of clusters is not easy, and mistakes can have undesirable consequences. As mentioned in [START_REF] Keogh | Towards parameter-free data mining[END_REF], incorrect parameter settings may cause an algorithm to fail in finding the true patterns and may lead the algorithm to report spurious patterns that do not really exist. Existing multi-view co-clustering methods have as many parameters as there are views in the data, plus one for the partition of objects.

In this paper we propose a new co-clustering formulation for high-order star-structured data that automatically determines the number of clusters at running time. This is obtained thanks to the employment of a multiobjective local search approach over a set of Goodman-Kruskal's τ objective functions, that will be introduced in Section 3.1. The local search approach aims at maximizing these functions, comparing two solutions on the basis of a Paretodominance relation over the objective functions. Due to the fact that those

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z functions have a defined upper limit that is independent of the numbers of clusters, the τ functions can be used to compare co-clusterings of different sizes. On the contrary, measures that are usually employed in co-clustering approaches, like the loss of mutual information [START_REF] Dhillon | Information-theoretic co-clustering[END_REF] or more generally the Bregman divergence [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], are optimized for a specific partition, usually the discrete one. Moreover, the stochastic local search algorithm we proposed enables the adaptation of the number of clusters during the optimization process. Last but not least, unlike most of the existing approaches, our algorithm is designed to fit any number of views: when a single feature space is considered, the algorithm is equivalent to the co-clustering approach presented in [START_REF] Robardet | Efficient local search in conceptual clustering[END_REF].

The remainder of this paper is organized as follows: Section 2 briefly explores the state of the art in star-structured co-clustering and multi-view clustering. The theoretic fundamental details are presented in Section 3 while the technical issues of our approach are provided in Section 4. In Section 5 we present the results of a comprehensive set of experiments on multi-view high-dimensional text and image data, as well as a qualitative evaluation on a specific result and a scalability analysis. Finally, Section 6 concludes.

Related Work

Co-clustering has been studied in many different application contexts including text mining [START_REF] Dhillon | Information-theoretic co-clustering[END_REF], gene expression analysis [START_REF] Cho | Minimum sum-squared residue co-clustering of gene expression data[END_REF][START_REF] Pensa | Constrained co-clustering of gene expression data[END_REF] and graph mining [START_REF] Chakrabarti | Fully automatic crossassociations[END_REF] where these methods have yielded an impressive improvement in performance over traditional clustering techniques. The methods differ primarily by the criterion they optimize, such as minimum loss in mutual information [START_REF] Dhillon | Information-theoretic co-clustering[END_REF], sum-squared distance [START_REF] Cho | Minimum sum-squared residue co-clustering of gene expression data[END_REF], minimum description length (MDL) [START_REF] Chakrabarti | Fully automatic crossassociations[END_REF], Bregman divergence [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF] and non-parametric association measures [START_REF] Robardet | Comparison of three objective functions for conceptual clustering[END_REF][START_REF] Ienco | Parameter-free hierarchical co-clustering by n-ary splits[END_REF]. Among these approaches, only those ones based on MDL and association measure are claimed to be parameter-free [START_REF] Keogh | Towards parameter-free data mining[END_REF]. However, methods based on MDL are strongly restricted by the fact they can only handle binary matrices. Association measures, such as Goodman and Kruskal τ , are internal measures of the quality of a co-clustering based on statistical considerations. They have also another advantage: they can deal with both binary and counting/frequency data [START_REF] Ienco | Parameter-free hierarchical co-clustering by n-ary splits[END_REF][START_REF] Robardet | Comparison of three objective functions for conceptual clustering[END_REF]. From an algorithmic point of view, the co-clustering problem has been shown to be NP-hard [START_REF] Anagnostopoulos | Approximation algorithms for coclustering[END_REF] when the number of row and column clusters are fixed. Therefore, proposed methods so far are based on heuristic approaches.

Several clustering and co-clustering methods for heterogeneous star-structured data have been proposed recently. Long et al. [START_REF] Long | Spectral clustering for multi-type relational data[END_REF] use spectral clustering to iteratively embed each type of data objects into low dimensional spaces in a way that takes advantage of the interactions among the different feature spaces. A partitional clustering algorithm like k-means is then employed to obtain the final clustering computed on the transformed spaces. In [START_REF] Long | A probabilistic framework for relational clustering[END_REF], a parametric probabilistic approach to cluster relational data is proposed. A Monte Carlo simulation method is used to learn the parameters and to assign objects to clusters. The problem of clustering images described by segments and captions is considered in [START_REF] Bekkerman | Multi-modal clustering for multimedia collections[END_REF]. The proposed algorithm is based on Markov random fields in which some of the nodes are random variables in the combinatorial problem.

The co-clustering problem on star-structured data was first considered in [START_REF] Gao | Star-structured high-order heterogeneous data coclustering based on consistent information theory[END_REF] where Gao et al. propose to adapt the Information Theory co-clustering approach [START_REF] Dhillon | Information-theoretic co-clustering[END_REF] to star-structured data. It consists in optimizing a weighted combination of mutual information evaluated over each feature space, where weights are chosen based on the supposed reliability/relevance of their correlation. Beyond the parameters inherited from the original algorithm, the weight involved in the linear combination also has to be fixed by the end-user. Another drawback of this approach is its complexity, that prevents its use on large-scale datasets. Greco et al. [START_REF] Greco | Co-clustering multiple heterogeneous domains: Linear combinations and agreements[END_REF] propose a similar approach based on the linear combination of mutual information evaluated on each feature space, where the parameter of the linear combination is automatically determined. Ramage et al. [START_REF] Ramage | Clustering the tagged web[END_REF], propose a generative clustering algorithm based on latent Dirichlet allocation to cluster documents using two different sources of information: document text and tags. Each source is modeled by a probability distribution and a weight value is used to weigh one vector space with respect to the other. During the learning step, the algorithm finds the distribution parameters, and models documents, words and tags. In addition to the weight parameter, the method has another drawback: it constrains the number of hidden topics in text and tag sources to be the same, which is a strong assumption on data that is not always true. Considering the multi-view clustering problem Cleuziou et al. [START_REF] Cleuziou | Cofkm: A centralized method for multiple-view clustering[END_REF] propose to find a consensus between the clusters from different views. Their approach merges information from each view by performing a fusion process that identifies the agreement between the views and solves the conflicts. To this purpose they formulate the problem using a fuzzy clustering approach and extend fuzzy k-means algorithm by introducing a penalty term. This term aims at reducing the disagreement between any pair of partitions coming from the different views. This approach is sensitive to the parameter setting because, apart from the number of clusters, the algorithm requires two other parameters that may influence the result. [START_REF] Chen | Semi-supervised document clustering with simultaneous text representation and categorization[END_REF] can be considered as the first attempt that performs multi-view co-clustering. The method is an extension of the Non-Negative Matrix Factorization approach to deal with multi-view data. It computes new word-document and document-category matrices by incorporating user provided constraints through simultaneous distance metric learning and modality selection. This method is shown to be effective, but its formulation is not flexible. In fact, in order to use more than two feature spaces, one needs to reformulate the whole process. Another problem of this approach is that the number of clusters for each feature space is given as a parameter. Furthermore, the number of parameters grows with the number of feature spaces.

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z 3 Measuring the quality of a co-clustering on star-structured data

In this section, we introduce the association measure optimized within our co-clustering approach. Originally designed for evaluating the dependence between two categorical variables, the Goodman and Kruskal τ measure can be used to evaluate the quality of a single-view co-clustering. We generalize its definition to a multi-view setting.

Evaluating the quality of a co-clustering

Association measures are statistical measures that evaluate the strength of the link between two or more categorical variables. Considering partitions of a co-clustering as categorical variables, these measures have been shown to be well adapted to determine co-clustering with high quality [START_REF] Robardet | Comparison of three objective functions for conceptual clustering[END_REF]. Goodman and Kruskal τ measure [START_REF] Goodman | Measures of association for cross classification[END_REF] is one of them that estimates the association between two categorical variables X and Y by the proportional reduction of the error in predicting X knowing or not the variable Y :

τ X|Y = e X -E[e X|Y ] e X Let say that variable X has m categories X 1 , • • • , X m with probabilities p 1 , • • • , p m , variable Y has n categories Y 1 , • • • , Y n with probabilities q 1 , • • • , q n and their joint probabilities are denoted r ij , for i = 1 • • • m and j = 1 • • • n. The error in predicting X can be evaluated by e X = m i=1 p i (1 -p i ) = 1 - m i=1 p 2
i . This is the probability that two independent observations from the marginal distribution of X fall in different categories. The error is minimum (equals to 0) when every p i 's is 0 or 1. The error is maximum when p i = 1 m for all i. E[e X|Y ] is the expectation of the conditional error taken with respect to the distribution of Y :

E[e X|Y ] = n j q j e X|Yj = n j q j m i r ij q j (1 - r ij q j ) = 1 - m i n j r 2 ij q j
The Goodman-Kruskal's τ X|Y association measure is then equal to:

τ X|Y = i j r 2 ij qj -i p 2 i 1 -i p 2 i (1)
Similarly, the proportional reduction of the error in predicting Y while X is known is given by:

τ Y |X = e Y -E[e Y |X ] e Y = i j r 2 ij pi -j q 2 j 1 -j q 2 j (2)
Let's now describe how this association measure can be employed to quantify the goodness of a co-clustering. Let O be a set of objects and F a set of Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z features. A dataset D associates with an object o u and a feature f v a value d uv that represents the frequency of feature v in the description of object u (or it could also be a Boolean value that stands for presence/absence of the feature). Figure 1(a) is an example of such a dataset. A co-clustering is made of one partition of objects, CO 1 , • • • , CO m , and one partition of features, CF 1 , • • • , CF n . To evaluate the quality of the co-clustering, we compute a contingency table that empirically estimates the joint probabilities of the two partitions and the marginal probabilities of each partition. It associates with each co-cluster (CO i , CF j ) the value t ij that aggregates the frequency values of features from cluster CF j on objects of CO i : In order to review the meaning of τ for the evaluation of a co-clustering consider the table in Figure 1(b) whose cell at the intersection of the row CO i with the column CF j contains the frequency of objects in cluster CO i having the features in cluster CF j . T Oi is the total counting for cluster CO i , T Fj is the sum of counts for cluster CF j , and T is the global total, i.e., T Oi = j t ij , T Fj = i t ij and T = i j t ij . Therefore, the probability that an object o u is in cluster CO i may be estimated by p i = TO i T , the relative frequency of the cluster. Similarly, the probability that a feature f v is in cluster CF j is estimated by q j = TF j T . Finally, the joint probability that an object o u and a feature f v are in co-cluster (CO i , CF j ) is r ij = tij T . Such a statistical estimate is represented by a two-dimensional contingency table or a two-way frequency table. We denote by X the random variable associated with the partition of objects, and by Y the random variable associated with the partition of features. τ X|Y is the proportional reduction in prediction error of partition

t ij = ou∈COi fv ∈CFj d uv (3)
T O 2 T F 1 T F 2 T (b)
{CO 1 , • • • , CO m } given partition {CF 1 , • • • , CF n }.
Example 1 Consider the dataset in Figure 2(a). It could represent one of the possible views on a document collection: the view on the document terms or

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z the view on the document tags that users adopted to annotate documents in a social network. Suppose the view is on the terms contained in the documents: each value in the matrix represents the number of occurrences of a specific word (indicated by the matrix column) in a specific document (indicated by the matrix row). We take into account two possible co-clusterings for this dataset:

C 1 = {{o 1 , o 2 , o 3 }, {o 4 , o 5 }}, {{f 1 , f 2 }, {f 3 , f 4 }} and C 2 = {{o 1 , o 2 }, {o 3 , o 4 , o 5 }}, {{f 1 , f 2 }, {f 3 , f 4 }}
Tables in Figure 2(b) and 2(c) represent the contingency tables for the first and second co-clustering respectively. The Goodman-Kruskal τ X|Y and τ Y |X for the first co-clustering are:

τ X|Y = τ Y |X = 0.5937
and for the second one are:

τ X|Y = 0.2158 and τ Y |X = 0.3375
We can observe that the first co-clustering better associates clusters of objects with clusters of features than the second one. In the first co-clustering, objects of CO 1 have mostly features of cluster CF 1 , whereas objects of CO 2 are characterized by features of CF 2 . In the second example, given an object having high frequency values on features of CF 1 , it is difficult to predict if it belongs to CO 1 or CO 2 . The τ measure corroborates these facts and in fact it has higher values for the first co-clustering. Analyzing the properties of τ , we can observe that it satisfies many desirable properties of a co-clustering measure. First, it is invariant by rows and

f 1 f 2 f 3 f 4 o 1 3 4
Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z columns permutation. Second, it takes values between [0, 1]: (1) τ X|Y is 0 iff knowledge of Y classification is of no help in predicting the X classification, i.e. there is independence between the object and feature clusters; (2) τ X|Y is 1 iff knowledge of an individual's Y cluster completely specifies his X class, i.e. each row of the contingency table contains at most one non zero cell. Third, it has an operational meaning: given an object, it is the relative reduction in the prediction error of the object's cluster given the partition of features, in a way that preserves the category distribution [START_REF] Goodman | Measures of association for cross classification[END_REF]. Finally, unlike many other association measures such as χ 2 , τ has a defined upper limit that is independent of the numbers of classes m and n. Therefore, τ can be used to compare co-clusterings of different sizes.

The star schema

In the context of high-order star-structured co-clustering, the same set of objects is represented in different feature spaces. Such data represent objects of a certain type, connected to other types of data, the features, so that the overall data schema forms a star structure of inter-relationships. The co-clustering task consists in clustering simultaneously the set of objects and the set of values in the different feature spaces. In this way we obtain a partition of the objects influenced by each of the feature spaces and at the same time a partition of each feature space. Considering N feature spaces that have to be clustered, in this paper we study a way to extend the τ measure to perform co-clustering of the data when the number of the involved dimensions is high.

In Figure 3 we show an example of a star-structured data schema. The set of objects is shown at the center of the star. For example it could be a set of documents. The set of objects is then described with multiple views by more features (Feature1, Feature2, Feature3). For instance, Feature1 could be the set of the vocabulary words and the values d 1 ij represent the number of occurrences of the word f 1 j in the document Oi; Feature2 could be the set of tags adopted by the users of a social community in their annotations and the values d 2 ij represent the number of users that used the tag f 2 j in the document Oi; Feature3 could be the set of authors who cited documents in the given set and the feature values represent the number of citations. For the complete description of each object we have to make reference to the different feature spaces. The values of each feature have their own distributions in the document set and are described by the different contingency tables. At this point both the document set and the feature spaces might be partitioned into co-clusters where all the dimensions contribute to determine the partitions.

Evaluating the quality of a star-structured co-clustering

Evaluating the quality of the partition of objects, given the partitions of features, is formalized as follows. The partition of objects is associated with the Feature 2

Fig. 3 The data star-schema.

dependent variable X, and the N partitions of the feature spaces are considered as many independent variables

Y = {Y 1 , . . . , Y N }. Each variable Y k ∈ Y has n k categories Y k 1 , • • • , Y k n k with probabilities q k 1 , . . . , q k n k and X has m cate- gories X 1 , • • • , X m . However, for each variable Y k , the m categories of X have different probabilities p k 1 , • • • , p k m , k = 1 • • • N . The joint probabilities between X and any Y k ∈ Y are denoted by r k ij , for i = 1 • • • m and j = 1 • • • n k .
The error in predicting X is the sum of the errors over the independent variables of Y:

e X = N k=1 m i=1 p k i (1 -p k i ) = N - N k=1 m i=1 (p k i ) 2 . E[e X|Y ]
is the expectation of the conditional error taken with respect to the distributions of all Y k ∈ Y:

E[e X|Y ] = N k n k j q k j e X|Y k j = N k n k j q k j m i r k ij q k j (1- r k ij q k j ) = N - N k m i n k j (r k ij ) 2 q k j
The generalized Goodman-Kruskal's τ X|Y association measure is then equal to:

τ X|Y = k i j (r k ij ) 2 q k j -k i (p k i ) 2 N -k i (p k i ) 2 (4)
If we consider Y k as a dependent variable, and X as an independent variable, the corresponding τ Y k |X is computed as described in the previous section (see

τ Y k |X = e Y k -E[e Y k |X ] e Y k = i j (r k ij ) 2 p k i -j (q k j ) 2 1 -j (q k j ) 2 (5)
We now describe how the generalized association measure can be employed for computing the quality of a high-order star-structured co-clustering. Let F = {F 1 , . . . , F N } be a set of N feature spaces. A dataset can be viewed under the different views given by the different feature spaces F k . Therefore, the view D k is associated with each feature space F k . We call star-structured dataset, denoted D = {D 1 , . . . , D N }, the collection of matrices built over O and F . The joint tables in Figure 4(b) are the contingency tables representing a co-clustering of the multi-space dataset of Figure 4(a). Symbol CO i represents i-th cluster on objects, while CF k j represents j-th cluster on the kth feature space (e.g.

(CO 1 = {o 1 , o 2 , o 5 }, CO 2 = {o 3 , o 4 }), CF 1 1 = {f 1 1 , f 1 2 }, CF 1 2 = {f 1 3 , f 1 4 } and CF 2 1 = {f 2 1 }, CF 2 2 = {f 2 2 , f 2 3 } . The association of co- cluster (CO i , CF k j
) is represented by the value t k ij , at the intersection of the i-th object cluster and the j-th cluster of the k-th feature space, and is equal to:

t k ij = ou∈COi f k v ∈CF k j d k uv (6)
T k Oi is the total counting for cluster CO i within the k-th feature space. T F k j is the total counting for the j-th cluster of the k-th feature space, and T k is the global total for the k-th feature space, i.e.,

T k Oi = j t k ij , T F k j = i t k ij and T k = i j t k ij .
Therefore, given a feature space k, the probability that an object o u is in cluster

CO i is p k i = T k O i T k . Similarly, the probability that a feature f k v is in cluster CF k j is q k j = T F k j T k .
Finally, the joint probability that an object o u and a feature f k v are in co-cluster

(CO i , CF k j ) is r k ij = t k ij T k .
Let X be the random variable associated with the partition on objects: its values range over all the possible partitions of O. Let Y be the set of N random variables Y k , k = 1 . . . N , associated with the feature partitions. Each random variable Y k has values that range over all the possible partitions of F k . Figure 4(b) shows that a possible value for X is the partition {CO 1 , CO 2 } and the values of the variables Y 1 and Y 2 are respectively the partitions {CF 1 1 , CF 1 2 } and {CF 2 1 , CF 2 2 }. τ X|Y is the proportional reduction in the prediction error of the object partition when the values of the feature partitions are known. The best coclustering solution will be the one that reduces the most the prediction error of the object partition when the object values on the feature partitions are known. As we will see, also the opposite prediction problem will be considered in our approach, by evaluation of τ Y|X : it considers the reduction in the prediction error of the feature partitions when the object partition is known.

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z (a)

f 1 1 f 1 2 f 1 3 f 1 4 f 2 1 f 2 2 f 2 3 o 1 d 1
CF 1 1 CF 1 2 CF 2 1 CF 2 2 CO 1 t 1 11 t 1 12 T 1 O 1 t 2 11 t 2 12 T 2 O 1 CO 2 t 1 21 t 1 22 T 1 O 2 t 2 21 t 2 22 T 2 O 2 T F 1 1 T F 1 2 T 1 T F 2 1 T F 2 2
T 2 (b) Fig. 4 An example of a star-structured dataset (a) and the contingency table associated with a related star-structured co-clustering (b).

Example 2 Example in Figure 4(a) has been instantiated in the dataset of Figure 5(a): it could represent a collection of documents. Suppose the first feature space is the space of the terms occurring in the documents, while the second feature space is the space of the tags associated by the users of a social network. Each counting in the first feature space (f

1 1 • • • f 1 4
) represents, for instance, the occurrences of a specific word in a specific document, while the counting in the second feature space (f

2 1 • • • f 2 3
) are the number of users who associated a specific tag with a specific document. We take into account two possible co-clusterings for this dataset: 5(b) and 5(c) represent the contingency tables for C 1 and C 2 respectively. For the first co-clustering, the Goodman-Kruskal τ X|Y and τ Y k |X measures are:

C 1 = {CO 1 = {o 1 , o 2 , o 3 }, CO 2 = {o 4 , o 5 }} × {{CF 1 1 = {f 1 1 , f 1 2 }, CF 1 2 = {f 1 3 , f 1 4 }} ∪ {CF 2 1 = {f 2 1 }, CF 2 2 = {f 2 2 , f 2 3 }}} C 2 = {CO 3 = {o 1 , o 3 , o 5 }, CO 4 = {o 2 , o 4 }} × {{CF 1 3 = {f 1 1 , f 1 3 }, CF 1 4 = {f 1 2 , f 1 4 }} ∪ {CF 2 1 = {f 2 1 }, CF 2 2 = {f 2 2 , f 2 3 }}} Tables in Figure
τ X|Y = 0.6390, τ Y 1 |X = 0.5937, τ Y 2 |X = 0.6890
τ values have been computed respectively on the contingency tables:

{CO 1 , CO 2 } × {CF 1 1 , CF 1 2 , CF 2 1 , CF 2 2 } {CF 1 1 , CF 1 2 } × {CO 1 , CO 2 } {CF 2 1 , CF 2 2 } × {CO 1 , CO 2 }
For the second co-clustering, they are:

τ X|Y = 0.0119, τ Y 1 |X = 0.0234, τ Y 2 |X = 0.0008
Since all the three τ values for the first co-clustering are greater than the second ones, we can conclude that the first co-clustering better captures the interactions between objects and features of each space.

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z

f 1 1 f 1 2 f 1 3 f 1 4 f 2 1 f 2 2 f 2 3 o 1 3 4 1 1 0 8 5 o 2 5 3 0 2 0 6 9 o 3 6 4 1 0 2 2 2 o 4 0 1 7 7 9 1 0 o 5 1 0 6 8 7 0 1 (a) CF 1 1 CF 1 2 CF 2 1 CF 2 2 CO 1 25

A stochastic local search approach for high-order star-structured co-clustering

We formulate our co-clustering approach for star-structured data as a multi-objective combinatorial optimization problem [START_REF] Jaszkiewicz | Genetic local search for multi-objective combinatorial optimization[END_REF] which aims at optimizing N + 1 objective functions based on Goodman-Kruskal's τ measure. Given a star-structured dataset D over O and F , the goal of the star-structured data co-clustering is to find a set of partitions Y = {Y 1 , . . . , Y k , . . . , Y N } over the feature set F = {F 1 , . . . , F k , . . . , F N }, and a partition X of the object set O such that

max Y∪X∈P τ (Y ∪ X) = τ Y 1 |X , . . . , τ Y k |X , . . . , τ Y N |X , τ X|Y ( 7
)
where P is the discrete set of candidate partitions and τ is a function from P to Z, where Z = [0, 1] N +1 is the set of vectors in the objective space, i.e., z i = τ Y i |X (see equation ( 5)), ∀i = 1, . . . , N and z N +1 = τ X|Y (see equation ( 4)).

To compare different candidate partitions of P, one can either use a scalarization function, that maps Z into R (such as the weighted sum), or one can employ a dominance-based approach, that induces a partial order over P. As for combinatorial problem there may exist many truly optimal partitions that are not optimal for any weighted sum scalarization function [START_REF] Liefooghe | On dominancebased multiobjective local search: design, implementation and experimental analysis on scheduling and traveling salesman problems[END_REF], we consider in the following a Pareto-dominance approach. The goal is to identify optimal partitions P opt ⊂ P, where optimality consists in the fact that no solution in P \ P opt is superior to any solution of P opt on every

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z objective function. This set of solutions is known as Pareto optimal set or non-dominated set. These concepts are formally defined below.

Definition 1 (Pareto-dominance) An objective vector z ∈ Z is said to dominate an objective vector z ′ ∈ Z iff ∀i ∈ {1, . . . , N + 1}, z i ≥ z ′ i and ∃j ∈ {1, . . . , N + 1}, such that z j > z ′ j . This relation is denoted z ≻ z ′ hereafter.

Definition 2 (Non-dominated objective vector and Pareto optimal solution) An objective vector z ∈ Z is said to be non-dominated iff there does not exist another objective vector z ′ ∈ Z such that z ′ ≻ z.

A solution p ∈ P is said to be Pareto optimal iff its mapping in the objective space (τ (p)) results in a non-dominated vector.

Thus, the star-structured data co-clustering problem is to seek for a Pareto optimal solution of equation [START_REF] Chen | Non-negative matrix factorization for semisupervised heterogeneous data coclustering[END_REF]. As this problem is more difficult than the co-clustering problem with fixed numbers of object and feature clusters, which is known to be NP-hard [START_REF] Anagnostopoulos | Approximation algorithms for coclustering[END_REF], we propose, in subsection 4.1, a Stochastic Local Search algorithm CoStar to solve it. This algorithm searches for a near optimal solution given a local knowledge provided by the definition of a neighborhood [START_REF] Paquete | Stochastic Local Search Algorithms for Multiobjective Combinatorial Optimization: Methods and Analysis[END_REF]. We demonstrate, in subsection 4.2, that CoStar outputs a Pareto local optimum solution of the problem. In subsection 4.3, we propose a way to evaluate a neighbor solution incrementally. Finally subsection 4.4 gives the overall complexity of CoStar and of its competitors.

Multiple feature space co-clustering algorithm

We build our algorithm on the basis of τ CoClust [START_REF] Robardet | Efficient local search in conceptual clustering[END_REF], whose goal is to find a partition X of objects and a partition Y of features such that Goodman-Kruskal's τ X|Y and τ Y |X are maximized. The algorithm locally optimizes the two coefficients by updating step by step the partitions. The basic updating procedure consists in running through the neighborhood of each partition and choose the partition that increases the most the measures. This neighborhood is made of partitions where a peculiar element is moved from one cluster to another. Let us first consider the improvement of partition X. After randomly picking up a cluster and randomly picking up an object in that cluster, the procedure assigns this object to the cluster (possibly an empty one) that mostly increases τ X|Y . In a following phase, it considers the improvement of the features partition Y : it randomly picks up a feature from a randomly obtained cluster and assigns this feature to the cluster that mostly increases τ Y |X . Starting from the discrete partitions, the algorithm alternatively updates the partitions X and Y until a convergence criterion is satisfied. Unlike other co-clustering approaches, that require the number of co-clusters as parameter, τ CoClust is able to automatically determine the most appropriate numbers of clusters for X and Y . This is due

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z to two characteristics: (1) the τ coefficient has a defined upper limit that is independent of the numbers of clusters (while other measures, such as the loss in mutual information [START_REF] Dhillon | Information-theoretic co-clustering[END_REF] have not) and thus enables comparison of co-clusterings of different sizes, and (2) the updating procedure can create or delete clusters at any time during the iterative process.

Let us now present algorithm CoStar as an extension of τ CoClust algorithm to multi-view co-clustering. CoStar is shown in Algorithm 1. It takes in input the dataset D and alternatively optimizes the partition of the objects, working on the random variable X, and each partition over each feature set F k , working on the random variables Y k , until a number of iterations N iter is reached.

Algorithm 1 CoStar(D, N iter ) 1: Initialize Y 1 , • • • , Y N , X with discrete partitions 2: i ← 0 3: T ← ∅ 4: for k = 1 to N do 5: T k ← ContingencyTable(X, Y k , D k ) 6:
T ← T T k 7: end for 8: while (i ≤ N iter ) do 9:

[X, T ] ← OptimizeMultiObjectCluster(X, Y, T ) 10:

for k = 1 to N do 11: [Y k , T k ] ← OptimizeFeatureCluster(X, Y k , T k ) 12:
end for 13:

i ← i + 1 14: end while 15: return Y 1 , • • • , Y N , X
CoStar first initializes the partitions Y 1 , • • • , Y N and X with discrete partitions (that is, partitions with a single element per cluster) and then computes the contingency tables T k between the partition X and alternatively each partition Y k over the feature set F k . The function Con-tingencyTable computes equation (6) (initialization step at lines 4-7). Then, from line 8 to 14, it performs a stochastic local search method by alternatively optimizing partition X (line 9) and each partition Y k (line 11). The candidate partitions considered by the local search algorithm belong to the neighborhood of Y ∪ X as defined in equation ( 8) page 19. Function OptimizeMultiObjectCluster substitutes the object partition X for the partition in the neighborhood of X that mostly increases τ X|Y (see Algorithm 2). This measure takes into account all the partitions Y k ∈ Y and their associated contingency tables T = k T k . Function Optimize-FeatureCluster substitutes the feature partition Y k for the partition in the neighborhood of Y k that mostly increases τ Y k |X (see Algorithm 4). These operations are alternated and repeated until a stopping condition is satisfied. This condition can be based on a convergence criterion of τ , but, if

(τ X ′ |Y (T ′ ) > maxτ X|Y ) then 8: maxτ X|Y ← τ X ′ |Y (T ′ ) 9: N D ← {xe} 10: else 11: if (τ X ′ |Y (T ′ ) = maxτ X|Y ) then 12: N D ← N D ∪ {xe} 13:
end if 14:

end if 15: end for 16: max X , max T ← FindParetoLocalOptimumObject(N D, X, T, x b , o) 17: return max X , max T Algorithm 3 FindParetoLocalOptimumObject(N D, X, T, x b , o) 1: while Card(N D) > 1 do 2: Take x i , x j ∈ N D s.t. x i = x j 3: nb + x i ← 0 4: nb + x j ← 0 5: [X i , T i ] ← Update(X, T, o, x b , x i ) 6:
[X j , T j ] ← Update(X, T, o, x b , x j ) 7:

for k = 1 to N do 8:

if

(τ Y k |X i (T i ) > τ Y k |X j (T j )) then 9: nb + x i ← nb + x i + 1 10: else 11: if (τ Y k |X i (T i ) < τ Y k |X j (T j )) then 12: nb + x j ← nb + x j + 1 13:
end if 14:

end if 15:

end for 16:

if nb + x i ≥ nb + x j then 17: N D ← N D \ {x j } 18: else 19: N D ← N D \ {x i } 20: end if 21: end while 22: return [X i , T i ] ← Update(X, T, o, x b , x i ), with x i ∈ N D
Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z

The stochastic local optimization of the partition X works as follows (see Algorithm 2). The algorithm starts by randomly picking up a cluster x b of X and an element o of x b (lines 1 and 2). Then, all the partitions of the neighborhood of X, constituted by the partitions where o is moved to another cluster of X (possibly the empty one), are considered. To retrieve the partition max X that increases at most τ X|Y among the partitions in the neighborhood of X and that is not dominated by any partition of this neighborhood, the algorithm uses a queue, N D, that gathers the possibly non-dominated partitions of this neighborhood. Each partition of the neighborhood of X is identified by the cluster to which o belongs. Thus, the queue N D is initialized with x b that represents the partition X (line 3). max τ X|Y is initialized with the current τ X|Y value. The loop from line 5 to 15 considers all the neighboring partitions of X (x e being the cluster where o is moved to, x e ∈ X ∪ ∅). The function Update (line 6) is called to modify the contingency tables T = k T k and the partition X with respect to the new cluster of o. To perform this step, it first moves the element o from the cluster x b to cluster x e and then performs one or more of the following three actions: 1. if cluster x b is empty after removing o, it deletes cluster x b and updates the contingency tables T ′ consequently; 2. if cluster x e is the empty cluster, it adds a new cluster to the partition and updates the contingency tables T ′ consequently; 3. if the two above mentioned cases do not apply, it simply updates the content of T ′ by modifying x b and x e rows.

Each modification of T ′ consists actually in modifying the N contingency tables T k , k = 1, . . . , N . From line 7 to 14, potentially non-dominated partitions of the neighborhood of X are stored in N D. τ X ′ |Y is evaluated using the contingency tables T ′ . If it is strictly greater than max τ X|Y then, given a partition of N D, either X ′ dominates it or neither of the partitions dominates the other (see proof of theorem 1 page 19). Therefore, N D and max τ X|Y are reinitialized respectively with X ′ and τ X ′ |Y . If τ X ′ |Y is equal to max τ X|Y then x e is added to the queue N D. Considering X ′ and any partition X of N D, all the Pareto dominance relation may exist between these two partitions: X ′ may dominate X, X ′ may be dominated by X or neither may dominate the other. Therefore, a partition that is not dominated by any partition of N D will be identified by the function FindParetoLocalOptimumObject.

When all the neighborhood has been explored, the function FindPare-toLocalOptimumObject (see Algorithm 3) is called. It processes the queue N D in order to retrieve a partition non dominated by any other partition of N D. To that end, the function compares two partitions X i and X j of N D by counting in nb the function removes X j from N D and continues until there is only one partition in N D.

Algorithm 4 OptimizeFeatureCluster(X, Y k , T k ) 1: Randomly choose a cluster y b ∈ Y k 2: Randomly choose a feature f ∈ y b 3: N D ← {y b } 4: maxτ Y k |X ← τ Y k |X (T k ) 5: for all ye ∈ {Y k ∪ ∅} s.t. y b = ye do 6: [Y k′ , T k′ ] ← Update(Y k , T k , f, y b , ye) 7: if (τ Y k′ |X (T k′ ) > maxτ Y k |X ) then 8: maxτ Y k |X ← τ Y k′ |X (T k′ ) 9: N D ← {ye} 10: else 11: if (τ Y k′ |X (T k′ ) = maxτ Y k |X ) then 12: N D ← N D ∪ {ye} 13: end if 14: end if 15: end for 16: max Y k , max T k ← FindParetoLocalOptimumFeature(N D, Y k , T k , y b , f ) 17: return max Y k , max

Local convergence of CoStar

Considering iterated local search algorithms for multi-objective optimization requires to extend the usual definitions of local optimum to that context [START_REF] Paquete | A study of stochastic local search algorithms for the biobjective qap with correlated flow matrices[END_REF]. Let N : P → 2 P be a neighborhood structure that associates to every solution p ∈ P a set of neighboring solutions N (p) ⊆ P. We define a Pareto local optimum with respect to N as follows.

Definition 3 (Pareto local optimum) A solution p ∈ P is a Pareto local optimum with respect to the set of solutions N (p) that constitutes the neighborhood of p, if and only if there is no r ∈ N (p) such that τ (r) ≻ τ (p) (see equation ( 7)).

In algorithm CoStar, the neighborhood function is defined as follows. Let p = Y ∪ X be an element of P. W i denotes either the partition Y i , if i ∈ {1, . . . , N }, or the partition X, if i = N + 1. Let w b be a cluster of W i and let u be an element of w b . The neighborhood N of p is determined by Take y i , y j ∈ N D s.t. y i = y j 3:

nb + y i ← 0 4: nb + y j ← 0 5: [Y k i , T k i ] ← Update(Y k , T k , f, y b , y i ) 6: [X k j , T k j ] ← Update(Y k , T k , f, y b , y j ) 7: Y i ← Y \ Y k ∪ Y k i 8: Y j ← Y \ Y k ∪ Y k j 9: if (τ X|Y i (T k i ) > τ X|Y j (T k j )) then 10: nb + y i ← nb + y i + 1 11: else 12: if (τ X|Y i (T k i ) < τ X|Y j (T k j ))
[Y k i , T k i ] ← Update(Y k , T k , f, y b , y i ), with y i ∈ N D
the movement of the element u from the cluster w b to the cluster w e in W i . It is given by the substitution in p of W i by W i′ :

N (i,w b ,u) (Y ∪ X) = {Y ∪ X \ W i ∪ W i′ } (8)
with Proof We consider a partition W i ∈ {Y ∪ X}. Algorithm CoStar (in functions OptimizeMultiObjectCluster and OptimizeFeatureCluster) randomly picks-up a cluster w b from partition W i and an element u ∈ w b . Then, it systematically explores the neighborhood N (i,w b ,u) (Y ∪ X) and retrieves the set of the corresponding neighboring solutions (co-clustering) N D defined by:

W i′ = W i \ w b \ w e ∪
N D = {s ∈ N (i,w b ,u) (Y ∪ X) | ∀t ∈ N (i,w b ,u) (Y ∪ X), z i (s) ≥ z i (t)}
Considering a solution s ∈ N D and a solution r ∈ N (i,w b ,u) (Y ∪ X), we have z i (s) ≥ z i (r) and the following situations can occur: 2. If z i (s) > z i (r) and ∃k ∈ {1, • • • , N + 1}, k = i, such that z k (s) < z k (r), then, from Definition 1, s ≻ r and r ≻ s.

1. If z i (s) > z i (r) and ∀k ∈ {1, . . . , N + 1}, k = i, z k (s) ≥ z k (r),
In both cases 1 and 2, r does not belong to N D since z i (s) > z i (r) (see functions OptimizeMultiObjectCluster or Optimize-FeatureCluster from lines 7 to 9) and s is a Pareto local optimum with respect to r.

If z

i (s) = z i (r) and ∀k ∈ {1, • • • , N + 1}, k = i, z k (s) ≥ z k (r)
, then, either s is equivalent to r (z(s) = z(r)) or s dominates r: s r.

Therefore, s is a Pareto local optimum with respect to r. From Algorithms OptimizeMultiObjectCluster or OptimizeFeatureCluster, both r and s belong to N D. When partitions r and s are compared in Algorithms FindParetoLocalOptimumObject or Find-ParetoLocalOptimumFeature, nb + s ≥ nb + r and the algorithms delete r from N D.

If z

i (s) = z i (r) and ∃k ∈ {1, • • • , N + 1}, k = i, such that z k (s) < z k (r)
and ∃ℓ ∈ {1, • • • , N + 1}, ℓ = i = k, such that z ℓ (s) > z ℓ (r), then s ≻ r and r ≻ s.

Both partitions r and s are non-dominated from each other and are Pareto local optima. They both belong to N D, and Algorithms FindParetoLocalOptimumObject or FindParetoLocalOpti-mumFeature will remove the partition that is strictly greater than the other on the smallest number of objective functions z k .

If z

i (s) = z i (r) and ∀k ∈ {1, • • • , N + 1}, z k (s) ≤ z k (r), and ∃ℓ ∈ {1, • • • , N + 1}, ℓ = i, such that z ℓ (s) < z ℓ (r), then r ≻ s.
z N +1 (r) > z N +1 (s) Indeed, in that case, z k (r) = z k (s), ∀k = 1, • • • , N , since the modi- fication of Y i does not change τ Y k |X , ∀k = i.
Author We can conclude that at the end of each call of OptimizeMultiObject-Cluster and OptimizeFeatureCluster, the output co-clustering is a Pareto local optimum with respect to the considered neighborhood. Finally, all the iterations are bounded and thus CoStar terminates in a finite number of steps. When moving a single element from a cluster to another, a large part of the computation of these functions can be saved, since a large part in the τ formula remains the same. To take advantage of this point, we only compute the variations of τ Y k |X and τ X|Y from one step to another as explained in the following sections 4.3.1 and 4.3.2. We evaluate the variation of each measure when the partition X is modified as well as it is Y k that is changed.

Computing the variation of τ Y k |X

When partition Y k is modified: Let's first consider the variation of τ Y k |X when one element f is moved from cluster y b to cluster y e of Y k . This variation is the result of some changes on the co-occurrence table T k . Moving f induces, for each cluster of objects x i , the transfer of a quantity λ i from t k ib to t k ie of T k . The overall vector of transferred quantities is called λ = [λ 1 , . . . , λ m ]. In the following, we call t k ij the elements of the co-occurrence table T k before the update operation, and s k ij the same elements after the operation. We call T k and S k the two co-occurrence tables, and

τ Y k |X (T k ), τ Y k |X (S k ) the values of τ Y k |X computed over T k and S k .
The following equations hold for t k ij and s k ij : Therefore, the variation of τ Y k |X is:

s k ij = t k ij , if j = b, j = e s k ib = t k ib -λ i s k ie = t k ie + λ i Author-
∆ τ Y k |X (T k , f, y b , y e ) = τ Y k |X (T k ) -τ Y k |X (S k ) = i j (t k ij ) 2 t k i. t k .. -j (t k .j ) 2 (t k .. ) 2 1 -j (t k .j ) 2 (t k .. ) 2 - i j (s k ij ) 2 s k i. s k .. -j (s k .j ) 2 (s k .. ) 2 1 -j (s k .j ) 2 (s k .. ) 2
where

t k i. = j t k ij , t k .j = i t k ij and t k .. = i j t k ij . Setting λ = i λ i , Ω = 1 -j (t k .j ) 2 (t k .. ) 2 and Γ = 1 -i j (t k ij ) 2 t k i. t k ..
, we obtain the following updating formula:

∆ τ Y k |X (T k , f, y b , y e ) = Ω i 2λi t k i. t k .. [t k ib -t k ie -λ i ] + Γ 2λ (t k .. ) 2 [t k .e -t k .b + λ] Ω 2 -2 (t k .. ) 2 λΩ
s k ij = t k ij , if i = b, i = e s k bj = t k bj -µ k j s k ej = t k ej + µ k j
Therefore, the variation of τ Y k |X is equal to:

∆ τ Y k |X (T k , o, x b , x e ) = τ Y k |X (T k ) -τ Y k |X (S k ) = i j (t k ij ) 2 t k i. t k .. -j (t k .j ) 2 (t k .. ) 2 1 -j (t k .j ) 2 (t k .. ) 2 - i j (s k ij ) 2 s k i. s k .. -j (s k .j ) 2 (s k .. ) 2 1 -j (s k .j ) 2 (s k .. ) 2
Author where

t k i. = j t k ij , t k .j = i t k ij and t k .. = i j t k ij . Setting µ k . = j µ k j , Ω = 1 -j (t k .j ) 2 (t k .. ) 2
, we obtain the following updating formula:

∆ τ Y k |X (T k , o, x b , x e ) = 1 Ω j (t k ej ) 2 t k e. t k .. - (t k ej + µ k j ) 2 (t k e. + µ k . )t k .. + (t k bj ) 2 t k b. t k .. - (t k bj -µ k j ) 2 (t k b. -µ k . )t k ..
Thanks to this approach, instead of computing τ Y k |X at lines 8 and 11 of Algorithm 3 with a complexity in O(m • n k ), we compute ∆ τ Y k |X (T k , o, x b , x e ) once at line 8 with a complexity in O(n k ) (that is in O(|F k |) in the worst case with the discrete partition F k ) and test if it is strictly positive in the first if condition or strictly negative in the second one. Computing Ω is in O(n k ) and is done only once at the beginning of Algorithm 3.

Computing the variation of τ X|Y

When partition X is modified: Let's now consider the variation of τ X|Y when one object o is moved from the cluster x b to the cluster x e of partition X. This operation induces, for each cluster of features y k j , the transfer of a quantity µ k j from t k bj to t k ej of T k . The overall vector of transferred quantities is called µ = [µ 1 1 , . . . , µ 1 n1 , . . . , µ N 1 , . . . , µ N nN ]. We call t k ij the elements of the k th contingency table T k before the update operation, and s k ij the same elements after the operation. We call T = {T 1 , . . . , T N } and S = {S 1 , . . . , S N } the two sets of contingency tables, and τ X|Y (T ) and τ X|Y (S) the values of τ X|Y computed over T and S. The following equations hold for t k ij and s k ij :

s k ij = t k ij , if i = b, i = e s k bj = t k bj -µ k j s k ej = t k ej + µ k j
Therefore, the variation of τ X|Y is:

∆ τ X|Y (T , o, x b , x e ) = τ X|Y (T ) -τ X|Y (S) = k i j (t k ij ) 2 t k .. ×t k .j -k i (t k i. ) 2 (t k .. ) 2 N -k i (t k i. ) 2 (t k .. ) 2 - k i j (s k ij ) 2 s k .. ×s k .j -k i (s k i. ) 2 (s k .. ) 2 N -k i (s k i. ) 2 (s k .. ) 2 Setting µ k . = j µ k j , Ω = N -k i (t k i. ) 2 (t k .. ) 2 , and Γ = N -k i j (t k ij ) 2 t k .. ×t k .j
, we obtain the following updating formula:

∆ τ X|Y (T , o, x b , x e ) = = Ω k j 2µ k j t k .. ×t k .j [t k bj -t k ej -µ k j ] +Γ k 2µ k . (t k .. ) 2 [t k e. -t k b. +µ k . ] Ω 2 -Ω k 2µ k . (t k .. ) 2 (t k e. -t k b. +µ k . )
Author When partition Y k is modified: Let's now consider the variation of τ X|Y when one feature f is moved from the cluster y b to the cluster y e of partition Y k . This variation is the result of some changes on the co-occurrence table T k . Moving f induces, for each cluster of objects x i , the transfer of a quantity λ i from t k ib to t k ie of T k , the other contingency tables been unchanged. The overall vector of transferred quantities is called λ = [λ 1 , . . . , λ m ]. In the following, we call t k ij the elements of the co-occurrence table T k before the update operation, and s k ij the same elements after the operation. We call T = {T 1 , . . . , T N } and S = {S 1 , . . . , S N } the two sets of contingency tables, and τ X|Y (T ) and τ X|Y (S) the values of τ X|Y computed over T and S.

The following equations hold for t k ij and s k ij :

s ℓ ij = t ℓ ij , if ℓ = k or ℓ = k, j = b, j = e s k ib = t k ib -λ i s k ie = t k ie + λ i
Therefore, the variation of τ X|Y is:

∆ τ X|Y (T , f, y b , y e ) = τ X|Y (T ) -τ X|Y (S) = k i j (t k ij ) 2 t k .. ×t k .j -k i (t k i. ) 2 (t k .. ) 2 N -k i (t k i. ) 2 (t k .. ) 2 - k i j (s k ij ) 2 s k .. ×s k .j -k i (s k i. ) 2 (s k .. ) 2 N -k i (s k i. ) 2 (s k .. ) 2
Setting Ω = Nk i (t k i. ) 2 (t k .. ) 2 and λ = i λ i , we obtain the following updating formula:

∆ τ X|Y (T , f, y b , y e ) = 1 Ω i (t k ie ) 2 t k .e t k .. - (t k ie + λ i ) 2 (t k .e + λ)t k .. + (t k ib ) 2 t k .b t k .. - (t k ib -λ i ) 2 (t k .b -λ)t k ..
Thanks to this approach, instead of computing τ X|Y at lines 9 and 12 of Algorithm 5 with a complexity in O(m • N ), we compute ∆ τ X|Y (T , f, y b , y e ) once at line 9 with a complexity in O(m) (that is in O(|O|) in the worst case with the discrete partitions X) and test if it is strictly positive in the first if condition or strictly negative in the second one. Computing Ω is in O(m • N ) and is done only once at the beginning of Algorithm 5.

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z

Overall complexity

In this subsection, we evaluate the overall complexity of CoStar and compare it to the complexity of competitors used in Section 5 (SRC [START_REF] Long | Spectral clustering for multi-type relational data[END_REF], NMF [START_REF] Chen | Semi-supervised document clustering with simultaneous text representation and categorization[END_REF] and ComRaf [START_REF] Bekkerman | Multi-modal clustering for multimedia collections[END_REF]).

CoStar starts by computing N contingency tables T k . This operation takes O(m • n k ) time for each feature space Let I denote the total number of iterations. The overall time complexity of CoStar is O(I • m • N ), i.e., it is linear in the number of objects and features.

F k , k = 1 • • • N .
We now consider the theoretical complexity of the competitors, as reported by the study [START_REF] Chen | Non-negative matrix factorization for semisupervised heterogeneous data coclustering[END_REF]. We use again I as the number of iterations, N as the number of feature spaces, and m as the number of objects, while n is the maximum feature dimension for all feature spaces and k is the maximum between the number of object clusters and the number of feature clusters. Hence, SRC has a complexity in O(I •N •(max(m, n) 3 +kmn)) and ComRaf has a complexity in O(I • N • max(m 3 , n 3 )). Notice that the complexity of these two algorithms is cubic in the biggest of the dimensions of objects and all the feature spaces. Complexity of NMF, instead, is similar to our method's: this approach takes time O(I

• N • k • m • n).
However, this method requires to define the number of clusters as a parameter.

Experiments

In this section, we empirically show the effectiveness of our algorithm CoStar in comparison with three state-of-the-art algorithms that deal with star-structured data. All of them are non deterministic, have been presented in recent years and have been summarized in Section 2. The first competitor is Comraf and it is presented in [START_REF] Bekkerman | Multi-modal clustering for multimedia collections[END_REF]. The second competitor is SRC and is described in [START_REF] Long | Spectral clustering for multi-type relational data[END_REF]. The third one is NMF and is introduced in [START_REF] Chen | Semi-supervised document clustering with simultaneous text representation and categorization[END_REF]. In particular, for obtaining all the experiment results of the competitors, we took care to use the same pre-and post-processing methods that have been adopted and described in the original papers. For instance, since the second approach uses the k-means algorithm as a post-processing method, we did the same.

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z All these methods require as input, at least, the number of row clusters. For any given dataset, we set this parameter equal to the number of preexisting classes. To validate our proposal, we use three different groups of datasets (details are given further in the text):

1. a publicly available synthetic dataset (BBC and BBCSport) for standard analysis, 2. three document-word-category datasets as illustrative cases of imbalanced problems, 3. an image-word-blob dataset (COREL Benchmark ) as illustrative cases of a balanced problem.

As any of these methods (the three competitors and CoStar) is non deterministic, for any given dataset, we run each algorithm 100 times and compute the average and standard deviation of the evaluation measures used (presented below). Since CoStar does not require to specify the number of clusters, we also study the distribution of the number of row clusters over the 100 trials for this algorithm. The number of iterations for CoStar was set to 10 × max(m, N ) (m and N being the number of objects and features respectively). As further information, CoStar is written in C++, Comraf in C, SRC in R, and NMF in MATLAB. Finally, all the experiments have been performed on a 2.6GHz Opteron processor PC, with 4GB RAM, running Linux.

External evaluation measures

We evaluate the performance of CoStar and the competitors using two external validation indices. We denote by C = {C 1 . . . C J } the partition built by the clustering algorithm on objects, and by P = {P 1 . . . P I } the partition inferred by the original classification. J and I are respectively the number of clusters |C| and the number of classes |P|. We denote by n the total number of objects. The first index is the Normalized Mutual Information (NMI). NMI provides an information that is impartial with respect to the number of clusters [START_REF] Strehl | Cluster ensembles -a knowledge reuse framework for combining multiple partitions[END_REF]. It measures how clustering results share the information with the true class assignment. NMI is computed as the average mutual information between every pair of clusters and classes:

NMI = I i=1 J j=1 x ij log nxij xixj I i=1 x i log xi n J j=1 x j log xj n
where x ij is the cardinality of the set of objects that belong to cluster C j and class P i ; x j is the number of objects in cluster C j ; x i is the number of objects in class P i . Its values range between 0 and 1.

The second measure is the Adjusted Rand index [START_REF] Hubert | Comparing partitions[END_REF]. Let a be the number of object pairs belonging to the same cluster in C and to the same class in P. This metric captures the deviation of a from its expected value corresponding

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z to the hypothetical value of a obtained when C and P are two random independent partitions. The expected value of a denoted by E[a] is computed as follows:

E[a] = π(C) • π(P ) n(n -1)/2
where π(C) and π(P ) denote respectively the number of object pairs that belong to the same cluster in C and to the same class in P. The maximum value for a is defined as:

max(a) = 1 2 (π(C) + π(P ))
The agreement between C and P can be estimated by the adjusted rand index as follows:

ARI(C, P) = a -E[a] max(a) -E[a]
Notice that this index can take negative values, and when ARI(C, P) = 1, we have identical partitions.

Datasets

Here, we describe in detail the datasets we employed for our experimental analysis. The first two datasets are multi-view text data from news corpora, downloaded from http://mlg.ucd.ie/datasets/segment.html: BBC and BBCSport. The multiple views were created by splitting the text corpora into related segments. Each segment is constituted by consecutive textual paragraphs as described in the cited website. Furthermore, as explained by the web site, the standard textual pre-processing steps have been performed on the texts: stemming, stop-words removal and removal of words occurring very few times (lower than 3 times). The BBCSport dataset consists in 2 views (bbcsport2), while the BBC dataset is available both with 2 and 4 views (bbc2 and bbc4). For each collection we selected only the documents that contain at least one segment for each view. The characteristics of the datasets are reported in Table 6. Object classes are provided with the datasets and consist in annotated topic labels provided for the news articles (with values such as business, entertainment, politics, sport and tech).

Also the second group of datasets comes from a textual domain and is described in [START_REF] Forman | An extensive empirical study of feature selection metrics for text classification[END_REF]. In particular, in our experiments we consider combinations of the following datasets:

-oh15 : is a sample from OHSUMED dataset. OHSUMED is a clinicallyoriented MEDLINE subset of abstracts or titles from 270 medical journals over a five-year period (1987-1991). -re0 : is a sample from Reuters-21578 dataset. This dataset is widely used as test collection for text categorization research. -wap: is a sample from the WebACE Project, consisting on web pages listed in the subject hierarchy of Yahoo!.

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, The first view consists in sets of terms contained in the documents, also available on the Weka Website 1 . Also for this data, we performed the same standard pre-processing used for the news corpora: stemming, stop-words removal and removal of words occurring less than 3 times. As a consequence of the preprocessing, we obtained a document-word matrix.

For the second view, we consider some categorization of the documents. To this purpose, we used the semi-supervised technique presented in [START_REF] Chen | Semi-supervised document clustering with simultaneous text representation and categorization[END_REF]. Each feature category of this view corresponds to a category provided by the documentset and depends on the dataset the document belongs to. Then, the value of each feature for each document is compiled as the probability that the document belongs to the related category: it is computed as the portion of the whole vocabulary associated to a specific category which is also present in a specific document. Each element of the document-category matrix has values in the range [0, 1]. The details of these datasets are reported in Figure 7 and in Figure 8. In particular, we generated three new datasets of increasing difficulty (named T1, T2 and T3) in which we created the concept of super-class by injecting an increased separation only between some of the classes. The aim is to see if the co-clustering is able to recognize the concept of super-class. We did the following: we took the object classes from the three datasets ( OHSUMED, re0 and wap) and mixed them as the information in Figure 7 suggests. For instance, for the generation of dataset T1 we placed the first three classes of objects indicated in the last column of Figure 7 The last group of datasets are extracted from the Corel Benchmark 2 . This benchmark has already been used in [START_REF] Bekkerman | Multi-modal clustering for multimedia collections[END_REF], and we use a similar pre-processing. It consists of 5000 images from 50 Corel Stock Photo CDs. Each CD contains 100 images on the same topic. The original dataset has 4500 training images and 500 test images. Here, we consider only the 4500 training images for our experiments. Every image has a caption and an annotation. The caption is a brief description of the scene and the annotation is a list of objects (entities, or segments) that appear in the image. On overall, 371 words and 42379 entities are used to describe the collection. We built two different views: the first view is based on the caption words, the second view is based on blobs. A blob is a portion of the image in which a user has recognized a particular entity. It is similar to a concept. In order to obtain blobs, we clustered (with a k-means algorithm) all the entities into 3000 groups (blobs). Each image is mapped onto the set of blobs or entities represented in the image itself. This model leads to a representation analogous to the bag-of-words model (BOW) often used in text processing. For our experiments, we built three different datasets (details are reported in Figure 9). Each of the generated dataset has balanced classes represented by an average of 90 objects.

Results and discussion

We now analyze in detail the different aspects of our experimental evaluation. First, we compare the behavior of our algorithm to the three competitors. Second, we assess the significance of the results statistically. Third, we study the variability of the number of detected clusters. Then, we perform an in-depth qualitative analysis of the specific results. Finally, we provide a scalability study.

CoStar against competitors

In Figure 10 and 11, we report the average and standard deviation of the performance indexes for each dataset and each clustering approach. We can observe that CoStar outperforms all other approaches on NMI index and achieves the best results on all the datasets except I3 when considering the ARI index. We recall that CoStar does not require any parameter. Despite this, it finds always high-quality partitions, while the other approaches do not take advantage of the additional information of the correct number of clusters. Notice that CoStar always outperforms SRC. On the other hand, Comraf achieves comparable results over the Images datasets. This is not surprising, since Comraf is suited for this specific application domain. But, when applying Comraf on the text datasets, it returns very poor results regarding both measures. NMF performs better than Comraf on text data, and has similar results on Corel Benchmarks.

Statistical significance of the results

To assess the statistical quality of our approach we use the Friedman statistic and the Nemenyi test [START_REF] Demsar | Statistical comparisons of classifiers over multiple data sets[END_REF]. These techniques are usually employed to evaluate the statistical relevance of results of different classifiers over multiple datasets. We briefly summarize the Friedman test:

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z 1. the performance of each method (in terms of a given evaluation measure such as ARI or NMI) is determined on each dataset; 2. the methods are ranked for each dataset according to these results, the best performing algorithm getting the rank of 1, the second best rank 2, and so on ; 3. for each method, its average rank R j w.r.t. the datasets is computed; 4. finally, the Friedman statistic is computed thanks to the following formula:

Q = 12N k(k + 1)   k j=1 R 2 j - k(k + 1) 2 4  
where N is the number of datasets and k is the number of methods that are compared.

On the null hypothesis, where is no difference among the k sets of measures, the average rank should approximate (k+1) 2 . When N and k are large (i.e. N > 15 and k > 4), the probability distribution of Q can be approximated by that of a chi-square distribution with k -1 degrees of freedom. In this case the p-value is given by P (χ 2 k-1 ≥ Q). If N or k is small, the approximation to chi-square becomes poor and the p-value should be obtained from tables of Q specially prepared for the Friedman test.

We compare ARI and NMI results of CoStar and its competitors (SRC, ComRaf and NMF) on the 9 datasets and compute the associated Q values: we obtain 21.63 and 22.73 for ARI and NMI respectively. As the critical value at significance level α = 0.01 is equal to 10.73 (see [START_REF] Zar | Biostatistical Analysis[END_REF][START_REF] Sheskin | Handbook of Parametric and Nonparametric Statistical Procedures[END_REF]) and is smaller than the Q values, we decide to reject the null hypothesis for both measures.

We can now proceed with the post-hoc Nemenyi test. According to this test, the performance of two algorithms is significantly different if the corresponding average ranks differ by at least the critical difference:

CD = q α k(k + 1) 6N
where α is the significance level (α = 0.10 in our case) and q α is the critical value for the two tailed Nemenyi test [START_REF] Demsar | Statistical comparisons of classifiers over multiple data sets[END_REF]. We compare the 4 methods at the critical value q 0.10 = 2.291. The ranking table is shown in Figure 12 (top). The critical difference between average ranks is CD = 1.394 (for ARI and NMI). Figure 12 (bottom) shows the difference between average ranks for both ARI and NMI. In the first case, CoStar is significantly different from ComRaf and SRC. In the second case, it performs significantly better than all other competitors.

Number of detected clusters

Since CoStar may find different numbers of clusters at each execution, we ran the algorithm 100 times and evaluate the variability of this result. In Figure 13, 14 and 15, we report the distributions of the number of row clusters for the three groups of datasets. On bbc2 and bbcsport2, CoStar produces most of the times partitions with the correct number of clusters (see Figure 13). However, the resulting partitions on bbc4 contain 5 or 6 clusters to the same extent. Now consider the three textual datasets (see Figure 14). By analyzing the results, we found that CoStar recognizes the "super-classes" corresponding to the originating dataset to which the document belongs to. For instance, for T1 there are two super-classes: one from oh15 and one from re0. This is reasonable because they represent two types of dictionary: a biology-centered one and a news-centered one. The same phenomenon happens for the other two textual datasets. Furthermore, the found number of clusters remains quite stable, in particular for T1 and T2.

Author
Figure 15 shows the distribution of the number of row clusters for the Corel benchmark. First, CoStar detects most of the time the correct number of clusters [START_REF] Chen | Non-negative matrix factorization for semisupervised heterogeneous data coclustering[END_REF] when applied to I1. When processing I2, the algorithm finds three different values of cluster numbers for the row partition: 5, 6 and 7. Notice that in this the correct number is 6, which is very close to the average of the detected number of clusters. For I3, which is a hard context for a clustering algorithm, the distribution of the number of detected clusters may be approximated with a normal distribution with a mean of 9, while the number of built-in classes is 13.

A qualitative evaluation

We now show an example to provide some insights about the kind of results our algorithm is able to provide. To this purpose, we choose a specific result of the I2 image dataset. In particular, among the 100 runs of CoStar, we choose the one providing the result with the highest value on the objective function computed over the row space (τ X|Y ).

To provide a readable description of each cluster, we perform a similar process to the one used in CoStar when optimizing the objective functions. The goal is to select a subset of elements within each cluster that represents and summarizes the content of the whole cluster. For any given cluster, we try to remove each of its elements and compute the ∆τ (see Section 4) that corresponds to this removal. In this case we use δτ as an evaluation function of the single cluster elements. We assign to each object the ∆τ (maximum in absolute value) that as a consequence of the movement causes the highest reduction of the objective function value. Finally, we rank all the elements by the obtained values of ∆τ , and select the top k objects as representatives. Notice that we use ∆ τ X|Y to rank images, and ∆ τ Y k |X to rank words (see Section 4). We report the cluster descriptions for the textual view of I2 in Figure 16. The reader can observe that the algorithm returns five well-separated clusters of words. The first cluster gathers words around polar animals; the second one is about bridges, coasts and waterways; the topic of the third one is around

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z animal with horns, the fourth one talks of horses and the last one is about birds and places where they live.

bear, polar, black, snow, antlers, grizzly, elk, ice, tundra, grass bridge, water, coast, arch, hills, sky, waves, beach, boats, steel deer, white-tailed, mule, horns, slope, fawn horses, foals, mare, field, fence, bush birds, nest, branch, fly, tree, leaf, wood, wings, stick, baby Using the top 5 ranked images, we distinguish five clusters: Clusters1 is about animals (reindeer and bears) that live in tundra or polar climate regions; Clusters2 contains pictures of birds; Cluster3 is about horses; Cluster4 is about coasts and bridges; finally, Cluster5 contains pictures of deer. Notice also that in this particular case it is possible to map each cluster of images onto a cluster of words and vice versa. Fig. 17 I2 dataset with 5 representative images for each cluster.

Scalability evaluation

We now analyze the time performances of CoStar. We did not perform comparisons with the competitors since all algorithms are written in different programming languages, but the reader may refer to Section 4.4 for a theoretical complexity comparison. Here, we focus our analysis on how the number of fea-Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z ture spaces and the number of objects influence the computational time (see Figure 18).

For the first experiment, we generate several synthetic boolean datasets varying the number of feature spaces. Each feature space contains 1000 features. We let the number of feature spaces vary from 2 to 10. The number of row objects varies from 500 to 2 000. Data are generated by considering a uniform distribution and an average matrix density of 5%. The related time curves are depicted in Figure 18(a). We observe that CoStar has reasonable time performances on high-order and high dimensional data. As expected, the behavior is linear w.r.t. the number of feature spaces.

The last experiment allows us to study the scalability of our algorithms when applied to large datasets. In this case, the number of spaces is fixed (2), while we let the number of objects vary from 1 000 to 10 000. Furthermore, we repeated this experiment varying the number of features per space from 500 to 2 500. As shown in Figure 18(b), even though the complexity of CoStar is quadratic in the number of objects/features, in practice, the computational time increases linearly with the number of objects.

Empirical Convergence over Synthetic Data

In this subsection we empirically evaluate the convergence of our method over a synthetic dataset where we know in advance the inner structure. In particular we consider a Boolean dataset D over O (1 000 objects) and F = {F 1 , F 2 } (1000 features for each views). The synthetic dataset has a perfect block nested structure with perfect partitions C O1 × C F 1 and C O2 × C F 2 . We have built C O1 × C F 1 to have two blocks structure (two perfect clusters over the objects and two perfect clusters over the features). In the same way C O2 × C F 2 has a four blocks structure (four perfect clusters over the objects and four perfect clusters over the features).

Let us consider a parameter p that is the probability to flip a Boolean value uniformly at random in D and denote by D p the resulting dataset. In our experiment we vary p over 5%, 10%, 15% and 20%.

In order to evaluate the quality of the non-dominated solutions, we follow the protocol proposed in [START_REF] Knowles | A Tutorial on the Performance Assessment of Stochastic Multiobjective Optimizers[END_REF] and [START_REF] Liefooghe | On dominancebased multiobjective local search: design, implementation and experimental analysis on scheduling and traveling salesman problems[END_REF]. Given a dataset D p we perform 30 runs. Then we define Z all as the union of the obtained solutions. Note that Z all could contain both dominated and non-dominated vectors, since a solution may dominate another. For this reason we define Z ⋆ as the set of non-dominated vectors of Z all . Second, we define z min = (z min ) denotes the lower (resp. upper) bound of the k th objective for all the points contained in Z all . In order to give a roughly equal range to the objective functions, the values are normalized with respect to z min and z max ( value-z min z max -z min ). Note that after the normalization step z min is equal to 0 N +1 and z max is equal to 1 N +1 .

Given a solution S, in order to measure the quality of S in comparison to Z ⋆ , we compute the difference between these two sets by using the unary hypervolume metric [START_REF] Knowles | A Tutorial on the Performance Assessment of Stochastic Multiobjective Optimizers[END_REF], where the vector z min = 0 N +1 is the reference point, as it corresponds to the lower bound of the N + 1 normalized objective functions. The hypervolume difference indicator (I H ) computes the portion of the objective space that is weakly-dominated by Z ⋆ and not by S. The closer this

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z measure to 0, the better the approximation S. At each iteration we compute the average of hypervolume difference indicator I H of the objective functions vectors S of the 30 runs. The maximum number of iterations is set to 10000 using the general rules employed in the experiments (10 × max(m, N ) where m and N are the number of objects and features respectively). Using this measure, we evaluate the volume of the dominated portion of the objective space and the convergence to a stable solution of CoStar. This process is repeated for each value of p. The results are reported in Figure 20.

In Figure 20, we observe that CoStar converges always to very similar solutions in all the runs with different values of p. As we state in Section 4.2, we can notice that at each iteration CoStar converges to a Pareto local optimum with respect the neighborhood. It is interesting to notice that, in this case, the algorithm converges to a stable solution very early w.r.t. the maximum number of iterations N iter foreseen in CoStar. We empirically observe the same phenomenon on the real datasets, too. In particular we have observed that the maximum number of iterations is always an upper bound of the real number of steps needed by CoStar to obtain the final stable solution.

Conclusion

Heterogeneous star-structured data are very common in real world applications. In this paper we presented CoStar, a novel algorithm that deals with this kind of data. Unlike previous approaches, CoStar is totally parameterless: this means that it does not require either the number of object clusters or the number of feature clusters for each of the feature spaces. Our approach employs a Pareto-dominance approach to optimize a set of measures for crossassociation in contingency tables which we extended to a high-dimensional setting. We assessed the performance of CoStar using objective evaluation measures on publicly available datasets over both the textual and image domains. We compared CoStar with three competitors, evaluate the statistical significance of the obtained results, studied the variability of the number of detected clusters and performed an in-depth qualitative analysis on an image dataset. Finally we provided a scalability study. The results show that our approach outperforms state-of-the-art methods for heterogeneous data (co-)clustering.

Fig. 1

 1 Fig. 1 An example dataset (a) and the contingency table associated with a related coclustering (b).

Fig. 2

 2 Fig. 2 An example dataset (a) and the contingency tables that correspond to C 1 (b) and C 2 (c).

Fig. 5

 5 Fig. 5 An example dataset (a) and two contingency tables associated with the starstructured co-clusterings C 1 (b) and C 2 (c).

1 :

 1 Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26[START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z for simplicity, we bound the number of iterations by N iter in the current version of the algorithm. Algorithm 2 OptimizeMultiObjectCluster(X, Y, T ) Randomly choose a cluster x b ∈ X 2: Randomly choose an object o ∈ x b 3: N D ← {x b } 4: maxτ X|Y ← τ X|Y (T ) 5: for all xe ∈ {X ∪ ∅} s.t. x b = xe do 6: [X ′ , T ′ ] ← Update(X, T, o, x b , xe) 7:

  T k Algorithms 4 and 5 optimize partition Y k and work in a similar way than Algorithms 2 and 3. However, it is simpler because (1) modifying partition Y k only impacts the contingency table T k′ , and (2), in Algorithm 5, the comparison of two partitions of the queue N D is only done on the basis of τ X|Y value, the others measures (τ Y k |X , k = 1, . . . N ) being equal.

  Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26(2), 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z Algorithm 5 FindParetoLocalOptimumFeature(N D, Y k , T k , y b , f ) 1: while Card(N D) > 1 do 2:

Theorem 1

 1 {w b \ {u}} ∪ {w e ∪ {u}} and w e ∈ {W i ∪ ∅}. The iterated local search algorithm CoStar terminates in a finite number of steps and outputs a Pareto local optimum with respect to N (i,w b ,u) .

Therefore, r

  is a Pareto local optimum with respect to s. Both s and r belong to N D and Algorithms FindParetoLocalOptimumObject or FindParetoLocalOptimumFeature will identify that s is dominated by r and will remove s from N D: (a) If i = N + 1, the set N D is processed by the function FindPare-toLocalOptimumObject, and s is removed from N D since N k=1 δ z k (r)>z k (s) > N k=1 δ z k (s)>z k (r) with δ C = 1 if C = T rue 0 otherwise (b) If i = N + 1, the set N D is processed by the function FindPare-toLocalOptimumFeature, and s is removed from N D since

4. 3

 3 Optimized computation of τ Y k |X and τ X|Y Algorithms 2 and 4 modify, at each iteration, partitions X and Y k , by evaluating functions τ X|Y and τ Y k |X respectively. Algorithms 3 and 5 find a non-dominated partition X and Y k , among a set of candidates, based on the computation of τ Y k |X and τ X|Y respectively. The computational complexity of these functions is in O(m • N ) and O(m • n k ) where N = k n k , |Y k | = n k and |X| = m. In the worst case, these complexities are in • k |F k |) and O(|O| • |F k |) where |O| represents the cardinality of the objects dimension. Moreover, during each iteration of algorithms 2 and 4, these operations are performed for each cluster (including the empty cluster). Thus the overall complexities are in O(m•m•N ) and O(m•n k •n k ).

  Using the notation introduced beforehand, computing all the contingency tables takes O(m•N ). The core of the algorithm is the iterative procedure which updates the object partition and the N feature partitions. Algorithm 2 has a complexity in O(m•N ) for computing Γ and O(m•N ) for computing Ω; the loop takes time O(m • N ) and Algorithm 3 takes O(m • n k ) in the worst case (|N D| = m). Therefore the whole complexity of the function is in O(m • N ). Concerning Algorithm 4, computing Γ and Ω takes time O(m • n k ) and O(m) respectively; the loop takes O(m • n k ) and Algorithm 5 takes time O(m • n k ) in the worst case (|N D| = n k ). The overall complexity of the function is in O(m • n k ).

Fig. 7

 7 Fig. 7 Text dataset characteristics.

  (Aden-Diph, Cell-Mov, Aluminium) from the original dataset OHSUMED (indicated in Figure 7 by symbol oh15) while the remaining 2 classes (cpi, money) have been taken from Reuters-21578 (indicated in Figure by re0). Analogously for the remaining datasets T2 and T3. Figure 8 reports instead the number of objects for each document class.

Fig. 9

 9 Fig. 9 Image dataset characteristics.

Fig. 13 Fig. 14

 1314 Fig. 13 Distribution of CoStar number of clusters on (a) bbc2 (b) bbcsport2 (c) bbc4.

Fig. 15

 15 Fig. 15 Distribution of the number of clusters of CoStar on (a) I1 (b) I2 (c) I3.

Fig. 16

 16 Fig. 16 Word clusters for image dataset 2.

Figure 17

 17 Figure17resumes the descriptions of the clusters obtained over the images. Using the top 5 ranked images, we distinguish five clusters: Clusters1 is about animals (reindeer and bears) that live in tundra or polar climate regions; Clusters2 contains pictures of birds; Cluster3 is about horses; Cluster4 is about coasts and bridges; finally, Cluster5 contains pictures of deer. Notice also that in this particular case it is possible to map each cluster of images onto a cluster of words and vice versa.

1 ,

 1 . . . , z min N +1 ) and z max = (z max 1 , . . . , z max N +1 ), where z min k (resp. z max k

Fig. 18

 18 Fig. 18 Running time of CoStar varying the number of feature spaces and objects

Figure 19

 19 illustrates this measure. The light gray area corresponds to the hypervolume difference indicator.

Fig. 19

 19 Fig. 19 Illustration of the hypervolume difference indicator.

Fig. 20

 20 Fig. 20 The average Hypervolume Difference Indicator at each iteration at different values of p. (a) p=5%, (b) p=10%, (c) p=15% and (d) p=20%

  Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26[START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z

						f 2 1	f 2 2	. . .	f 2 n2
					O1 d 2 11 d 2 12	. . .	. . .
					O2 d 2 21
					.	.
					.	.
					.	.
					OM d 2 M1 . . .	. . .	d 2 Mn2
			Feature 3						Feature 1
		f 3 1	f 3 2	. . .	f 3 n3		Objects		f 1 1	f 1 2	. . .	f 1 n1
	O1 d 3 11	d 3 12	. . .	. . .		O1		O1 d 1 11	d 1 12	. . .	. . .
	O2 d 3 21							O2 d 1 21
	. .	. .					O2		. .	. .
	.	.					...		.	.
	OM d 3 M1 . . .	. . .	d 3 Mn3		OM		OM d 1 M1 . . .	. . .	d 1 Mn1

  + xi the number of times τ Y k |Xi is strictly greater than τ Y k |Xj , and in nb + xj the number of times τ Y k |Xj is strictly greater than τ Y k |Xi , for k = 1, . . . , N . If nb + xi ≥ nb + xj then X j is either dominated by X i or neither X i nor X j dominate the other (see proof of theorem 1). Thus, Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26(2), 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z

  then s dominates r: s ≻ r.

	Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26(2), 217-254.
	The original publication is available at http://www.springer.com
	DOI : 10.1007/s10618-012-0248-z

  -produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26(2), 217-254.

	The original publication is available at http://www.springer.com
	DOI : 10.1007/s10618-012-0248-z

  produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26(2), 217-254.

	The original publication is available at http://www.springer.com
	DOI : 10.1007/s10618-012-0248-z

  Thanks to this approach, instead of computing τ Y k |X at lines 7 and 11 of Algorithm 4 with a complexity in O(m• n k ), we compute ∆ τ Y k |X (T k , f,y b , y e ) once at line 7 with a complexity in O(m) (that is in O(|O|) in the worst case with the discrete partition X) and test if it is strictly positive in the first if condition or null in the second one. Computing Γ is in O(m • n k ) and Ω in O(m) and is done only once at the beginning of Algorithm 4.When partition X is modified: Let's now consider the variation of τ Y k |X when one element o is moved from cluster x b to cluster x e of partition X. This operation induces, for each cluster of features y k j , the transfer of a quantity µ k j from t k bj to t k ej of T k . The overall vector of transferred quantities is calledµ k = [µ k 1 , . . . , µ k n k ].We call t k ij the elements of the k th contingency table T k before the update operation, and s k ij the same elements after the operation. τ Y k |X (T k ) and τ Y k |X (S k ) are the values of τ Y k |X computed over T k and S k . The following equations hold for t k ij and s k ij :

	(t k .e -t k .b + λ)

  -produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26(2), 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z

  -produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26(2), 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z Thanks to this approach, instead of computing τ X|Y at lines 7 and 11 of Algorithm 2 with a complexity in O(m•N ), we compute ∆ τ X|Y (T , o, x b , x e ) once at line 7 with a complexity in O(N ) (that is in O( k |F k |) in the worst case with the discrete partitions Y k ) and test if it is strictly positive in the first if condition or null in the second one. Computing Γ is in O(m • N ) and Ω in O(m • N ) and is done only once at the beginning of Algorithm 2.

  , 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z Synthetic dataset characteristics.

		Name	# obj. # obj. classes	View	# features
		bbc2	2012		5		1st view 2nd view	6838 6790
		bbcsport2	544		5		1st view 2st view	3183 3203
								1st view	4659
		bbc4		685		5		2nd view 3rd view	4633 4665
								4th view	4684
	Fig. 6 Name	Datasets	# terms	# cat.	# obj.	# classes	Classes
	T1	oh15 re0	3987		2	833	5	{Aden-Diph, Cell-Mov, Aluminium} {cpi, money }
	T2	oh15 re0	3197		2	461	5	{Blood-Coag, Enzyme-Act, Staph-Inf} {jobs, reserves}
		wap						{ Film, Television, Health}
	T3	oh15	8282		3	2129	9	{Aden-Diph, Cell-Mov, Enzyme-Act}
		re0						{ interest, trade, money}

  Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26[START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z Classes distributions for text datasets.

			Class Name # obj. Class Name	# obj.
			Aden-Diph	56	Cpi	60
			Cell-Mov	106	Money	608
			Aluminium	53	Interest	219
			Blood-Coag	69	Trade	319
			Enzyme-Act	154	Film	196
			Staph-Inf	157	Television	130
			Jobs		39	Health	341
			Reserves	42
	Fig. 8 Name	# words	# obj.	# classes		Classes
						sunsets, tigers, train,
	I1	114	630	7		swimmers, formula One car,
						skyscrapers, war airplanes
	I2	116	541	6		bears, deers, horses, cliffs, birds, bridges
						sunsets, tigers, train, swimmers,
	I3	170	1171	13	formula One car, skyscrapers, war airplanes, bridges,
						bears, deers, horses, cliffs, birds,
	1 http://www.cs.waikato.ac.nz/ml/weka/

  Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26[START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF], 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z

		CoStar	Comraf	SRC	NMF
	bbc2	0.68 ± 0.02	0.58 ± 0.04	0.31 ± 0.11	0.5 ± 0.04
	bbcsport2	0.69 ± 0.06	0.09 ± 0.05	0.24 ± 0	0.57 ± 0.01
	bbc4	0.68 ± 0.03	0.41 ± 0.03	0.33 ± 0.04	0.45 ± 0.0
	I1	0.88 ± 0.04	0.85 ± 0.05	0.80 ± 0.11	0.86 ± 0.06
	I2	0.75 ± 0.04	0.66 ± 0.05	0.61 ± 0.09	0.74 ± 0.04
	I3	0.76 ± 0.03	0.75 ± 0.04	0.61 ± 0.11	0.73 ± 0.06
	T1	0.72 ± 0	0.50 ± 0.03	0.23 ± 0.05	0.49 ± 0.01
	T2	0.73 ± 0.1	0.35 ± 0.03	0.33 ± 0.08	0.58 ± 0.02
	T3	0.71 ± 0.02	0.66 ± 0	0.55 ± 0.06	0.63 ± 0.02
	Fig. 10 Normalized Mutual Information results.
		CoStar	Comraf	SRC	NMF
	bbc2	0.67 ± 0.03	0.42 ± 0.05	0.18 ± 0.08	0.48 ± 0.06
	bbcsport2	0.58 ± 0.14	0.23 ± 0.04	0.14 ± 0.09	0.56 ± 0.02
	bbc4	0.56 ± 0.2	0.34 ± 0.05	0.16 ± 0.06	0.41 ± 0.0
	I1	0.83 ± 0.09	0.80 ± 0.08	0.69 ± 0.11	0.79 ± 0.1
	I2	0.72 ± 0.07	0.55 ± 0.06	0.51 ± 0.09	0.67 ± 0.07
	I3	0.50 ± 0.11	0.62 ± 0.06	0.43 ± 0.13	0.61 ± 0.09
	T1	0.73 ± 0	0.23 ± 0.05	0.09 ± 0.09	0.23 ± 0.02
	T2	0.68 ± 0	0.19 ± 0.04	0.18± 0.09	0.48 ± 0.01
	T3	0.44 ± 0.02	0.42 ± 0.03	0.36 ± 0.08	0.44 ± 0.04
	Fig. 11 Adjusted Rand Index results.
	2 http://kobus.ca/research/data/eccv_2002

  -produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26(2), 217-254. The original publication is available at http://www.springer.com DOI : 10.1007/s10618-012-0248-z

					Average ranking
			Methods	ARI	NMI
			CoStar		1.2	1
			NMF		2.2	2.4
			ComRaf		2.6	2.7
			SRC		4	3.9
	ARI	NMF ComRaf	SRC	NMI	NMF ComRaf	SRC
	CoStar	1	1.4	2.8		CoStar	1.4	1.7	2.9
	NMF	-	0.4	1.8		NMF	-	0.3	1.5
	ComRaf	-	-	1.4		ComRaf	-	-	1.2
	SRC	-	-	0		SRC	-	-	0
	Fig. 12 Average Ranking results for the Nemenyi Test (top); Difference between average
	ranks for ARI and NMI evaluation measures (bottom).

Acknowledgments

Céline Robardet is partially funded by the French Research Agency (ANR) through the FOSTER project (ANR-2010-COSI-012-02).

Author-produced version of the article published in Data Mining and Knowledge Discovery, 2013, 26 [START_REF] Banerjee | A generalized maximum entropy approach to bregman co-clustering and matrix approximation[END_REF]