Neighborhood MinMax Projections

1 Introduction

Linear dimensionality reduction is an important method when facing with high-dimensional data. Many algorithms have been proposed during the past years. Among these algorithms, Principal Component Analysis (PCA) [Jolliffe, 2002] and Linear Discriminant Analysis (LDA) [Fukunaga, 1990] are two of the most widely used methods. PCA is an unsupervised method, which does not take the class information into account. LDA is one of the most popular supervised dimensionality reduction techniques for classification. However, there exist several drawbacks in it. One drawback is that it often suffers from the Small Sample Size problem when dealing with high dimensional data. In this case, the withinclass scatter matrix Sw may become singular, which makes LDA difficult to be performed. Many approaches have been proposed to address this problem [Belhumeur et al., 1997;

∗This paper is Funded by Basic Research Foundation of Tsinghua National Laboratory for Information Science and Technology (TNList).

Chen et al., 2000; Yu and Yang, 2001]. However, these variants of LDA discard a subspace and thus some important discriminative information may be lost. Another drawback in LDA is its distribution assumption. LDA is optimal in the case that the data distribution of each class is Gaussian, which can not always be satisfied in real world applications. When the class distribution is more complex than Gaussian, LDA may fail to find the optimal discriminative directions. Moreover, the number of available projection directions in LDA is smaller than the class number [Duda. et al., 2000], but it may be insufficient for many complex problems, especially when the number of class is small.

For the distance metric based classification methods, such as the nearest neighbor classifier, learning an appropriate distance metric plays a vital role. Recently, a number of methods have been proposed to learn a Mahalanobis distance metric [Xing et al., 2003; Goldberger et al., 2005; Weinberger et al., 2006]. Linear dimensionality reduction can be viewed as a special case of learning a Mahalanobis distance metric(see section 5). This viewpoint can give a reasonable interpretation for the fact that the performance of nearest neighbor classifier can always be improved after performing linear dimensionality reduction.

In this paper, we propose a new supervised linear dimensionality reduction method, Neighborhood MinMax Projections (NMMP). The method is largely inspired by the classical supervised linear dimensionality reduction method, i.e., LDA, and the recent proposed distance metric learning method, large margin nearest neighbor (LMNN) classification [Weinberger et al., 2006]. In our method, we focus only on the pairwise points where the two points are neighbors of each other. After the transformation, we try to pull the considered pairwise points within the same class as close as possible, and take those between different classes apart. This goal can be achieved by formulating the task as a constrained optimization problem, in which the global optimum can be effectively and efficiently obtained. Compared with LDA, our method avoids the three drawbacks in LDA discussed in above. Compared with the LMNN method, our method is computationally much more efficient. The performance on several data sets demonstrates the effectiveness of our method.

2 Problem Formulation

Given the data matrix X = [x1, x2, ..., xn], xi ∈ Rd , our goal is to learn a linear transformationW : Rd → Rm, where W ∈ Rd×m and WTW = I. I is m×m identity matrix. Then the original high-dimensional data x is transformed into a low-dimensional vector:y =WT x (1)

Let each data point of class i have two kinds of neighborhood: within-class neighborhood Nw(i) and between-class neighborhood Nb(i), where Nw(i) is the set of the data's kw(i) nearest neighbors in the same class i and Nb(i) is the set of the data's kb(i) nearest neighbors in the class other than i. Obviously, 1 ≤ kw(i) ≤ ni - 1, and 1 ≤ kb(i) ≤ n - ni, where ni is the data number of class i.

Here, we focus only on the pairwise points where the two points are neighbors of each other. After the transformation, we hope that the distance of the considered pairwise points within the same class will be minimized, while the distance of those between different classes will be maximized. (see Figure 1).

After the transformationW, the sum of the Euclidean distances of the pairwise points within the same class can be formulated as:

where tr(·) denotes the trace operator of matrix, and

Here, Ci denote the class label of xi, and Cj denote the class label of xj. Obviously, Ci = Cj, and ˜Sw is positive semi-definite.

Similarly, the sum of the Euclidean distances of the pairwise points between different classes is:

where

Here, Ci ̸= Cj, and ˜Sb is positive semi-definite too.

To achieve our goal, we should maximize sb while minimize sw. The following two function can be used as objective:

As it is difficult to determine a suitable weight λ for the former objective function, we select the latter as our objective to optimization. In fact, as we will see, the latter is just a special case of the former, where the weight λ is automatically determined.

Therefore, we formulate the problem as a constrained optimization problem:

Fortunately, the globally optimal solution of this problem can be efficiently calculated. In the next section, we will describe the details for solving this constrained optimization problem.

3 The Constrained Optimization Problem

We address the above optimization problem in a more general form which is described as follows:

The constrained optimization problem: Given the real symmetric matrix A ∈ Rd×d and the positive semi-definite matrix B ∈ Rd×d, rank(B) = r ≤ d. Find a matrix W ∈ Rd×m that maximize the following objective function with the constraint of WT W = I:

At first, we propose Lemma 1, which shows that when WT W = I and m > d - r, the value of tr(WT BW) will not be equal to zero.

Lemma1. Suppose W ∈ Rd×m, WT W = I, B ∈ Rd×d is a positive semi-definite matrix, and rank(B) = r ≤ d, m > d- r, then it holds that tr(WT BW) > 0.

proof. According to the result of Rayleigh quotient [Golub and van Loan, 1996], min WT W=I tr(WT BW) = ∑mi=1 βi, where β1, β2, . . . , βm are the first m smallest eigenvalues of B. As B is positive semi-definite, rank(B) = r, and m> d-r, then ∑mi=1 βi > 0. Therefore, with the constraint of WT W = I, tr(WT BW) ≥min tr(WT BW) > 0. Thus we discuss this optimization problem in two cases.

Lemma 1 ensures that the optimal value is finite in this case. Suppose the optimal value is λ∗ , Guo [2003] has derived that max WT W=I tr(WT (A - λ∗B)W) = 0.

Note that tr(WT BW) > 0, so we can easy to see, max WT W=I tr(WT (A - λB)W) < 0 ⇒ λ > λ∗ , and max WT W=I tr(WT (A - λB)W) > 0 ⇒ λ < λ∗ .

On the other hand, max

WT W=I

tr(WT

(A - λB)W) =

γ,

where γ is the sum of the first m largest eigenvalues of AλB. Given a value λ, if γ = 0,then λ is just the optimal value, otherwise γ > 0 implies λ is smaller than the optimal value and vice versa. Thus the global optimal value of the problem can be obtained by an iterative algorithm. Subsequently, in order to give a suitable value λ, we need to determine the possible bound of the optimal value. Theorem 1 is proposed to solve this problem.

Theorem 1. Given the real symmetric matrix A ∈ Rd×d and the positive semi-definite matrix B ∈ Rd×d, rank(B) = r ≤ d. If W1 ∈ Rd×m1 , W2 ∈ Rd×m2 and m1 > m2 > d- r, then max WT 1 W1=I tr(WT 1 AW1) tr(WT 1 BW1) ≤ max WT 2 W2=I tr(WT 2 AW2) tr(WT 2 BW2) .

The proof of Theorem 1 is based on the following lemma: Lemma 2. If ∀i, ai ≥ 0, bi > 0 and a1

Proof. Let ak

b1 ≤ a2 b2 ≤ · · · ≤ ak bk , then a1+a2+···+ak b1+b2+···+bk ≤ ak bk .

bk =

q. So ∀i, ai ≥ 0, bi >

0, we have ai ≤

qbi. Therefore a1+a2+···+ak

b1+b2+···+bk ≤ ak

bk

Now we give the proof of Theorem 1 in the following.

Proof of theorem 1. Suppose W1

and W∗

1 = arg max

1 W1=I

WT

tr(WT

∗

= [w1,w2, ...,wm1 ]

1 AW1)

tr(WT

m1 =

1 BW1) . Let Cm2

without loss of generality, we suppose tr(WT p(1)AWp(1)) tr(WT p(1)BWp(1)) ≤

where Wp(i) ∈ Rd×m2 is the i-th combination of w1,w2, ...,wm1 with m2 elements(note that m1 > m2), so the number of combinations is h.

Let Cm2-1 m1-1 = l, note that each of wj(1 ≤ j ≤ m1) occurs l times in {Wp(1),Wp(2), ...,Wp(h)}. According to Lemma 2, we have

According to Theorem 1 we know, with the reduced dimension m increases, the optimal value is decreased monotonously. When m = d, the optimal value is equal to tr(A) tr(B) . So max WT W=I tr(WT AW) tr(WT BW) ≥ tr(A) tr(B) .

On the other hand, max WT W=I tr(WT AW) = ∑mi=1 αi, and min WT W=I tr(WT BW) = ∑mi=1 βi, where α1, α2, . . . , αm are the first m largest eigenvalues of A, and β1, β2, . . . , βm are the first m smallest eigenvalues of B. Therefore, max WT W=I tr(WT AW) tr(WT BW) ≤ α1+α2+···+αm β1+β2+···+βm .

As a result, the bound of the optimal value is given by tr(A) tr(B) ≤ max WT W=I tr(WT AW) tr(WT BW) ≤ α1+α2+···+αm β1+β2+···+βm

h, Now, we obtain an iterative algorithm for obtaining the optimal solution, which is described in Table 1. From the algorithm we can see, only a few iterative steps are needed to obtain a precise solution. Note that the algorithm need not calculate the inverse of B, and thus the singularity problem does not exist in it naturally.

In this case, when W lies in the null space of matrix B, then tr(WT BW) = 0, the value of the objective function becomes infinite. Therefore, we can reasonably replace the optimization problem with V∗ = arg max VT V=I tr(VT (ZT AZ)V), where V ∈ R(d-r)×m, and Z = [z1, z2, ..., zd-r] are the eigenvectors corresponding to d - r zero eigenvalues of B.

We know that V∗ = [μ1 ,μ2 , ...,μm], where μ1,μ2, ...,μm are the first m largest eigenvectors of ZT AZ. So, in this case, the final solution isW∗ = Z · V∗

Input:

The real symmetric matrix A ∈ Rd×d and the positive semi-definite matrix B ∈ Rd×d, rank(B) = r ≤ d. The error constant ϵ.

Output:

Projection matrix W∗ , where W∗ ∈ Rd×m and W∗T W∗ = I.

In the case of : m> d- r.

1.λ1 ← tr(A) tr(B) , λ2 ← α1+α2+···+αm β1+β2+···+βm , λ ← λ1+λ2 2 , where α1, α2, ..., αm are the first m largest eigenvalues of A, β1, β2, ..., βm are the first m smallest eigenvalues of B.

1. Calculate γ, where γ is the sum of the first m largest eigenvalues of A - λB.
2. If γ > 0, then λ1 ← λ, else λ2 ← λ.
3. λ ← λ1+λ2 2 .

End while.

W∗ = [ν1, ν2, ..., νm], where ν1, ν2, ..., νm are the first m largest eigenvectors of A - λB.

In the case of : m≤ d- r.

W∗ = Z · [μ1, μ2, ..., μm], where μ1, μ2, ..., μm are the first m largest eigenvectors of ZT AZ, and Z = [z1, z2, ..., zd-r] are the eigenvectors corresponding to d - r zero eigenvalues of B.

Table 1: The algorithm for the optimization problem

4 NeighborhoodMinMax Projections

The method of Neighborhood MinMax Projections(NMMP) is described in Table 2 . In order to speed up, PCA can be used as a preprocessing step before performing NMMP.

Denote the covariance matrix of data by St, and denote the null space of St by φ, the orthogonal complement of φ by φ⊥.

1. Preprocessing: eliminate the null space of the covariance matrix of data, and obtain new data X = [x1, x2, ..., xn] ∈ Rd×n, where rank(X) = d
2. Input:

1. calculate ˜Sw and ˜Sb according to Eq.(3) and Eq.(5)
2. calculate W using the algorithm described in Table 1
3. Output:

Table 2: Algorithm of NMMP

It is well known that the null space of St can be eliminated without lose of any information. In fact, it can be easy to prove that the null space of St comprises the null space of ˜Sw and the null space of ˜Sb defined in Section 2. Suppose w ∈ φ⊥ and ξ ∈ φ, then

Eq.(10) demonstrates that eliminating the null space of the covariance matrix of data will not affect the result of the proposed method. Thus we use PCA to eliminate the null space of the covariance matrix of data.

5 Discussion

Our method is closely connected with LDA. Both of them are supervised dimensionality reduction methods, and the goals are also similar. They both try to maximize the scatter between different classes, and minimize the scatter within the same class. The matrix ˜Sb defined in Eq.(5) and ˜Sw defined in Eq.(3) are parallel to the between-class scatter matrix Sb and within-class scatter matrix Sw in LDA respectively. In fact, when the number of neighbors reaches the number of the total available neighbors(kw(i) = ni - 1, and kb(i) = n-ni, where ni is the data number of class i, and n is the number of total data), we have ˜Sb + ˜Sw = n2St, which is similar to Sb +Sw = St in LDA.

However, in comparison with LDA, we do not impose the faraway pairwise points within the same class to be close to each other, which makes us focus more on the improvement of the discriminability of local structure. This property is especially useful when the distribution of class data is more complex than Gaussian. We give a toy example to illustrate it(Figure 2).

The toy data set consists of three classes(shown by different shapes). In the first two dimensions, the classes are distributed in concentric circles, while the other eight dimensions are all Gaussian noise with large variance. Figure 2 shows the two-dimensional subspace learned by PCA, LDA and NMMP, respectively. It illustrates that NMMP can find a low-dimensional transformation preserving manifold structure with more discriminability.

Moreover, compared with LDA, our method is able to extract more discriminative features and the singularity problem existing in LDA will not occur naturally.

Distance metric learning is an important problem for the distance based classification method. Learning a Mahalanobis distance metric is to learn a positive semidefinite matrix M, and using the Mahalanobis distance metric (xi - xj)T M(xi - xj) to replace the Euclidean distance metric (xi - xj)T (xi - xj), where M ∈ Rd×d, xi, xj ∈ Rd . Note that M is positive semidefinite, with the eigen-decomposition, M = VVT , where V = [σ1v1, σ2v2, ..., σdvd], σ and v are eigenvalues and eigenvectors of M. Therefore, the Mahalanobis distance metric can be formulated as (VT xi -VT xj )T (VT xi -VT xj). In this form, we can see that Learning a Mahalanobis distance metric is to learn a weighted orthogonal linear transformation. NMMP learns a linear transformationWwith the constraint of WT W = I. So it can be viewed as a special case of learning a Mahalanobis distance metric, where the weight value σi is either 0 or 1. Note that directly learning the matrix M is a very difficult problem and it is usually formulated as a semidefinite programming (SDP) problem, where the computation burden is extremely heavy. However, if we learn the transformationV instead of learning the matrixM, the problem will become much easier to solve.

6 Experimental Results

We evaluated the proposed NMMP algorithm on several data sets, and compared it with LDA and LMNN method. The data sets we used belong to different fields, a brief description of these data sets is list on Table 3.

We use PCA as the preprocessing step to eliminate the null space of data covariance matrix St. For LDA, due to the singularity problem existing in it, we further reduce the dimension of data such that the within-class scatter matrix Sw is nonsingular.

In each experiment, we randomly select several samples per class for training and the remaining samples for testing. the average results and standard deviations are reported over 50 random splits. The classification is based on k-nearest neighbor classifier(k = 3 in these experiments).

It is worth noting that the parameters in our method are not sensitive. In fact, in each experiment, we simply set kb(i) to 10, and set kw(i) to ni/2+2 for each class i, where ni is the training number of class i.

The experimental results are reported in Table 4. We use the recognition result directly performed after the preprocessing by PCA as the baseline.

In the following we describe the details of each experiment.

The UCI data sets

In this experiment, we perform on two small data sets, Iris and Balance, taken from the UCI Machine Learning Repository 1 . As the class distributions of this two data sets are not very complex, LDA works well, and our method also demonstrates the competitive performance.

Face recognition

The AT&T face database (formerly the ORL database) includes 40 distinct individuals and each individual has 10 different images. Some images were taken at different times, and have variations [Samaria and Harter, 1994] including expression and facial details. Each image in the database is of size 112 × 92 and with 256 gray-levels.

1 Available at http://www.ics.uci.edu/ mlearn/MLRepository.html

In this experiment, no other preprocessings are performed except the PCA preprocessing step. The result of our method is much better than those of LDA and the baseline. LMNN have a good performance too, but the computation burden is extremely heavy.

We also perform the experiments on many other face databases, and obtain the similar results, say, our method demonstrates the much better performances uniformly.

Object recognition

The COIL-20 database [Nene et al., 1996] consists of images of 20 objects viewed from varying angles at the interval of five degrees, resulting in 72 images per object.

In this experiment, each image is down-sampled to the size of 16 × 16 for saving the computation time.

Similar to the face recognition experiments, the results of our method and LMNN are much better than those of LDA

and the baseline. Note that both face images and object images distribute on an underlying manifold, the experiments verify that NMMP can preserve manifold structure with more discriminability than LDA.

Digit recognition

In this experiment, we focus on the digit recognition task using the USPS handwritten 16 × 16 digits data set 2 . The digits 1,2,3,and 4 are used in this experiment as the four classes. There are 1269, 929, 824 and 852 examples for each class, with a total of 3874.

On this data set, the baseline already works well, and our method still makes a little improvement. LDA fails in this case, which demonstrates that the available projection number of LDA may be insufficient when the data distributions are more complex than Gaussian.

Text categorization

In this experiment, we investigated the task of text categorization using the 20-newsgroups data set 3 . The topic rec which contains autos, motorcycles, baseball, and hockey was chosen from the version 20-news-18828. The articles were preprocessed with the same procedure as in [Zhou et al., 2004]. This results in 3970 document vectors in a 8014dimensional space. Finally the documents were normalized into TFIDF representation.

Our method and LMNN both bring significant improvements comparing with the baseline on this data set. In comparison, the performance of LDA is limited in that the available projection number of it is only 3, which is insufficient for this complex task.

7 Conclusion

In this paper, we propose a new method, Neighborhood MinMax Projections (NMMP), for supervised dimensionality reduction. NMMP focuses only on the pairwise points where the two points are neighbors of each other. After the dimensionality reduction, NMMP minimizes the distance of the considered pairwise points within the same class, and maximizes the distance of those between different classes. In comparison with LDA, NMMP focuses more on the improvement of the discriminability of local structure. This property is especially useful when the distribution of class data is more complex than Gaussian. Toy example and real world experiments are presented to validate it. Moreover, other disadvantages of LDA i.e., the singularity problem of Sw and the limitation of the available number of dimension are also avoided in our method.

As a linear dimensionality reduction method, NMMP can be viewed as a special case of learning a Mahalanobis distance metric. Usually, the computation burden of learning a Mahalanobis distance metric is extremely heavy. Our method formulates the problem as a constrained optimization problem, and the global optimum can be effectively and efficiently obtained. Experiments demonstrate that our method has a competitive performance compared with the recent proposed Mahalanobis distance metric learning method, LMNN, but the computation cost is much lower.

2 Available at http://www.kernel-machines.org/data

3 Available at http://people.csail.mit.edu/jrennie/20Newsgroups/