A Scalable Kernel-Based Algorithm for Semi-Supervised Metric Learning

1 Introduction

1.1 Semi-Supervised Learning

In supervised learning, we are given a training sample in the form of input-output pairs. The learning task is to find a functional relationship that maps any input to an output such that disagreement with future input-output observations is minimized. Classification and regression problems are the most common supervised learning problems for discrete-valued and continuous-valued outputs, respectively. In unsupervised learning, we are given a training sample of objects with no output values, with the aim of extracting some structure from them to obtain a concise representation or to gain some understanding of the process that generated the data. Clustering, density estimation, and novelty detection problems are common unsupervised learning problems.

Over the past decade or so, there has been growing interest in exploring new learning problems between the supervised and unsupervised learning extremes. These methods are generally referred to as semi-supervised learning methods, although there exist large variations in both the problem formulation and the approach to solve the problem. The semisupervised learning literature is too large to do a comprehensive review here. Interested readers are referred to some good surveys, e.g., [Zhu, 2006].

∗This research has been supported by Competitive Earmarked Research Grant 617404 from the Research Grants Council of the Hong Kong Special Administrative Region, China.

One way to categorize many, though not all, semisupervised learning methods is to consider the type of supervisory information available for learning. Unlike unsupervised learning tasks, supervisory information is available in semi-supervised learning tasks. However, the information is in a form that is weaker than that available in typical supervised learning tasks. One type of (weak) supervisory information assumes that only part (usually a limited part) of the training data are labeled. This scenario is commonly encountered in many real-world applications. One example is the automatic classification of web pages into semantic categories. Since labeling web pages is very labor-intensive and hence costly, unlabeled web pages are far more plentiful on the web. It would be desirable if a classification algorithm can take advantage of the unlabeled data in increasing the classification accuracy. Many semi-supervised classification methods [Zhu, 2006] belong to this category.

Another type of supervisory information is even weaker in that it only assumes the existence of some pairwise constraints indicating similarity or dissimilarity relationships between training examples. In video indexing applications, for example, temporal continuity in the data can be naturally used to impose pairwise similarity constraints between successive frames in video sequences. Another example is in proteomic analysis, where protein-protein interactions can naturally be represented as pairwise constraints for the study of proteins encoded by the genes (e.g., Database of Interacting Proteins (DIP), http://dip.doe-mbi.ucla.edu/). Yet another example is the anti-spam problem. Recent studies show that more than half of the e-mail in the Internet today is spam, or unsolicited commercial e-mail. One recent approach to spam detection is based on the trustworthiness of social networks [Boykin and Roychowdhury, 2005]. Such social networks can naturally be represented as graphs with edges between nodes representing pairwise relationships in the social networks.

While supervisory information in the form of limited labeled data can be transformed into pairwise similarity and dissimilarity constraints, inverse transformation is in general not possible except for the special case of two-class problems.

In this sense, the second type of supervisory information is of a weaker form and hence the corresponding learning problem is more difficult to solve. The focus of our paper is on this category of semi-supervised learning problems.

1.2 Semi-SupervisedMetric Learning Based on Pairwise Constraints

While many semi-supervised learning methods assume the existence of limited labeled data [Zhu, 2006], there are much fewer methods that can work with pairwise constraints only. We survey the most representative methods in this subsection.

Very often, the pairwise constraints simply state whether two examples belong to the same class or different classes. [Wagstaff and Cardie, 2000] first used such pairwise information for semi-supervised clustering tasks by modifying the standard k-means clustering algorithm to take into account the pairwise similarity and dissimilarity constraints. Extensions have also been made to model-based clustering based on the expectation-maximization (EM) algorithm for Gaussian mixture models [Shental et al., 2004; Lu and Leen, 2005]. However, these methods do not explicitly learn a distance function but seek to satisfy the constraints, typically for clustering tasks. Hence, they are sometimes referred to as constraint-based clustering methods.

As a different category, some methods have been proposed to learn a Mahalanobis metric or some other distance function based on pairwise constraints. [Xing et al., 2003] formulated a convex optimization problem with inequality constraints to learn a Mahalanobis metric and demonstrated performance improvement in a subsequent clustering task. Solving a similar problem to learn a Mahalanobis metric, the relevant component analysis (RCA) algorithm [Bar-Hillel et al., 2003; 2005] was proposed as a simpler and more efficient algorithm than that of [Xing et al., 2003]. However, RCA can make use of similarity constraints only. [Hertz et al., 2004] proposed a distance function learning method called DistBoost. However, there is no guarantee that the distance function learned is a metric. [Bilenko et al., 2004] explored the possibility of integrating constraint-based clustering and semi-supervised clustering based on distance function learning. The distance function learning methods reviewed above either learn a Mahalanobis metric that corresponds to linear transformation only, or learn a distance function that is not a metric. In our previous work [Chang and Yeung, 2004], we proposed a metric learning method that corresponds to nonlinear transformation. While this method is more powerful than the linear methods, the optimization problem is non-convex and hence is more complicated to solve.

In this paper, we focus on the distance function learning approach because distance functions are central to many learning models and algorithms. This also makes it easier to achieve out-of-sample generalization. Moreover, we focus on learning metrics because this allows us to formulate the metric learning problem based on the kernel approach [Sch¨olkopf and Smola, 2002], which provides a disciplined, computationally appealing approach to nonlinear metric learning.

1.3 Organization of Paper

In Section 2, we propose a simple and efficient kernel-based metric learning method based on pairwise similarity constraints. Like many kernel methods, a limitation of this method is that it does not scale up well with the sample size. In Section 3, we address the scalability issue by applying low-rank approximation to the kernel matrix. This extended method can also naturally give rise to out-of-sample generalization, which is addressed in Section 4. We present some experimental results in Section 5 to demonstrate the effectiveness of our metric learning algorithm. Finally, Section 6 concludes the paper.

2 Kernel-BasedMetric Learning

2.1 Problem Setup

Let {xi}ni=1 be a set of n data points in some input space X. Suppose we have a Mercer kernel ˆ k which induces a nonlinear featuremap ˆ φ from X to some reproducing kernel Hilbert space H [Sch¨olkopf and Smola, 2002]. The corresponding set of feature vectors in H is { ˆ φ(xi)}ni=1 and the kernel matrix is ˆK = [ˆ k(xi, xj)]n×n = [〈ˆ φ(xi), ˆ φ(xj)〉]n×n. Choices for ˆ k include the Gaussian RBF kernel and the polynomial kernel. We apply a centering transform such that the feature vectors in H have zero mean. The resulting kernel matrix K can be computed as K = [k(xi, xj )]n×n = [〈φ(xi), φ(xj )〉]n×n = H ˆKH, where the centering matrix H = I - (1/n)11T with I being the n×n identity matrix and 1 being the n×1 vector of ones.

We consider a type of semi-supervised metric learning in which the supervisory information is given in the form of pairwise similarity constraints. 1 Specifically, we are given a set of point pairs, which is a subset of X × X, as S = {(xi, xj) | xi and xj belong to the same class}. Our goal is to make use of S to learn a better metric through modifying the kernel so that the performance of some subsequent task (e.g., clustering, classification) based on the metric is improved after kernel learning.

2.2 Kernel Learning

Since the kernel matrix K is symmetric and positive semidefinite, we can express it as K = ∑pr=1 λrvrvT r = ∑pr=1 λrKr, where λ1 ≥ · · · ≥ λp > 0 are the p ≤ n positive eigenvalues of K, v1, . . . , vp are the corresponding normalized eigenvectors, and Kr = vrvT r .

We consider a restricted form of kernel matrix learning by modifying K through changing the λr's while keeping all Kr's fixed. To ensure that the eigenvalues are nonnegative, we rewrite K as Kβ = [kβ(xi, xj)]n×n = [〈φβ(xi), φβ(xj)〉]n×n = ∑pr=1 β2 rKr, which represents a family of kernel matrices parameterized by β = (β1, . . . , βp)T .

1 As we will see later in the formulation of the optimization problem for kernel learning, these constraints are "soft" constraints rather than "hard" constraints in that they are only preferred, not enforced. This makes it easy to handle noisy constraints, i.e., erroneous supervisory information, if we so wish.

We perform kernel learning such that the mean squared Euclidean distance induced by Kβ between feature vectors in H corresponding to point pairs in S is reduced. Thus the criterion function for optimization is

where bi is the ith column of the p × p identity matrix 2 and DS is a p × p diagonal matrix with diagonal entries

To prevent β from degenerating to the zero vector 0 and to eliminate the scaling factor, we minimize the convex function JS(β) subject to the linear constraint 1T β = c for some constant c > 0. This is a simple convex optimization problem with a quadratic objective function and a linear equality constraint. We introduce a Lagrange multiplier ρ to minimize the following Lagrangian:

The optimization problem can be solved easily to give the optimal value of β as the following closed-form solution:

Note that D-1 S exists as long as all the diagonal entries of DS are positive, which is usually true. 3 We set the constant c = ∑pr=1 √λr .

3 Scalable Kernel Learning

The kernel learning method described above requires performing eigendecomposition on K. In case n is very large and hence K is a large matrix, operations such as eigendecomposition on K are computationally demanding. In this section, we apply low-rank approximation to extend the kernel learning method so that it scales up well with n.

2 This indicator variable will be "overloaded" later to refer to a column of any identity matrix whose size is clear from the context.

3 In case DS is really singular (though a rare case), a common way to make it invertible is to add a term ϵI to DS where ϵ is a small positive constant.

3.1 Low-Rank Approximation

We apply low-rank approximation to approximate K by another n × n symmetric and positive semi-definite matrix ˜K:

where W ∈ Rn×m is an n × m matrix and L ∈ Rm×m is an m × m symmetric and positive semi-definite matrix, with m ≪ n.

There are different ways to construct L for low-rank approximation. We consider one way which constructs L using a subset of the n data points. We refer to these points as landmarks [de Silva and Tenenbaum, 2003; Weinberger et al., 2005; Silva et al., 2006]. Without loss of generality, we assume that the n points are ordered in such a way that the first m points form the set of landmarks {xi}mi=1 . We should ensure that all points involved in S are chosen as landmarks. Other landmarks are randomly sampled from all the data points. Similar to K in the previous section, L is obtained here by applying the centering transform to ˆL, as L = HˆLH, where ˆL = [ˆ k(xi, xj)]m×m = [〈ˆ φ(xi), ˆ φ(xj)〉]m×m is the upper-left m×m submatrix of ˆK.

We apply eigendecomposition on L and express it as

where μ1 ≥ · · · ≥ μq > 0 are the q ≤ m positive eigenvalues of L, α1, . . . , αq are the corresponding normalized eigenvectors, Dμ = diag(μ1, . . . , μq), and Vα = [α1, . . . , αq]. Substituting (6) into (5), we can rewrite ˜K as ˜K = ∑qr=1 μr ˜Kr , where ˜Kr = (Wαr)(Wαr)T .

3.2 Kernel Learning

We apply low-rank approximation to devise a scalable kernel learning algorithm which can be seen as an extension of the algorithm described in Section 2. We use ˜K to approximate K and define the following parameterized family of kernel matrices: ˜Kβ = ∑qr=1 β2 r ˜Kr. Note, however, that β is now a q × 1 vector rather than a p × 1 vector.

The optimal value of β has the same form as (4), except that the constant c is set to ∑qr=1 √μr andDS is now a q×q diagonal matrix with diagonal entries

3.3 Computing the EmbeddingWeights

A question that remains to be answered is how to obtainW for low-rank approximation. We use a method that is similar to locally linear embedding (LLE) [Roweis and Saul, 2000; Saul and Roweis, 2003], with two differences. First, we only use the first part of LLE to obtain the weights for locally linear fitting. Second, we perform locally linear fitting in the kernelinduced feature space H rather than the input space X .

LetW= [wij]n×m and wi = (wi1, . . . ,wim)T . If xi is a landmark, i.e., 1 ≤ i ≤ m, then

If xi is not a landmark, then we minimize the following function to obtain wi:

where Ni is the set of K nearest landmarks of φ(xi) in H, subject to the constraints ∑φ(xj)∈Ni wij = 1T wi = 1 and wij = 0 for all φ(xj) ̸∈ Ni. We can rewrite E(wi) as

where

is the local Grammatrix of xi in H.

To prevent wi from degenerating to 0, we minimize E(wi) subject to the constraints ∑φ(xj)∈Ni wij = 1T wi = 1 and wij = 0 for all φ(xj) ̸∈ Ni. As above for the kernel learning problem, we solve a convex optimization problem with a quadratic objective function and a linear equality constraint. The Lagrangian with a Lagrange multiplier α is as follows:

The closed-form solution for this optimization problem is given by wi = (G-1 i 1)/(1T G-1 i 1) if G-1 i exists. 4

Instead of performing matrix inversion, a more efficient way of finding the solution is to solve the linear system of equations Gi ˆwi = 1 for ˆwi and then compute wi as wi = ˆwi/(1T ˆwi) to ensure that the equality constraint 1T wi = 1 is satisfied.

We assume above that the neighborhood relationships between points in H and the local Gram matrices remain fixed during the kernel learning process. A simple extension of this basic algorithm is to repeat the above procedure iteratively using the learned kernel at each iteration. Thus the basic algorithm is just a special case of this iterative extension when the number of iterations is equal to one.

4 Out-of-Sample Generalization

The exact form of out-of-sample generalization depends on the operation we want to perform. For example, the n given points are first clustered into C ≥ 2 classes after kernel learning, and then a new data point x is classified into one of the C classes. We are interested in the case where both the clustering of the n points and the classification of new points are based on the same Euclidean metric in H.

The key idea of our out-of-sample generalization scheme rests on the observation that kernel principal component analysis (KPCA) [Sch¨olkopf et al., 1998] can be performed on {xi}ni=1 to obtain an embedding in a q-dimensional subspace Y of H, so that the non-landmark points {xi}ni=m+1 and any out-of-sample point x can be embedded into Y in the same way.

4 Similar to DS above, we may add ϵI to make sure that Gi is invertible.

Let {u1, . . . , uq} be an orthonormal basis with each ur ∈ H being a unit vector along the direction of the rth principal component. We define U = [u1, . . . , uq]. Then each landmark xi (i = 1, . . . ,m) can be embedded into Y to give a q-dimensional vector yi = UT φ(xi). Since ur = 1 √μr ∑mj=1 αjrφ(xj), we can express yi as yi = D-1/2 μ VTαLbi = D1/2 μ VTαbi. Let Ym = [y1, . . . , ym]. So

For the non-landmark points {xi}ni=m+1 , from (9), we use ˜ φ(xi) = ∑φ(xj)∈Ni wijφ(xj) to approximate φ(xi) and ˜yi to denote the embedding of ˜ φ(xi) in Y. Thus we have

(13) Similarly, for any out-of-sample example x, we use ˜ φ(x) = ∑φ(xj)∈Ni wjφ(xj) to approximate φ(x) and ˜y to denote the embedding of ˜ φ(x) in Y. The embedding weights w = (w1, . . . ,wm)T can be determined as in Section 3.3. Similar to the non-landmark points {xi}ni=m+1 , we can obtain

Based on (13) and (14), the squared Euclidean distance between ˜yi and ˜y in Y before kernel learning can be expressed

The squared Euclidean distance after kernel learning is

5 Experimental Results

In this section, we present some experiments we have performed based on both artificial and real-world data.

5.1 Experimental Setup

We compare two versions of our kernel-based metric learning method described in Sections 2 and 3 with RCA [Bar-Hillel et al., 2003; 2005], which is a promising linear metric learning method that usually performs equally well as other computationally more demanding methods. For baseline comparison, we also include two metrics without metric learning. They are the Euclidean metric in the input space and the Euclidean metric in the feature space induced by a Gaussian RBF kernel.

For each data set, we randomly generate 10 different S sets. For small data sets, we can learn Kβ without low-rank approximation. For large data sets, in addition to the data points involved in S, we also randomly select some other points as landmarks for learning ˜Kβ. We use the iterative extension of the scalable kernel learning algorithm with the number of iterations equal to 3. We also measure the change in metric learning performance as the number of landmarks increases.

5.2 An Illustrative Example

Figure 1 uses the XOR data set to compare the performance of different metrics, some with metric learning. Figure 1(a) shows 360 data points in the input space. The points with the same color and mark belong to the same class. The randomly generated point pairs corresponding to the similarity constraints in S are shown in solid lines. Figure 1(b) shows 40 new data points in the input space. The results obtained by RCA, RBF kernel and two versions of our method are shown in Figure 1(c)-(f). RCA performs metric learning directly in the input space, which can be generalized to new data using the learned linear transformation. For the kernel methods, we apply KPCA using the (learned) kernel matrix to embed the data points to a 2-dimensional space, as shown in Figure 1(d)-(f). New data points can be embedded using the generalization method described in Section 4. As expected, RCA does not perform satisfactorily for the XOR data set since it can only perform linear transformation. On the other hand, our kernel-based metric learning method can group the points according to their class membership. The result of the scalable kernel learning method with low-rank approximation based on 40 landmarks is almost as good as that of the basic algorithm based on all 400 data points. Moreover, the result on embedding of new data points verifies the effectiveness of our out-of-sample generalization method.

5.3 Quantitative Performance Comparison

Let {yi}ni=1 be the set of true class labels of the n data points. We define the following performance measure:

where db = (1/nb) ∑yi̸=yj ‖xi - xj‖ is the mean betweenclass distance with nb being the number of point pairs with different class labels, and dw = (1/nw) ∑yi=yj ‖xi -xj‖ is the mean within-class distance with nw being the number of point pairs with the same class label. Note that J is closely related to the optimization criterion JS in (1), except that JS is defined for the labeled data (i.e., data involved in the pairwise constraints) only while J is for all data assuming the existence of true class labels. For kernel methods, we use φ(xi) or ˜ φ(xi) in place of xi and apply the kernel trick to compute the mean distances and hence J. A larger value of J corresponds to a better metric due to its higher class separability.

We first perform some experiments on a much larger XOR data set with 8,000 data points. We randomly select 50 similarity constraints to form S and measure the metric learning performance in terms of the J value for an increasing number of landmarks. Table 1 shows the results for different metrics. For the metric learning methods (i.e., RCA and our method), we show for each trial the mean (upper) and standard deviation (lower) over 10 random runs corresponding to different S sets. From the results, we can see that our method outperforms the other methods significantly. Moreover, using more landmarks generally gives better results.

We further perform some experiments on real-world data sets. One of them is the Isolet data set from the UCI Machine Learning Repository which contains 7,797 isolated spoken English letters belonging to 26 classes, with each letter represented as a 617-dimensional vector. Other data sets are handwritten digits from the MNIST database. 5 The digits in the database have been size-normalized and centered to 28×28

5 http://yann.lecun.com/exdb/mnist/

gray-level images. Hence the dimensionality of the input space is 784. In our experiments, we randomly choose 2,000 images for each digit from a total of 60,000 digit images in the MNIST training set. We use 50 similarity constraints and 100 landmarks in the experiments. The results for Isolet and different digit data sets are shown in Table 2. From the results, we can again see that the metric learned by our method is the best in terms of the J measure.

6 Concluding Remarks

We have presented a simple and efficient kernel-based semisupervised metric learning method based on supervisory information in the form of pairwise similarity constraints. Not only does it scale up well with the data set size, it can also naturally lead to out-of-sample generalization. Although previous studies by other researchers showed that pairwise dissimilarity constraints usually cannot help much in many realworld applications, there are situations when incorporating them may still be helpful and hence we plan to extend our method to incorporate dissimilarity constraints as well. In our low-rank approximation scheme, besides including those points involved in S, we also randomly sample some other points as landmarks. A recent study [Silva et al., 2006] shows that non-uniform sampling of landmarks for manifold learning can give parsimonious approximations using only very few landmarks. We will pursue research along this direction in our future work.