Local Learning Regularized Nonnegative Matrix Factorization

1 Introduction

Nonnegative Matrix Factorization (NMF) [Lee and Seung, 2000] has been widely used in machine learning and data mining. It aims to find two nonnegative matrices whose product can well approximate the nonnegative data matrix, which naturally lead to parts-based representation. Recent years, many variants of NMF have been proposed. [Li et al., 2001] proposed a local NMF (LNMF) which imposes a spatially localized constraint on the bases. [Hoyer, 2004] proposed a NMF with sparseness constraint. [Ding et al., 2008] proposed a semi-NMF approach which relaxes the nonnegative constraint on the base matrix. [Ding et al., 2006] proposed a nonnegative matrix tri-factorization for co-clustering. All the methods mentioned above are unsupervised, while [Wang et al., 2004] and [Zafeiriou et al., 2006] proposed independently a discriminative NMF (DNMF), which adds an additional constraint seeking to maximize the between-class scatter and minimize the within-class scatter in the subspace spanned by the bases.

Recent studies show that many real world data are actually sampled from a nonlinear low dimensional manifold which is embedded in the high dimensional ambient space [Roweis and Saul, 2000] [Niyogi, 2003]. Yet NMF does not exploit the geometric structure of the data. In other word, it assumes that the data points are sampled from a Euclidean space. This greatly limits the application of NMF for the data lying on manifold. To address this problem, [Cai et al., 2008] proposed a graph regularized NMF (GNMF), which assumes that the nearby data points are likely to be in the same cluster, i.e. cluster assumption [Chapelle et al., 2006].

In this paper, we present a novel nonnegative matrix factorization method. It is based on the assumption that the cluster label of each point can be predicted by the data points in its neighborhood, i.e. local learning assumption, which is the philosophy of local learning algorithm [Bottou and Vapnik, 1992]. This assumption is embodied by a local learning regularization, which exploits both the discriminative information and the geometric structure. We constrain the NMF with local learning regularization, resulting in a local learning regularized NMF (LLNMF). LLNMF not only inherits the advantages of NMF, e.g. nonnegativity, but also overcomes the shortcomings of NMF, i.e. Euclidean assumption based and does not take into account discriminative information. We will show that it can be optimized via an iterative multiplicative updating algorithm and its convergence is theoretically guaranteed. Experiments on many benchmark data sets demonstrate that the proposed method outperforms NMF and its variants as well as many other state of the art clustering algorithms.

The remainder of this paper is organized as follows. In Section 2 we will briefly review NMF. In Section 3, we present LLNMF, followed with its optimization algorithm along with the proof of the convergence of the proposed algorithm. Experiments on many benchmark data sets are demonstrated in Section 4. Finally, we draw a conclusion in Section 5.

2 AReview of NMF

In this section, we will briefly review NMF [Lee and Seung, 2000]. Given a nonnegative data matrix X = [x1, . . . , xn] ∈ Rd×n + , each column of X is a data point. NMF aims to find two nonnegative matrices U ∈ Rd×m + and V ∈ Rm×n + which minimize the following objective

where || · ||F is Frobenius norm. To optimize the objective, [Lee and Seung, 2000] presented an iterative multiplicative updating algorithm as follows

In the clustering setting of NMF [Xu et al., 2003], V ∈ Rm×n + is the cluster assignment matrix where m is the prescribed number of clusters. And the cluster label yi of data point xi can be calculated as

In the rest of this paper, we denote V = [v1 , . . . , vm]T where vl = [vl1 , . . . , vln ]T ∈ Rn+, 1 ≤ l ≤ m is the row vector of V, and vl i is the l-th cluster assignment of xi.

3 The ProposedMethod

In this section, we first introduce local learning regularization. Then we will present Local learning Regularized Nonnegative Matrix Factorization, followed with its optimization algorithm. The convergence of the proposed algorithm is also proved.

3.1 Local Learning Regularization

According to [Bottou and Vapnik, 1992], selecting a good predictor f in a global way might not be a good strategy because the function set f(x) may not contain a good predictor for the entire input space. However, it is much easier for the set to contain some functions that are capable of producing good predictions on some specified regions of the input space. Therefore, if we split the whole input space into many local regions, then it is usually more effective to minimize predicting cost for each region. This inspires the work [Wu and Sch¨olkopf, 2006], which adopted supervised learning idea for unsupervised learning problem. In the following, we will introduce how to construct the local predictors, and derive a local learning regularization.

For each data point xi, we use N(xi) to denote its neighborhood. We further construct predicting function f l i (x), 1 ≤ l ≤ m in N(xi), to predict the cluster label of {xj}xj∈N(xi). Note that for f l i (x), the superscript l indicates that it is for the l-th cluster, while the subscript i means that it is trained within the neighborhood of xi.

To fit the predictor f l i (x), we use regularized ridge regression [Hastie et al., 2001], which minimize the following loss function

where ni is the cardinality of N(xi), vlj is the l-th cluster assignment of xj, ||f l i ||I measures the smoothness of f l i with respect to the intrinsic data manifold, λi > 0 is the regularization parameter. In this paper, we assume λ1 = . . . = λn = λ and n1 = n2 = . . . = nn = k, for simplicity.

According to Representor Theorem [Sch¨olkopf and Smola, 2002], we have

where K : X × X → R is a positive definite kernel function, and βl ij are the expansion coefficients. Therefore we have

where f l i ∈ Rni denotes the vector [f l i (xj )]T , xj ∈ N(xi), Ki ∈ Rni×ni is the kernel matrix defined on the neighborhood of xi, i.e. N(xi), and βli = [βi1, . . . , βini ]T ∈ Rni×1 is the expansion coefficient vector. Bringing Eq.(5) back into Eq.(4), we can derive the following loss function

where vli ∈ Rni denotes the vector [vlj ]T , xj ∈ N(xi). By setting ∂J l i ∂βli = 0, we can get that

where I ∈ Rni×ni is the identity matrix. Substituting Eq.(8) into Eq.(5), we have

where ki ∈ Rni denotes the vector [k(xi, xj)]T , xj ∈ N(xi), and αi = [αi1, . . . , αini ]T ∈ Rni×1 .

After the local predictors are constructed, we will combine them together by minimizing the sum of their prediction errors

Substitute Eq.(9) into Eq.(10), we obtain

where vl = [vl1 , . . . , vln ]T ∈ Rn is the l-th row of V, I ∈ Rn×n is identity matrix, L = (G-I)(G-I) and G ∈ Rn×n is defined as follows

Eq.(11) is called as Local Learning Regularization. The better the cluster label of each point is predicted by the data points in its neighborhood, the smaller the local learning regularizer will be.

3.2 Local Learning Regularized NMF

Our assumption is that the cluster label of each point can be predicted by the data points in its neighborhood. To apply this idea for NMF, we constrain NMF in Eq.(1) with local learning regularization in Eq.(11) as follows

where μ is a positive regularization parameter controlling the contribution of the additional constraint. We call Eq.(13) as Local Learning Regularized Nonnegative Matrix Factorization (LLNMF). Letting μ = 0, Eq.(13) degenerates to the original NMF. To make the objective in Eq.(13) lower bounded, we use L2 normalization on rows of V in the optimization, and compensate the norms of V to U.

In the following, we will give the solution to Eq.(13).

Since U ≥ 0, V ≥ 0, we introduce the Lagrangian multiplier γ ∈ Rd×m and η ∈ Rm×n , thus, the Lagrangian function is

Setting ∂L(U,V) ∂U = 0 and ∂L(U,V) ∂V = 0, we obtain

Using the Karush-Kuhn-Tucker condition [Boyd and Vandenberghe, 2004] γijUij = 0 and ηijVij = 0, we get

Introduce

where L+ij = (|Lij | + Lij )/2 and L-ij = (|Lij | - Lij )/2. Substitute Eq.(17) into Eq.(16), we obtain

Eq.(18) leads to the following updating formula

3.3 Convergence Analysis

In this section, we will investigate the convergence of the updating formula in Eq.(19). We use the auxiliary function approach [Lee and Seung, 2000] to prove the convergence. Here we first introduce the definition of auxiliary function [Lee and Seung, 2000].

Definition 3.1. [Lee and Seung, 2000] Z (h, h′ ) is an auxiliary function for F (h) if the conditions

are satisfied.

Lemma 3.2. [Lee and Seung, 2000] If Z is an auxiliary function for F, then F is non-increasing under the update

Lemma 3.3. [Ding et al., 2008] For any nonnegative matrices A ∈ Rn×n , B ∈ Rk×k , S ∈ Rn×k ,S′ ∈ Rn×k , and A, B are symmetric, then the following inequality holds

Theorem 3.4. Let

Then the following function

is an auxiliary function for J(U). Furthermore, it is a convex function in U and its global minimum is

Proof. See Appendix A

Theorem 3.5. UpdatingUusing Eq.(19) will monotonically decrease the value of the objective in Eq.(13), hence it converges.

Proof. By Lemma 3.2 and Theorem 3.4, we can get that J (U0 ) = Z (U0 , U0 ) ≥ Z (U1 , U0 ) ≥ J(U1 ) ≥ . . . So J(U) is monotonically decreasing. Since J(U) is obviously bounded below, we prove this theorem.

Theorem 3.6. Let

Then the following function

is an auxiliary function for J(V). Furthermore, it is a convex function in V and its global minimum is

Theorem 3.7. Updating V using Eq.(19) will monotonically decrease the value of the objective in Eq.(13), hence it converges.

Proof. By Lemma 3.2 and Theorem 3.6, we can get that J (V0 ) = Z (V0 , V0 ) ≥ Z (V1 , V0 ) ≥ J(V1 ) ≥ . . . So J(V) is monotonically decreasing. Since J(V) is obviously bounded below, we prove this theorem.

4 Experiments

In this section, we evaluate the performance of the proposed method. We compare our method with Kmeans, Normalized Cut (NC) [Shi and Malik, 2000], Local Learning Clustering (LLC) [Wu and Sch¨olkopf, 2006], NMF [Lee and Seung, 2000], Semi-NMF [Ding et al., 2008], Orthogonal Nonnegative Matrix Tri-Factorization (ONMTF) [Ding et al., 2006] and GNMF [Cai et al., 2008].

4.1 EvaluationMetrics

To evaluate the clustering results, we adopt the performance measures used in [Cai et al., 2008]. These performance measures are the standard measures widely used for clustering.

Clustering Accuracy Clustering Accuracy discovers the one-to-one relationship between clusters and classes and measures the extent to which each cluster contained data points from the corresponding class. Clustering Accuracy is defined as follows:

where ri denotes the cluster label of xi, and li denotes the true class label, n is the total number of documents, δ(x, y) is the delta function that equals one if x = y and equals zero otherwise, and map(ri) is the permutation mapping function that maps each cluster label ri to the equivalent label from the data set.

Normalized Mutual Information The second measure is the Normalized Mutual Information (NMI), which is used for determining the quality of clusters. Given a clustering result, the NMI is estimated by

where nk denotes the number of data contained in the cluster Dk(1 ≤ k ≤ C), ˆnm is the number of data belonging to the Lm(1 ≤ m ≤ C), and nk,m denotes the number of data that are in the intersection between the cluster Dk and the class Lm. The larger the NMI is, the better the clustering result will be.

4.2 Data Sets

In our experiment, we use 6 data sets which are widely used as benchmark data sets in clustering literature [Cai et al., 2008] [Ding et al., 2006].

Coil20 1 This data set contains 32 × 32 gray scale images of 20 3D objects viewed from varying angles. For each object there are 72 images.

PIE The CMU PIE face database [Sim et al., 2003] contains 68 individuals with 41368 face images as a whole. The face images were captured by 13 synchronized cameras and 21 flashes, under varying pose, illumination and expression. All the images were also resized to 32 × 32.

CSTR This is the data set of the abstracts of technical reports published in the Department of Computer Science at a university. The data set contained 476 abstracts, which were divided into four research areas: Natural Language Processing (NLP), Robotics/Vision, Systems and Theory.

Newsgroup4 The Newsgroup4 data set used in our experiments is selected from the famous 20-newsgroups data set 2 . The topic rec containing autos, motorcycles, baseball and hockey was selected from the version 20news-18828. The Newsgroup4 data set contains 3970 documents.

WebKB4 The WebKB dataset contains webpages gathered from university computer science departments. There are about 8280 documents and they are divided into 7 categories: student, faculty, staff, course, project, department and other, among which student, faculty, course and project are four most populous entity-representing categories.

WebACE The data set contains 2340 documents consisting of news articles from Reuters new service via the Web in October 1997. These documents are divided into 20 classes.

Table.1 summarizes the characteristics of the real world data sets used in this experiment.

4.3 Parameter Settings

Since many clustering algorithms have one or more parameters to be tuned, under each parameter setting, we repeat clustering 20 times, and the mean result is computed. We report the best mean result for each method to compare with each other. We set the number of clusters equal to the true number of classes for all the data sets and clustering algorithms.

For NC [Shi and Malik, 2000], the scale parameter of Gaussian kernel for constructing adjacency matrix is set by the grid {10-3 , 10-2 , 10-1 , 1, 10, 102 , 103 }.

1 http://www1.cs.columbia.edu/CAVE/software/softlib/coil-

20.php

2 http://people.csail.mit.edu/jrennie/20Newsgroups/

For LLC, the neighborhood size for computing the local learning regularization is determined by the grid {5, 10, 20, 30, 40, 50, 80} according to [Wu and Sch¨olkopf, 2006]. For the image datasets (e.g. Coil20, PIE), Gaussian kernel is used with the scale parameter tuned the same as in NC, while for the text datasets (e.g. CSTR, Newsgroup4, WebKB4, WebACE), the cosine kernel is adopted.

For ONMTF, the number of word clusters is set to be the same as the number of document clusters, i.e. the true number of classes in our experiment, according to [Ding et al., 2006].

For GNMF, the neighborhood size to construct the graph is set by search the grid {1, 2, 3, . . . , 10} according to [Cai et al., 2008], and the regularization parameter is set by the grid {0.1, 1, 10, 100, 500, 1000}.

For the proposed algorithm, the neighborhood size k for computing the local learning regularization is determined by the grid {5, 10, 20, 30, 40, 50, 80}, and the regularization parameter μ is set by search the grid {0.1, 1, 10, 100, 500, 1000}. Kernel selection is the same as that in LLC.

Note that no parameter selection is needed for Kmeans, NMF and Semi-NMF, given the number of clusters.

4.4 Clustering Results

The clustering results are shown in Table 2 and Table 3. Table 2 shows the clustering accuracy of all the algorithms on all the data sets, while Table 3 shows the normalized mutual information.

We can see that our method outperforms the other clustering methods on all the data sets. The superiority of our method may arise in the following two aspects: (1) the local learning assumption, which is usually more powerful than cluster assumption [Chapelle et al., 2006] [Cai et al., 2008] for clustering data on manifold. (2) the nonnegativity, inheriting from NMF, which is suitable for nonnegative data, e.g. image data and text data.

5 Conclusion

In this paper, we present a local learning regularized nonnegative matrix factorization (LLNMF) for clustering, which considers both the discriminative information and the geometric structure. The convergence of the proposed algorithm is proved theoretically. Experiments on many benchmark data sets demonstrate that the proposed method outperforms NMF as well as many state of the art clustering methods.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (No.60721003, No.60673106 and No.60573062) and the Specialized Research Fund for the Doctoral Program of Higher Education. We thank the anonymous reviewers for their helpful comments.

A Proof of Theorem 3.4

Proof. We rewrite Eq.(20) as

By applying Lemma 3.3, we have

To obtain the lower bound for the remaining terms, we use the inequality that

Then

By summing over all the bounds, we can get Z(U,U′ ), which obviously satisfies (1) Z(U,U′ ) ≥ JLLNMF(U); (2)Z(U, U) = JLLNMF (U)

To find the minimum of Z(U,U′ ), we take the Hessian matrix of Z (U, U′ )

which is a diagonal matrix with positive diagonal elements. So Z (U, U′ ) is a convex function of U, and we can obtain the global minimum of Z (U, U′ ) by setting ∂Z(U,U′ ) ∂Uij = 0 and solving for U, from which we can get Eq.(21).

B Proof of Theorem 3.6

Proof. We rewrite Eq.(22) as

By applying Lemma 3.3, we have

By the inequality in Eq.(27), we have

By summing over all the bounds, we can get Z(V,V′ ), which obviously satisfies (1) Z(V,V′ ) ≥ JLLNMF(V); (2)Z(V, V) = JLLNMF (V)

To find the minimum of Z(V,V′ ), we take the Hessian matrix of Z (V, V′ )

which is a diagonal matrix with positive diagonal elements. So Z (V, V′ ) is a convex function of V, and we can obtain the global minimum of Z (V, V′ ) by setting ∂Z(V,V′ ) ∂Vij = 0 and solving for V, from which we can get Eq.(23).