Discriminative Feature Selection via A Structured Sparse Subspace Learning Module

1 Introduction

Data in many areas are represented by high-dimensional features, and a natural question emerges out: what are the most discriminative features in the high-dimensional data? Feature selection plays a crucial role in machine learning which endeavors to select a subset of features from the highdimensional data for improving data compactness and reducing noisy features so that the over-fitting, high computational consumption and low performance issues can be alleviated in many real-world applications.

Feature selection approaches can be roughly divided into three categories, i.e., filter methods [Gu et al., 2012], wrapper methods [Maldonado and Weber, 2009] and embedded methods [Xiang et al., 2012; Ming and Ding, 2019; Tian et al., 2019]. Wherein, filter methods focus on evaluating the correlation of features with respect to class label of data, which results in that the correlated features are redundant. Wrapper methods measure the importance of features according to classification performance, so as to the computational cost of wapper model is very high and the representability of learned features may be poor in other tasks. In contrast, embedded methods gain more attentions since it incorporates feature selection and classification model into a unified optimization framework by learning sparse structural projections. Besides, some feature selection approaches based on neural network [Zhang et al., 2020] and autoencoder model [Han et al., 2018; Abid et al., 2019] pay overmuch attentions on minimizing the reconstruct errors, which possible not benefit for classification task and tuning parameters may affect its efficiency of practical applications.

∗Corresponding Author

Most of embedded methods commonly employ ℓ2,1-norm regularization to improve the row-sparsity in learned subspace. Concretely, Nie et al. propose robust feature selection model (RFS) [Nie et al., 2010] which simultaneously incorporates ℓ2,1-norm into loss function and regularization term to overcome outliers issues and make feature selection come true. Similar, [He et al., 2012] proposes to use Correntropy Robust loss function in Feature Selection (CRFS) to enhance model's robustness. Additionally, [Xiang et al., 2012] develops a feature selection effort joints discriminative least square regression (DLSRFS) model and sparse ℓ2,1-norm regularization into an unified framework. Recently, RJFWL [Yan et al., 2016] aims to learn the ranking of all features by imposing ℓ2,1-norm constraint and non-negative constraint on learned feature weights matrix. Regrettably, all aforementioned approaches are always unable to escape the dilemma of tuning the trade-off parameter between loss term and sparse regularization term. Moreover, the compactness of sparse row is sensitive to the trade-off parameter, namely, the uniqueness of optimum is difficult to guarantee.

An illustrative example is depicted in Figure 1, from which we can figure out that ℓ2,1-norm calculates the score of features integrally, then ranks them resorting to the learned scores. By doing so, the feature ranking results may be changed dramatically when the feature score vary subtly. In contrast, ℓ2,0-norm obtains the score matrix only depend upon the k features, and no further ranking operation is required. Based on these, Cai et al. [Cai et al., 2013] questioning that whether the method based on convex problem is always better than that based on non-convex problem. They exploit top-k features selection method (RPMFS) based on least square regression model with ℓ2,0-norm equality constraint. Pang et al. [Pang et al., 2018] develop an Efficient Sparse Feature Selection (ESFS) via ℓ2,0-norm constraint based on least square regression model as well. In [Du et al., 2018], Unsupervised Group Feature Selection (UGFS) algorithm selects a group of features initially then update the selection until a better group appears by using ℓ2,0-norm constraint. Nevertheless, both of above methods optimize the ℓ2,0-norm constraint problem by using Augmented Lagrangian Multiplier (ALM) and gradient descent optimization algorithms, which may result in that the solution is sensitive to the initialization and easy to stuck in local optimum. Ultimately, the subsequent performance of classification is fluctuant.

Facing these obstacles that hinder the development of feature selection, we entail solving ℓ2,0-norm constraint problem along with rigorous theoretical guarantee firstly in this paper. Then, we propose a new discriminative feature selection model via orthogonal ℓ2,0-norm constraint. Finally, we analyze the performance of proposed method from two aspects, i.e., discriminability and efficiency. Our contributions can be summarized as follows:

• We propose a Structured Sparse Subspace Learning (S3L) module to solve the subspace sparsity issue with theoretical guarantee, which is first presented in our paper according to our best knowledge.
• We elaborately design a discriminative feature selection model, named Subspace Sparsity Discriminant Feature Selection (S2DFS) that integrates trace ratio formulated objective and structured sparse subspace constraint for acquiring more informative feature subset.
• Optimizing the proposed S2DFS model is a NP-hard problem. We provide an alternative iterative optimization algorithm to solve trace ratio problem, in which the proposed S3L module is employed to acquire a closedformed solution rather than an approximate one, so as to ensure the stability of selection results.
• Experimental results show the effectiveness of proposed optimization algorithm in two perspectives, i.e., performance: our method outperforms other related SOTA

feature selection methods in terms of classification on several real-world datasets; convergent speed: our algorithm reaches convergence within few iterations.

2 RelatedWorks

In this section, we briefly review constrained Linear Discriminative Analysis (LDA) in trace ratio formulation, then the subspace sparsity issue will be introduced.

2.1 Constrained Trace Ratio LDA

Given the training data X = {x1, x2, ..., xn} ∈ Rd×n , where xi ∈ Rd×1 and the label vector is y = {y1, y2, ..., yn}, where yi ∈ {1, 2, ..., c}, and c denotes the total class number. LDA aims to find a transformation matrix W ∈ Rd×m to map original high-dimensional sample x ∈ Rd×1 into lowdimensional subspace with m dimensions where m << d. For improving the discriminative power of model, the samples within the same class are pulled together as well as those points between different classes are pushed far away in learned subspace ideally. There are many kinds of objective function of LDA [Bishop, 2006], one of the most discriminative objectives is the following constrained trace ratio one

where Sb = ∑ci=1 ni (µi -µ) (µi -µ)T and Sw = ∑ci=1 ∑xj∈πi (xj -µi) (xj -µi)T denote between-class and within-class scatter matrix respectively. ni is the number of samples belong to i-th class, µi denotes the mean of samples in i-th class, µ is the mean of total samples and πi denotes the set of i-th class points. Orthogonal constraint is used to avoid trivial solution [Nie et al., 2019]. Note that, optimizing such model in Eq.(1) is a non-convex optimization problem that does not have a closed-form solution, and several attempts [Wang et al., 2007; Nie et al., 2019] have tried to solve it. In what follows, we will provide an alternative iterative optimization algorithm to get the solution.

2.2 Subspace Sparsity Problem

Subspace sparsity issue is originally derived from the sparse principal subspace estimation [Vu et al., 2013] which can be generally formulated as

where m ≤ k ≤ d and A ∈ Rd×d is a positive semidefinite matrix. ‖W‖2,0 = ∑di=1 ∥∥∥∑mj=1 w2 ij ∥∥∥0 denotes the number of non-zero rows in matrix W, and the constraint ‖W‖2,0 = k forces the learned subspace is row sparse, in which the index of non-zero rows equals to the index of selected features. Since both the orthogonal and row sparsity constraints are the non-convex constraint, optimizing problem (2) is challenging. In the following Section 3.1, we will propose the S3L module to solve it with theoretical guarantee.

3 ProposedMethod

3.1 Structured Sparse Subspace LearningModule

In this section, we propose a novel S3L module to solve problem (2) that has been tentatively solved in [Wang et al., 2014; Yang and Xu, 2015], but none of them provides deterministic guarantee on global convergence. We first develop a straightforward strategy for solving the case that rank (A) ≤ m. Then, an iterative optimization Algorithm 1 is derived to solve the general case that rank (A) > m in which a lowrank proxy covariance P is designed to approximate A. For readability, we define two definitions that will be used in following deducing:

Definition 3.1 Given an indices I whose elements span from 1 to d. Then, the row extraction matrix Qkd ∈ Rd×k is defined as:

Define Qd,k(I) = Qkd as an operator that inputs indices I and outputs row extraction matrix Qkd ∈ Rd×k .

Case 1: rank (A) ≤ m

We first consider the simplest problem, i.e., k = m, thus, the problem (2) reduces to

Due to the constraint WT W = Im×m and ‖W‖2,0 = m, W can be rewritten as W = Qmd V, where V ∈ Rm×m and satisfies VT V = Im×m. Additionally, Qmd = Qd,m(I) ∈ Rd×m is the row extraction matrix that has been defined. Then, problem (3) has been reduced to solve

where ˜A = Qmd T AQmd . So far, the row extraction matrix Qmd is still unknown. As a consequence, considering V ∈ Rm×m is a square matrix and VT V = VVT = I, we can infer the fact Tr(VT ˜AV) = Tr(˜AVVT ). Thus, problem (4) can be reduced to

which has a globally optimal solution that m largest diagonal elements of A, and Qmd is generated by operator Qd,m(I) with input I that equals to the indices of m largest diagonal elements of A.

Next, we consider the case that k ̸= m, and problem in Case 1 becomes to

We use similar technique to solve problem (6) as

which can be easily solved globally by selecting k largest diagonal elements of A, and Qkd = Qd,k(I) where I denotes the indices of k largest diagonal elements of A. Throughout above analysis, our S3L module is able to efficiently (non-iteratively) obtain globally optimal solution of NPhard problem (2) when rank(A) ≤ m.

Nevertheless, the Case 1 is not usually occurred in practice, since PCA often used as a data-preprocessing technique to remove null space of data, which cause a general Case 2, i.e., rank(A) = d > m. Next, we develop an iterative S3L module to solve the problem (2) in the following Case 2.

Case 2: rank (A) > m

We solve the problem (2) in case 2 based on the Minimization-Majorization framework (MM) [Sun et al., 2016] whose key insight is to maximize a lowerbound function of target problem. Therefore, how to design a suitable lowerbound function is the key to the success of MM. Fortunately, through a lot of attempts, we elaborately design an ideal lowerbound function of problem (2) as follow:

where K = {W ∈ Rd×m|WT W = I, ‖W‖2,0 = k} is the structured sparse subspace constraint, the surrogate variable Γt = AWt (WT t AWt )† WT t A and (·)† denotes the Moore-Penrose pseudoinverse. As a lowerbound function, a deserved property should be invariably held:

which will be proved in Section 3.2.

According to MM framework, solving problem (2) in Case 2 becomes to solve

which is still a difficult problem that can not be optimized directly. Now, we discuss some properties of designed surrogate variable Γt through the following Theorem 1.

Theorem 1 In the t-th iteration (t ≥ 1), the conditions rank (Γt) ≤ m and Γt ⪰ 0 will always hold.

Proof 1 According to the property of matrix rank rank(AB) ≤ min{rank(A), rank(B)}, we have rank(Γt) ≤ rank(Wt) = m, then the first condition has been proved. Subsequently, we suppose that Ω = AWt (WT t AWt )† WT t A 1 2 , then according to the fact that A ⪰ 0, A is a symmetric matrix and B† = B†BB† , we can obtain

(11)

which completes the proof of Theorem 1.

According to Theorem 1, we wondrously discover an interesting phenomenon that the designed surrogate variable Γt satisfies all conditions of A in case 1, i.e, rank(A) ≤ m and Γt ⪰ 0. Those facts enlighten a thought that problem (10) can be solved by using S3L module iteratively, which is summarized in the following Algorithm 1.

Algorithm 1 Solve problem (2), when rank(A) > m

Input: A, k, m, d andW0;

Initialization: W0 = rand(d, m), t = 0;

Repeat:

1. Γt ←AWt (WT t AWt )† WT t A;
2. It ←sort ( diag(Γt), 'descend', k );
3. Qt kd ←Qd,k(It);

1. Wt+1 ← Qt kdVt;
2. t = t +1;

Until convergence

Output: W∗ ∈ Rd×m.

3.2 AlgorithmAnalysis

We prove the convergence of the proposed Algorithm 1 through the following Theorem 2.

Theorem 2 Algorithm 1 increases the objective function value of Eq.(2) in each iteration, when rank(A) > m.

Proof 2 Note the fact that B = BB†B, we can infer that

where rank (Γt) ≤ m. Then, we can use the S3L module to maximize the above Eq.(12) as

where ̂WT t+1 is generated according to Eq.(4) with input A = Γt. Thanks to the commutative law in algebraic operation of matrix trace and A is a symmetric matrix, Eq.(13) can be decomposed as

where two symmetric matrices are respectively defined as Φ = A 12Wt (WT t AWt )† WT t A 1 2 ∈ Rd×d , Θ = A 12 ̂Wt+1 ̂WT t+1A 12 ∈ Rd×d . Considering Theorem 4.3.53 in [Horn and Johnson, 2012], we have

which provides an evidence of correctness for Eq.(9). Moreover, note that rank(Θ) ≤ rank(̂Wt+1) = m, then we have ∑mi=1 λi (Θ) = Tr (Θ). Combining Eq.(12), Eq.(13) and Eq.(15), we can conclude that

According to the previous definitions, i.e., ̂Wt+1 = Qt kd ̂Vt and Wt+1 = Qt kdVt, we can finally obtain that

then the proof of Theorem 2 has been completed.

3.3 Subspace Sparsity Discriminant Feature Selection

In this part, we propose to introduce a novel feature selection method named S2DFS, and optimize it by using S3L module.

Model Formulation. In order to improve model's discriminative power, S2DFS leverages trace ratio LDA model in Eq.(1) for pulling samples within same class together and pushing those samples between different classes far away in learned subspace. Besides, we use the structured sparse subspace constraint to directly facilitate the exploration and inference of important features of data. In light of these points, the objective function of S2DFS model can be formulated as

where m ≤ k < d, m is the number of selected features and k means the number of non-zero rows in W. Due to the nonlinear constraints and trace ratio objective, optimizing such a maximization ratio problem in Eq.(18) is challenging.

Optimization Algorithm Before solving problem (18), we first consider to solve a general maximization ratio problem as follow:

where x ∈ C is arbitrary constraint on x and as the denominator, g(x) > 0. In Algorithm 2, we summarize an alternative optimization algorithm to solve the general problem (19), and a series of theorems and proofs are presented to guarantee its convergence.

Algorithm 2 Algorithm to solve the general maximization ratio problem (19).

1. Initialization: x ∈ C, t = 1.

While not converge do

1. Calculate λt = f(xt) g(xt) .
2. Calculate xt+1 = argmax x∈C f(x) -λtg(x).
3. t = t +1. end while

Theorem 3 The global solution of the general maximization ratio problem (19) is equivalent to the root of following function:

Proof 3 Supposing x∗ is the global solution of problem (19) and its corresponding maximal objective function value is λ∗ , then the following equation holds: f(x∗ ) g(x∗ ) = λ∗ . As a result,

g(x) > 0, we can infer that f(x)-λ∗ g(x) ≤ 0 which means:

Consequently, the global maximal objective function value λ∗ of problem (19) is the root of function h(λ). Thus, the proof of Theorem 3 is completed.

□

Theorem 4 Algorithm 2 increases the objective function value of problem (19) in each iteration until it reaches convergence.

Proof 4 According to the step 2 in Algorithm 2, we have f(xt) - λtg(xt) = 0, and from step 3, we can infer that f (xt+1) - λtg(xt+1) ≥ f (xt) - λtg(xt). Combining above two inequalities, we can obtain f(xt+1) - λtg(xt+1) ≥ 0 that equals to f(xt+1) g(xt+1) ≥ λt = f(xt) g(xt) . That is, Algorithm 2 increases the objective function value in Eq.(19) in each iteration and the proof of Theorem 4 has been finished.

□

Systematically, it is worth noting that Theorem 3-4 provide a complete framework to optimize the general maximization ratio problem in Eq.(19) with rigorously convergent proof.

According to the Algorithm 2, we can infer that the key of solving proposed maximization ratio problem in Eq.(18) is to solve the following problem:

(22)

max W∈Rd×m T r (WT BW) s.t. WT W= I, ‖W‖2,0 = k, where B = Sb - λSw + ηI is a positive semi-definite matrix if η is large enough. In theory, the relaxation parameter η is easily to be calculated as the largest eigenvalue of Sb -λSw. The first way that can be thought to acquire η is eigenvalue decomposition, which however is time-consuming in dealing with high-dimensional data. Here, we employ a power iteration method [Nie et al., 2017] to solve above issue, which can obtain the relaxation parameter η faster than eigenvalue decomposition and improves the applicability of our method.

Note that problem (22) is a special case of problem (2), which can be directly addressed by using proposed S3L module and obtain a closed-form solution.

4 Experiments

In this section, for verifying the effectiveness of proposed method, we first present an visualization experiment on Mnist dataset. Then, we evaluate its classification performance on nine publicly available high-dimensional datasets. Finally,

we exhibit the convergent speed of our method.

4.1 Preliminary

Datasets. We evaluate the performance of proposedmethod on several high-dimensional real-world datasets, and more details about them are shown in Table ??. For the color image datasets, i.e., Pubfig [Xu et al., 2018], we firstly downsample each image into a suitable scale and extract 100 LOMO features [Liao et al., 2015] for data representation.

1 http://qwone.com/∼jason/20Newsgroups/

Experimental setting. We compare our method to several SOTA feature selection methods including DLSRFS [Xiang et al., 2012], RFS [Nie et al., 2010], CRFS [He et al., 2012], RJFWL [Yan et al., 2016], infFS [Roffo et al., 2015] and ESFS [Pang et al., 2018] in the classification task. We use the k-nearest neighbor algorithm as the classifier, all experiments are repeatedly conducted 10 times and the average recognition rate and standard deviation are recorded as the measurement of performance for all competitors.

4.2 Visualization Experiment

In order to intuitively show the performance of feature selection algorithm, we provide an visualization experiment conducted on Mnist digits dataset. Figure 2b exhibits the two types of selected features, i.e., digital pixels (white dots) and background pixels (red dots) respectively. In the Figure 2c,

we show the selected features lied on the Mnist digits, from

which we can observe that the shape of digits can be dis-

tinctly reconstructed by selected features generated by pro-

posed method. Experimental results verify the discriminabil-

ity of selected features that learned from our method.

(a)

(b)

4.3 Classification Experiment

Figure 3 shows the statistical estimation of classification results on six high-dimensional real-world datasets, which demonstrates the discriminability of feature subsets selected by our method. Concretely, in general, comparing to all competitors, our method always achieves comparable or even highest mean accuracies on almost all datasets with different number of selected features. Specifically, due to the nonparameter property of our method, it achieves stable performance on all datasets. Furthermore, our method outperforms ESFS in different degrees on various datasets, which results from that ESFS applies approximate gradient descent algorithm to obtain a suboptimal solution, while our method obtains a closed-form solution, which guarantees stability of our method's performance.

2 http://www.cs.cmu.edu/afs/cs/project/theo-20/www/data/

3 http://www.escience.cn/system/file?fileId=82035

4.4 Convergence Analysis

Figure 4 plots the convergence curves of our iterative algorithm with different number of selected features m and model's sparsity k on 20News and WebKB datasets. Overall, our algorithm reaches convergent within 15 iterations in most cases or even less than 5 iterations on 20News dataset. In detail, a subfigure embedded in each figure is used to facilitate observation of the curve convergence when k = m. We discover the fact that the curve reaches convergent with one iteration on all datasets, which verifies the high-efficiency of proposed optimization algorithm.

5 Conclusion

In this paper, we have proposed explicit discriminative feature selection algorithm via employing structured sparse subspace constraint, which is a NP-hard optimization problem. To our best knowledge, there are no techniques able to obtain closed-form solution and proves convergence. As the major contribution of this paper, we provide an ideal iterative optimization algorithm that can directly solve proposed NP-hard problem and reaches convergent very fast in both theory and practice. Experimental results demonstrate the effectiveness of proposed method. In future works, we intend to develop nonlinear feature selection model via proposed S3L module for future improving discriminability of selected features.

Acknowledgements

This work was supported in part by the National Key Research and Development Program of China under Grant 2018AAA0101902, in part by the National Natural Science Foundation of China under Grant 61936014, Grant 61772427 and Grant 61751202, and in part by the Fundamental Research Funds for the Central Universities under Grant G2019KY0501.