{'fraction_non_alphanumeric': 0.0359023077227566, 'fraction_numerical': 0.01092506828167676, 'mean_word_length': 4.156970810369463, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 30, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Large-Scale Clustering Based on Data Compression

This paper considers the clustering problem for large data sets. We propose an approach based on distributed optimization. The clustering problem is formulated as an optimization problem of maximizing the classification gain. We show that the optimization problem can be reformulated and decomposed into small-scale sub optimization problems by using the Dantzig-Wolfe decomposition method. Generally speaking, the Dantzig-Wolfe method can only be used for convex optimization problems, where the duality gaps are zero. Even though, the considered optimization problem in this paper is non-convex, we prove that the duality gap goes to zero, as the problem size goes to infinity. Therefore, the Dantzig-Wolfe method can be applied here. In the proposed approach, the clustering problem is iteratively solved by a group of computers coordinated by one center processor, where each computer solves one independent small-scale sub optimization problem during each iteration, and only a small amount of data communication is needed between the computers and center processor. Numerical results show that the proposed approach is effective and efficient.


I. INTRODUCTION
In the recent years, due to the rapid progress of data acquisition and communication technologies, it has become readily easy to collect and store large amounts of data. Large databases of scientific measurements at the scale of terabyte or even petabyte can be frequently observed in high energy physics, astronomy, space exploration and human genome projects. Large databases of financial data and sale transactions at the scale of terabyte or petabyte can also be frequently observed. These huge amounts of data usually contain valuable scientific and business information. For example, a large collection of sale transaction data may contain important information of consumer behaviors and market trends. However, the data analysis on such large databases presents many technique challenges. The database size is usually far larger than the memory size of any single computer. Many existing centralized data analysis algorithms fail for these instances. In fact, most data analysis problems for large databases are currently open or not well-solved.
In this paper, we consider one important data analysis problem, the clustering problem for large databases. The clustering problem is the problem that a set of given data samples are classified into different groups, so that, the data samples within each group are similar according to certain metrics. Clustering is a fundamental problem in data analysis. It has many applications in pattern recognition, machine learning, data mining, computer vision, and signal processing. For example, clustering is usually an important step in many data mining algorithms.
Many algorithms for clustering problems have been previously discussed in the literature, see for example [1] and references therein. These algorithms range from heuristic algorithms to statistical modeling based algorithms. Among the previous algorithms, the statistical modeling based methods generally have better clustering performance compared with other types of algorithms, especially when the data clusters are not well separated. The Expectation-Maximization (EM) algorithms with mixture Gaussian modeling [2] [3] are the major state-of-the-art statistical modeling based clustering algorithms. The EM algorithms can be considered as iterative algorithms for computing the maximum likelihood estimation. It has been proven that the likelihood functions do not decrease during iterations.
However, it is well-known that the EM algorithms have certain limitations. First, according to previous experimental results, the EM algorithms may convergence very slowly [4], [5]. It is shown in [6], that the EM algorithms are first-order optimization algorithms, which provides a theoretical explanation for the slow convergence speeds. In fact, it has been a long-standing open problem that super-linear and secondorder methods should be found and preferred for the clustering problems [7]. Second, the EM algorithms do not converge and have numerical difficulties for certain types of instances [4], [8]. For example, the EM algorithms do not converge, when the covariance matrices are singular. The EM algorithms also do not converge, when the numbers of components in the mixture modeling are greater than the actual numbers of data clusters.
In addition, the standard EM algorithms require memory spaces proportional to the database size, therefore, do not scale well. Various scaling-up versions of the standard EM algorithms have been proposed in the literature [9], [10]. However, these previous approaches are approximation algorithms. The accuracy of the obtained results decreases as the ratio between the database size and main processor memory space size increases.
In this paper, we propose a new clustering algorithm for large databases based on data compression principles and mixture Gaussian modeling. Following the approaches in [11], we formulate the clustering problems as optimization problems. Instead of using a centralized approach, we propose a distributed algorithm to solve the global optimization problems. In our approach, the global optimization problem is decomposed into small-scale sub optimization problems using the Dantzig-Wolfe decomposition method [12]. Generally speaking, the Dantzig-wolfe method can only be used in the convex optimization case, where the duality gaps are zero. Even though, the considered problem in this paper is nonconvex, we show that the duality gap goes to zeros as the problem size goes to infinity. Therefore, the Dantzig-Wolfe method can be applied here. Our algorithm is especially suitable for the cases of distributed databases, where data are stored at multiple hosts or even at different geographical locations. The global optimal solutions can be computed with only intra-host computations, intra-host local database queries and a small amount of inter-host communications. Unlike many clustering algorithms for large databases, which compute approximate solutions, our algorithm computes exact solutions. Numerical results show that the proposed algorithm does not have any numerical difficulties for the case that the covariance matrices are singular. Numerical results also show that the algorithm has fast convergence speeds.
The rest of this paper is organized as follows. We present the proposed algorithm in Section II. We prove that the duality gap is vanishing for sufficiently large databases in Section III. Numerical results are presented in Section IV. We present the conclusion remark in Section V.
Notation: We use bold face lower-case letters and bold face capital letters to denote the column vectors and matrices respectively. For example, we use a a a to denote a column vector a a a. We use a a a(d) to denote the d-th element of the vector a a a. We use A A A t to denote the transpose of the matrix A A A. We use H(p 1 , . . . , p J ) to denote the entropy function, We use log(x) to denote the natural logarithm of the number x. We use det(A A A) to denote the determinant of the matrix A A A.

II. CLUSTERING ALGORITHM
In this paper, we consider a data set consisting of N data samples, x x x 1 , x x x 2 , . . . , x x x N , where each data sample is a D dimensional vector. We assume that the data samples are randomly distributed with a mixture Gaussian distribution. That is, Alternatively, we may consider x x x 1 , . . . , x x x n , . . . , x x x N as a mixture of data samples from J information sources, where each information source is Gaussian distributed. The considered problem is therefore estimating the membership of each data sample to one of the J information sources, and also the probability distribution of each information source.
In this paper, we propose a distributed algorithm for the above clustering problem. Our algorithm is efficient for the case that the data set contains a large amount of data samples. The data samples can be stored at multiple computers or database hosts. The proposed algorithm formulates the clustering problem as an optimization problem and decomposes the optimization problem into multiple small-scale sub optimization problems. Each sub optimization problem can be solved at one database host using only locally stored data samples. A center processor coordinates the computation at the database hosts. The final solution is obtained from the sub optimization results. A diagram of the system is shown in Fig. 1 The algorithm in this paper is built up on the data compression based algorithm for clustering in [11]. The main idea behind the algorithm is that optimal data clustering should induce optimal adaptive data compression. That is, if we partition the data set into several clusters and use one data compression encoder for each cluster, then the optimal compression performance is achieved if each cluster contains only the data samples from one information source. The algorithm in [11] then formulates the data cluster problem as an optimization problem, where the classification gain is maximized. The classification gain is a measure of data compression efficiency previously proposed in the data compression literature [13].
If the covariance matrices of all clusters are not singular, then the classification gain is inversely proportional to the following function, where, p i is the fraction of data samples in the i-th cluster, and Σ Σ Σ i is the covariance matrix of the i-th cluster. The above function is the objective function in our optimization formulation. In the sequel, we will always assume that the covariance matrices of all clusters are not singular without the loss of generality. Because, if any covariance matrix is singular, we can minimize the following function in the algorithm instead, where, σ 2 n is a sufficiently small positive number, and I I I D is the D dimensional identity matrix. This is equivalent to adding white noise with covariance matrix σ 2 n I I I D to the data samples and clustering the noise corrupted data samples instead. The optimality of the final obtained clustering results is not much affected, if σ 2 n is small enough.
The proposed algorithm formulates the clustering problem as an optimization problem. We introduce a variable a ni for each n, i, 1 ≤ n ≤ N , and 1 ≤ i ≤ J. The variable a ni is a likelihood that the n-th data sample belongs to the ith information source. The mean µ µ µ i , covariance matrix Σ Σ Σ i , and occurrence probability p i are functions of the likelihood variables a ni , The formulated optimization problem is therefore, Subject to: a a a ∈ Ω where, a a a is a vector obtained by stacking all the variables a ni , The final estimation results can be obtained by randomly rounding the optimal solution a * ni of the above optimization problem as in [11]. The near-optimality of this optimization based approach has been shown in [11] and [14].
In the sequel, we show that the optimization problem in Eqn. 8 can be reduced into sub optimization problems that can be locally solved at each database host. The reduction and reformulation procedure consists of four steps.
In the first step of reformulating the problem, we adopt an approach of first solving the restricted optimization problems with p i being fixed, Subject to: a a a ∈ Ω, and And then, we optimize over p 1 , . . . , p J to find the overall optimization solution, min p1,..., pJ g( p 1 , . . . , p J ) * , Subject to: The problem in Eqn. 11 can be easily solved by using the gradient descent approach. The main problem is therefore reduced to the optimization problem in Eqn. 10.
In the second step of reformulating the problem, we introduce auxiliary unitary matrices It can be shown that the optimization problem in Eqn. 10 is equivalent to the following optimization problem.
where, σ 2 id is the d-th diagonal element of the matrix B B B i . The two optimization problems are equivalent, because We solve the optimization problem in Eqn. 12 by an alternating optimization approach. That is, we iteratively first fix A A A 1 , . . . , A A A J and optimize over a a a, and then fix a a a and optimize over A A A 1 , . . . , A A A J . The latter optimization problem is easy to solve, because the optimal A A A 1 , . . . , A A A J are clearly the matrices, such that B B B i becomes diagonal. The main optimization problem is therefore reduced to Subject to: a a a ∈ Ω, where, A A A 1 , . . . , A A A J are fixed and given.
In the third step of reformulating the problem, we use an iterative upper bounding and minimizing approach to solve the optimization problem in Eqn. 14. Let σ 2 id [t] denote the solution obtained in the t-th iteration. Note that the objective function in Eqn 14 can be upper bounded as follows, due to the fact that the objective function is concave with respect to σ 2 id .
In the (t + 1)-th iteration, we find a solution a a a, such that the corresponding σ 2 id minimizes the above upper bound. It can be seen clearly that the objective function never increase during iterations. Therefore, the main optimization problem is reduced to the following optimization problem.
In the fourth and final step of reformulating the problem, we decompose the optimization problem in Eqn. 16 into sub optimization problems by using the Dantzig-Wolfe decomposition method. Each sub optimization problem can be locally solved at each database host. The Dantzig-Wolfe decomposition method is introduced initially for linear programming problems [12]. The method has been then generalized to the convex optimization cases, where the duality gaps are zero, (see for example [16] and references therein). For non-convex optimization problems, the decomposition method generally can not be applied due to the non-zero duality gaps. Even though the optimization problem in Eqn. 16 is not convex, we show in Theorem 3.6 that the duality gap goes to zeros as the number of data samples N goes to infinity. Therefore, the decomposition method can be applied here.
Let us assume that the data samples x x x 1 , . . . , x x x n , . . . , x x x N are stored at K database hosts. Let N k denote the set of the indexes of the data samples stored at the k-th host. We use A A A i x x x n (d) to denote the d-th element of the vector A A A i x x x n . The optimization problem in Eqn. 16 is equivalent to the following optimization problem. min a a a, µ µ µ Subject to: where, µ µ µ is the vector obtained by stacking all the variables µ ikd . The real number µ ikd can be considered as a local guess or estimation of the mean of A A A i x x x n (d) at the k-th database host. If all the local guesses are equal, then the above objective function is equal to the objective function in Eqn. 16. Because the duality gap is approximately zero as proven in Theorem 3.6, the optimization problem in Eqn. 17 is approximately equivalent to its Lagrangian dual problem as follows.
where, λ λ λ denotes the vector obtained by stacking all variables λ µikd and λ pi . The above optimization problem is separable and can be rewritten as, where, each f * k is the optimization result of one sub optimization problem. Let a a a k denote the vector obtained by stacking all variables a ni with n ∈ N k . Let µ µ µ k denote the vector obtained by stacking all parameters µ ikd , i = 1, . . . , J, d = 1, . . . , D.
It can be clearly checked that each f * k can be solved locally at each database host using only information about local data samples x x x n , n ∈ N k with given parameters β id , λ λ λ, and A A A 1 , . . . , A A A J . Therefore, the proposed algorithm iteratively computes the clustering result. During each iteration, each database host solves one local small-scale optimization problem as in Eqn. 20. The center processor then solves the global optimization problem as in Eqn. 19 using the local optimization results. The global optimization problem can be solved by using, for example, the subgradient method [16, Section 6.3.1].

III. VANISHING DUALITY GAP
In this section, we prove that the duality gap between the primal optimization problem in Eqn. 17 and the dual optimization problem in Eqn. 18 goes to zero as the problem size N goes to infinity. We need the Azuma inequality in our discussion. A proof of the inequality can be found, for example in [17] [18]. Lemma 3.1: (Azuma Inequality) Let Z 1 , . . . , Z N be independent random variables, with Z k taking values in a set Λ k . Assume that a (measurable) function f : Λ 1 ×Λ 2 ×· · ·×Λ N → R satisfies the following Lipschitz condition (L).
Then, the random variable X = f (Z 1 , . . . , Z N ) satisfies, for any t ≥ 0, The basic idea is to use randomization. Randomization has been used previously in establishing stronger duality theories. We refer interested readers to [19] and references therein. Let p(a a a, µ µ µ) denote the probability distribution of a a a and µ µ µ, where the range of a a a is Ω, and We introduce the following randomized primal optimization problem.
min p(a a a, µ µ µ) Subject to: The corresponding Lagrangian randomized dual problem is max λ λ λ min p(a a a, µ µ µ) Let us denote the optimal solutions of the primal optimization problem in Eqn. 17, randomized primal optimization problem in Eqn. 24, dual optimization problem in Eqn. 18, and randomized dual optimization problem in Eqn. 25 by P * , P R * , D * , and DR * respectively. We have the following lemmas.

Lemma 3.2:
P R * ≤ P * (26) Proof: The lemma follows from the fact that each deterministic variable can be considered as a random variable with a singleton probability distribution.
(29) Subject to: It can be check that P R * ǫ ≤ P R * , and P R * ǫ → P R * , as ǫ → 0. The dual of the P R ǫ problem DR ǫ is max λ λ λ − ,λ λ λ + min p(a a a, µ µ µ) It can be also checked that DR * ǫ → DR * , as ǫ → 0. Now we show that P R ǫ is a convex optimization problem. Let p 1 (a a a, µ µ µ), p 2 (a a a, µ µ µ) be two probability distributions satisfying all the constraints in the P R ǫ problem. Let p(a a a, µ µ µ) = αp 1 (a a a, µ µ µ) + (1 − α)p 2 (a a a, µ µ µ), where, 0 ≤ α ≤ 1. Equivalently, we may introduce a random variable z, P(z = 1) = α, P(z = 2) = 1 − α; p(a a a, µ µ µ) = p 1 (a a a, µ µ µ), if z = 1, and p(a a a, µ µ µ) = p 2 (a a a, µ µ µ), if z = 2. We can show that p(a a a, µ µ µ) satisfies the constraint in Eqn. 30 as follows. Similarly, We can also show that p(a a a, µ µ µ) satisfies the constraint in Eqn. 31 by using the fact that the expectation is a linear functional. Finally, the objective function in Eqn. 29 is also convex, because the expectation is a linear functional. Therefore, the optimization problem P R ǫ is a convex optimization problem. Because, P R ǫ is a convex optimization problem and the Slater condition holds, P R * ǫ = DR * ǫ according to the strong duality theorem [20, Thm. 6.7]. Therefore, P R * = DR * . Lemma 3.5: Assume max n,m ||x x x n −x x x m || 2 ≤ V , for a fixed upper bound V , where ||·|| 2 denotes the Euclidean norm. Then P R * → P * , as N goes to infinity.
Proof: Let p * (a a a, µ µ µ) denote the optimal solution of the randomized primal problem. We can construct a probability distribution p(a a a, µ µ µ) as follows.
p(a a a, µ µ µ) = p * ( µ µ µ) N n=1 p * (a n1 , . . . , a nJ | µ µ µ), where, the probability distributions at the right hand are marginal distributions. It can be checked that the probability p(a a a, µ µ µ) achieves the exactly same objective function and constraint function values in the randomized primal problem as the probability distribution p * (a a a, µ µ µ). Therefore, we can assume that p * (a a a, µ µ µ) takes the form in Eqn. 36 without the loss of generality.
We define the typical set T (ǫ) as The probability that (a a a, µ µ µ) is not in the typical set T (ǫ) can be upper bounded by using the Azuma inequality and the union bound as follows.
Due to the fact that the objective function is non-negative, the average achieved objective function values by (a a a, µ µ µ) in the typical set, Also by the above discussions, Therefore, we have that the average of the objective function in the typical set is bounded by There must exist one ( a a a,μ μ µ) in the typical set, such that the corresponding objective function is less than or equal to the above average. We can further modify the above a a a into a certain a a a ∈ Ω, a a a = (. . . , a ni , . . .), such that N n=1 and the corresponding objective function is raised by at most (J − 1) max{β id }V 2 ǫ. We can now set Clearly, a ni and µ ikd satisfy all the constraints in the primal problem. Therefore, where, (a) follows from the fact that µ ikd are the minimizer of the above quadratic function. The lemma then follows from the fact that P R * ≤ P * . Theorem 3.6: The duality gap P * − D * between the primal problem and dual problem goes to zero as the data sample number N goes to infinity.

IV. NUMERICAL RESULTS
In this section, we present numerical results for the proposed clustering algorithm. In Fig. 2, we depict the result of the proposed algorithm for the case of two overlapping clusters in a two dimensional space. Both the two clusters have zero mean. Their covariance matrices are as follows.
The total data sample number is 2048 and each cluster contains 1024 data samples. We assume that the data samples can be observed by two database hosts, where the first database host can only observe the 1024 data samples from the first cluster, and the second database host can only observe the 1024 data samples from the second cluster. After the clustering result is obtained, we randomly select 128 data samples from the first cluster and 128 data samples from the second cluster and plot these data samples in the figure. The data samples classified into one cluster are plotted as red circles and the data sample classified into the other cluster are plotted as blue squares. The percentage of missed classified data samples is 5.32%.
The clustering errors mainly occur at the regions where the two clusters overlap. The algorithm starts with two randomly selected unitary matrices A A A 1 , and A A A 2 . We observe that these matrices converge quickly. We also experiment with the cases that each database host observes a mixture of data samples from the two clusters with various percentages. The obtained results are not significantly different from the result in Fig. 2. In Fig. 3, we depict the result of the proposed algorithm for the case of two overlapping clusters with one cluster having a singular covariance matrix. Both the two clusters have zero mean. Their covariance matrices are as follows.

52000 52000 35600
, The total data sample number is 2048 and each cluster contains 1024 data samples. There are two database hosts, and the first database host can only observe the 1024 data samples from the first cluster, and the second database host can only observe the 1024 data samples from the second cluster. In the formulated optimization problem, a term σ 2 n I I I 2 , σ 2 n = 0.5, is added to the objective function. The clustering results of randomly selected 256 data samples are shown in the figure. The percentage of missed classified data samples is 1.71%. The results for the cases that each database host observes a mixture of data samples from the two clusters with various percentages are not significantly different from the result in the figure. The proposed clustering algorithm does not have any numerical or convergence difficulties for these cases. In Fig. 4, we depict the result of the proposed algorithm for the case of two clusters with different means. The first cluster has zero mean and covariance matrix 80000 52000 52000 35600 .
The total data sample number is 2048 and each cluster contains 1024 data samples. There are two database hosts, the first database host can only observe the 1024 data samples from the first cluster, and the second database host can only observe the 1024 data samples from the second cluster. The percentage of missed classified data samples is 2.29%. The results for the cases that each database host observes a mixture of data samples from the two clusters with various percentages are not significantly different from the result in the figure. In summary, we find that the proposed clustering algorithm has low missed classification probability and fast convergence speeds. The algorithm does not have numerical or convergence difficulties for the case of singular covariance matrices. The proposed algorithm is a promising approach for future largescale data analysis.

V. CONCLUSION
This paper proposes a large-scale data clustering algorithm based on distributed optimization. We show that the duality gap of the considered optimization problem goes to zero as the problem size goes to infinity. Therefore, the global optimization problem can be decomposed into small-scale sub optimization problems by using the Dantzig-Wolfe method. The small-scale sub optimization problems can be solved using a group of computers coordinated by one center processor. Numerical results show that the proposed algorithm is effective, efficient and does not have numerical or convergence difficulties.