models, which are easier to interpret. Such models can then be analyzed by domain experts and are easier to validate. Getting more interpretable models is also a key concern nowadays and even considered by many as a requirement when deployed in the medical domain.
Feature selection has been already largely studied. Yet, current methods are still widely unsatisfactory mainly because of the typical instability they exhibit. Instability here refers to the fact that the selected features may be drastically different for similar data, even though the true underlying processes (explaining the target variable) are essentially constant. Such instability is a key issue as it reduces the interpretability of the predictive models as well as the trust of domain experts towards the selected feature subsets. We address this problem here by designing methods balancing between the classification performance and the selection stability of the well-known Recursive Feature Elimination (RFE) algorithm. Our approach allows domain experts to explicitly control the trade-off and to select Pareto-optimal compromises based on their personal preferences.
In the rest of this section, two distinct stability problems that are tackled in this paper are introduced. 1.1
The Stability Problems Single Task Stability (1) Feature selection methods are often inherently unstable, i.e. they return highly different feature sets when the training data is slightly modified. Figure 1a illustrates such an instability. The initial dataset is perturbed1 to form different datasets. Instability arises when little overlap of the selected features occurs. This prevents a correct and sound interpretation of the selected features and strongly impacts their further validation by domain experts. Unlike optimizing the accuracy of predictive models, optimizing selection stability may look trivial since an algorithm always returning an arbitrary but fixed set of features would be stable by design. Yet, such an algorithm is not expected to select informative and predictive features. This illustrates that optimizing stability is only well-posed jointly with predictive accuracy, and possibly additional criteria such as minimal model size or sparsity.
Transfer Learning Selection Stability (2) Multi-task feature selection aims at discovering variables that are relevant for several similar, yet distinct, classification tasks. Different feature subsets can be returned for each task. In this paper, we focus on the case where all learning tasks are not directly available. Information from the tasks arising first can be propagated to subsequent tasks, via transfer learning. Stability has to be encouraged from the domain expert point of view as features that are relevant for different data sources are likely to be particularly interesting to study. The accuracy-stability trade-off on such a learning problem (represented in Figure 1b) can take two extreme values. With complete disregard to stability, each feature set could be selected on a given task 1 Here by bootstrapping which is often used to measure such instability, but it could be any small perturbation. independently of the others, with no control on the across task stability. On the contrary, maximum stability can trivially be reached by returning the feature set computed for the first task, for all subsequent tasks. However, this is expected to reduce the accuracy of the models built on these subsequent tasks as previously learned features might turn out to be less informative for them. This would be the case if the different tasks are obtained by gradually enriching or correcting the data as features learned on the error-corrected data are expected to be more relevant.
(a) Single-task (b) Transfer learning
In section 2, feature selection methods and propositions to increase stability are reviewed. Section 3 introduces a metric to assess the performance of methods compromising between feature selection stability and classification performance. Then a biased variant of the RFE algorithm is proposed in section 4. Section 5 demonstrates the ability of this biased RFE to tackle the previously mentioned stability problems. 2
Feature selection techniques are generally split into three categories: filters,
wrappers and embedded methods. Filters evaluate the relevance of features
independently of the final model, most commonly a classifier, and remove low ranked
features. Simple filters ( e.g. t-test or ANOVA) are univariate, which is
computationally efficient and tends to produce a relatively stable selection but they
E4xplicit VC.oHntarmolero,fPF.eDautpuornetRelevance and Stability
plainly ignore the possible dependencies between various features. Information
theoretic methods, such as MRMR [
Some works specifically study the causes of selection instability. Results show
that it is mostly caused by the small sample/feature ratio [
Looking for a stable feature selection also requires a proper way to quantify
stability itself and lots of measures have already been proposed: the Kuncheva
index [
1 Pd 2
d f=1 sf
φ = 1 − k k )
d ∗ (1 − d
(1)
with sf2 = MM− 1 pˆf (1 − pˆf ) the estimator of the variance of the selection of the
fth feature over the M selected subsets and k the mean number of features
selected from the original d features.2 This measure is the only existing measure
satisfying the 5 (good) properties described in [
Several authors already proposed different methods to increase stability. For
instance, instance-weighting for variance reduction [
Multi-task feature selection has already been largely studied [
A Multi-Objective Evaluation Framework Through Pareto Optimality In this section, we propose to use a classical evaluation framework in multiobjective optimization to assess the efficiency of methods balancing between classification performance and selection stability. An (accuracy,stability) pair 3 (a1, φ 1) dominates another pair (a2, φ 2) iff a1 ≥ a2 ∧ φ 1 ≥ φ 2 and at least one of the inequalities is strict (>). A given method m is able to generate some pairs Pm in the space of all possible pairs4 P = { (a, φ ) : 0 ≤ a, φ ≤ 1} . From the set of generated pairs Pm, the set of pairs that are not dominated by any other pair, 3 Common alternatives to the classification accuracy, such as specificity/sensitivity or
AU C, can also be used.
4 The careful reader may remember that the stability measure φ formally lies in the
[
P am, can be found. This set, called the Pareto set, defines a subspace where no point dominates any other point. A domain expert would then choose his favorite pair based on his personal preference towards classification performance and feature selection stability.
As performance metric, we propose the widely used hypervolume measure
[
An example of the DA metric can be seen in Figure 2. The blue method
starts from the left with a higher accuracy. It thus gains some area over the
red method. Nonetheless, the red method can reach higher stabilities without
dropping the accuracy as much as the blue one. Overall, the red method has a
larger DA. Note that this DA is also equal to the fraction of the total area that
is dominated by the method, or 1 minus the fraction of area that dominates the
method. Its value thus lies in the [
As noticed by [
DAP =
Z 1 Z 1
0
0
w(a, φ )fP (a, φ )dadφ
(2)
with w the weighting function and fP the attainment function which is equal
to 1 if (a, φ ) is dominated and 0 otherwise. To preserve the [
In order to evaluate the pair (a, φ ) ∈ P corresponding, for instance, to a set of meta-parameters, some data has to be used to learn the features and some data to evaluate them on independent examples. This can be done via standard crossvalidation or by bootstrapping. Each set of meta-parameters produces different pairs in P and their average value is reported. Concretely, each point in Figure 2 comes with an uncertainty linked to the sampling of the data. In the following, we define a confidence interval on the true value of DA based on the derivation of confidence regions for each Pareto-optimal pair.
Let A be the random variable representing the accuracy value measured on a given subsampling of the data and Φ be the corresponding stability value. Let P = (A, Φ ) be the multivariate random variable with the accuracy and stability as dimensions. Let us assume that the evaluation protocol produces B measurements of P for each Pareto-optimal point5, represented by the vector p. The Hotelling distribution T 2 is the multivariate counterpart of the Student’s t distribution, with which we can define confidence (here 2-dimensional) regions.
T 2 = B(p¯ − μ (p))0C− 1(p¯ − μ (p)) ∼ 2(B − 1)
B − 2 ∗ F 2,B− 2 with C the sample covariance matrix. It can be shown that T 2 is distributed like a Fisher distribution F2,B− 2. Thus,
P (p¯ − μ (p))0C− 1(p¯ − μ (p) ≤ 2 ∗ (B − 1)
B − 2 ∗ F 2,B− 2(α ) = 1 − α The inequality defines an ellipsoidal region, that is likely to cover μ (p). The center of the ellipsoid is p¯. The length of the axis and their angle can be found by computing the eigenvalues and eigenvectors of the sample covariance matrix C. To compute a confidence interval on the DA, the most dominant and dominated point of each ellipse are found and used to compute the upper and lower bound of the confidence interval (see Figure 3c and 3d for concrete examples in our experiments). (3) (4) 4
In this section, we propose a simple method to balance between the classification performance and selection stability of a logistic RFE algorithm. The RFE 5 B could be e.g. the number of bootstrap samples. E8xplicit VC.oHntarmolero,fPF.eDautpuornetRelevance and Stability algorithm was originally introduced with a hinge loss. We prefer here the logistic variant for an expected smoother control of the trade-off under study. RFE iteratively drops the least relevant features until the desired number of features k is reached. We opt here to drop a fixed fraction (20%) of the features at each iteration. The loss function that a logistic RFE minimizes for a binary classification task is the following, with n the number of samples, xi sample number i made of d features as dimensions, and yi its label.
n L = X log(1 + e− yi∗ (w∗ xi)) + λ | w| 2 i=1 (5) The weight vector w contains a weight assigned to each feature. The features are then ranked based on the absolute value of their weight, which represents the importance of the feature in the final prediction. The term λ | w| 2 of the loss function is a regularization term, preventing coefficients of the model to take too high values, which would most likely result in overfitting. In the classical approach, every feature is regularized by the same amount λ .
We propose to extend equation (5), such that the regularization term
becomes λ β | w| 2. The function of the vector β is to bias the selection towards
certain features via differential shrinkage. A feature f with a small β f is less
regularized and vice-versa. Its selection in the model is less penalized than a
feature with a higher β f . The search is thus biased towards features with small β f .
A similar differential shrinkage has already been applied to the `1-AROM and
`2-AROM methods [
Biased RFE for Single Task Feature Selection By varying the distribution of β , the accuracy-stability trade-off of the biased RFE can be controlled. The biased RFE is equivalent to a standard RFE when β = 1. Otherwise, the selection is biased towards features with a small β f . This is expected to increase stability at the possible cost of some classification performance, as uninformative features could be prioritized. In this initial approach, we decide to favor some features non-uniformly at random, following a gamma distribution.
β f ∼ Γ (α, 1) with α the shape of the gamma distribution, which controls the trade-off. All β f are post centered such that μ (β ) = 1. As α → ∞ , the gamma distribution tends to a Dirac delta, δ (α ). All features have then the same weight (equal to μ (β ) = 1) and no bias is put in the selection. As α → 0, the distribution of β departs from δ (α ) which increases the bias. Domains experts can thus play with the α values and therefore explicitly tune the trade-off between selection stability and prediction accuracy. Using Prior Knowledge The biased RFE can take advantage of available prior knowledge. If a ranking of the features is known, then the β f can be assigned such that this ranking is respected. If the prior knowledge is meaningfull, the selection is no longer biased towards arbitrary features, but towards features that are high in the ranking, and thus potentially informative. Another type of prior knowledge could be an unordered set of features that are suspected to take part in the process of interest. The lowest β f could then systematically be assigned to those features.
Biased RFE for Transfer Learning We are now interested in the across task stability that can be obtained via transfer learning. Tasks are thus ranked such that information from previous tasks can be used in the selection of features for subsequent tasks.6 In task i, features that have been returned for tasks 0..i − 1 should be prioritized over the rest, such that the feature stability is increased. Given the definition of stability used here (equation (1)), it is actually possible to compute the drop/gain in stability that the selection of a feature would cause. Intuitively, we propose to bias the selection, through a specific choice of the vector β , towards features that would cause the highest gain/lowest drop in stability if they were to be selected. Constant terms put aside7, each feature influences (negatively) the total stability by its variance in the selection sf2 ∝ pf (1 − pf ). Feature f is given an attractiveness score scf , expressed in equation (6).
scf = 2 with sf,no the selection variance of feature f assuming f is not selected in the 2 current task and sf,yes its selection variance if it were to be selected. N is the number of tasks for which feature sets have already been selected. This score is thus proportional to the difference of stability between the cases where the feature is selected for a given task and not. This is illustrated in table 1a where the current task is T 4. For instance, the selection of the feature F 2 in task T 4 would make its mean selection, pf , equal to 0.75. If it were not to be selected, pf would be equal to 0.5. The attractiveness score of F 2, scF 2 is actually positive, meaning that the selection of F 2 in T 4 would increase the measured stability.
The (N+1)2 factor of equation (6) is there to correct a downards tendency of
N scf when the index of the considered task increases. This is illustrated on table 1b. If feature f is selected in each task, scf would actually decrease which would decrease the bias. It can be shown that including the correction term leads to scf = 2 ∗ pf − 1 with pf the proportion of the selections of feature f in the past N tasks. Let S be the sum of the such selections. By definition, pf = NS , 6 Tasks can be ranked naturally from their chronological order or by the domain expert. 7 We purposely drop the M/(M − 1) term, for convenience. Also, the denominator kd (1 − kd ) is constant if the number k of selected features is fixed.
sf,no = N+1 ∗ (1 − NS+1 ) and sf2,yes = NS++11 ∗ (1 − NS++11 ). Thus, 2 S scf = to bias the selection towards previously selected features.8 With α t = 0, features are learned independently on each task. On the contrary, an increasing α t raises the bias towards features that were already selected in past tasks. Domain experts can thus tune the α t values to control the accuracy-stability trade-off in such a transfer learning setting. 5
In this section, we evaluate to what extent an actual compromise between
prediction accuracy and selection stability can be made with the proposed approaches.
Experiments are performed on two distinct tasks, prostate cancer diagnosis from
microarray data and handwritten digit recognition. The Prostate dataset
contains 12600-dimensional (microarray) gene expression data from 52 patients with
prostate tumors and 50 healthy patients [
Evaluation Methodology To obtain the results presented in the next sections, the following methodology has been used. Each (a, φ ) pair is obtained with a different α (problem 1) or α t (problem 2). For the single task stability problem, the β f are first sampled from the gamma distribution. Then, M bootstrap samples are built. k features are then selected using the proposed biased RFE on each bootstrap sample. For the transfer learning stability problem, a single bootstrap sample for each task is created. Features are selected from it, then β for the next task is computed according to equation (7). The final prediction model is learned by minimizing the classical, unbiased, logistic loss with a L2 regularization (see equation (5)) with a non-strongly fitted 9 regularization parameter λ . Every model is evaluated on its out-of-bag examples. The mean accuracy as well as the stability of the selected features are computed. As these values are obviously dependent on the sampling of β (problem 1) or the features learned on the first few tasks (problem 2), this procedure is repeated B times and the mean values are reported. The 95% confidence regions of the expected value of the accuracy-stability trade-off are computed as well as the confidence interval on DA described in section 3. Stability of the feature selection (x-axis on Figures 2, 3 and 4) has not to be confused with its corresponding uncertainty which is the width of the confidence regions along the x-axis. 5.2 Single Task Selection Stability The λ meta-parameter of the RFE formulation (equation (5)) has not been strongly optimized. A value of 0.1 which provides a good accuracy has been used for both tasks. To obtain the below graphs, the methodology detailed in section 5.1 has been used with M = 30, k = 20 and B = 100.
Results on both data sets can be seen in Figure 3. The blue curves are obtained without any prior knowledge. The top-left point of each subgraph corresponds to the (accuracy,stability) trade-off obtained with the classical logistic RFE method. Following Pareto lines from left to right, the shape α of the gamma distribution decreases. This makes the biased logistic RFE departs from its unbiased version which raises stability but reduces classification performance. As the method fails to reach maximum stability, it was extended with the trivial point (arand, 1), obtained by always returning the same arbitrary feature subset.10 It 9 Values used for λ are 0.1 for problem 1 and 1 for problem 2. The final classification algorithm does not influence the selection stability. It can thus be optimized to maximize the predictive accuracy only. 10 It is actually impossible to reach a maximum stability of 1 for a finite regularization parameter λ . In such a case, even with no regularization, a feature is not guaranteed to be always selected.
seems that, while it is possible to increase signicfiantly the stability without degrading too much the accuracy on the Gisette dataset (Figure 3b), it is not the case for the Prostate dataset (Figure 3a) where the accuracy drops directly.
To measure the effect of prior knowledge, N = 10 examples are sampled randomly. The 100 features with the highest variance are selected as part of the prior knowledge, here representing a set of potentially relevant features. As can be seen on Figure 3a and 3b, even such a small prior knowledge improves the Dominance Area considerably.
Figures 3c and 3d have been obtained with a small subset of the Pareto points. The ellipses are the 95% confidence regions of the expected value (on the β sampling) of the accuracy-stability trade-off. For large α values, the importance of β is reduced, and thus the uncertainty limited. As α decreases, this confidence region grows. The ellipses are also all inclined towards the right. This represents the covariance between the accuracy and stability for a single β sampling. If the sampling appears to be bad, .i.e. poor features are prioritized, poor accuracy and poor stability are obtained. The opposite is true for a good sampling. By using the top-right and bottom-left point of each ellipse, it is possible to derive a confidence interval on the true DA of the method on these datasets. 5.3
Multi-task Selection Stability via Transfer Learning
To generate different, yet similar, classification tasks, normally distributed noise has been added to the Prostate dataset. This noise is centered on 0 and has a specific standard deviation for every couple of feature and task, such that features relevant in some task, could be irrelevant in others. Yet, tasks are expected to share common informative features. 8 tasks are considered here, with an arbitrary order between them. Results with k = 10, B = 500, λ = 1 are shown on Figure 4a. The blue area is obtained by combining two trivial options. First, the features learned on the first task can be selected for all subsequent tasks, achieving a stability of 1. Or features can be learned independently from each other (equivalent to α t = 0). This strategy offers poor selection stability, but also a sub-optimal classification performance. Knowledge from previous tasks can be used to guide the search towards potentially good features for subsequent tasks. This increases both the accuracy and stability at first. Then, the accuracy starts to decrease, as the selection of features is forced too much. This tendency is
better illustrated in Figure 4b, which contains some non-Pareto optimal points.
This result is consistent with the conclusion drawn by the analysis of the
GroupLasso with `1/`p regularization [
The typical instability of standard feature selection methods is a key concern nowadays as it reduces the interpretability of the predictive models as well as the trust of domain experts towards the selected feature subsets. Such experts would often prefer a more stable feature selection algorithm over an unstable and slightly more accurate one. In this paper, the compromise between feature relevance and selection stability is made explicit by biasing the selection towards some features through differential shrinkage of the Recursive Feature Elimination algorithm. Domain experts are given the opportunity to select any Paretooptimal trade-off of accuracy and selection stability based on their preferences. We propose the use of the hypervolume metric to assess the performance of methods realizing such a compromise. An associated confidence interval, based on the derivation of confidence regions of the accuracy-stability trade-off, is derived.
Results on prostate cancer diagnosis and handwritten recognition tasks show that the selection stability can be increased at will, often with a cost of classification performance. When some prior knowledge is available, far better compromises can be made. The design and evaluation of hybrid methods, learning the prior knowledge from the data, and using it to stabilize the selection is part of our future work.
Motivated by the needs of domain experts, across tasks feature stability is also studied in a transfer learning setting (i.e. when tasks are ordered). A biasing scheme that takes the stability measure explicitly into account is proposed. For similar, yet different, tasks, we show on microarray data that some bias is at first beneficial for both the accuracy and the stability. A too strong bias continues to increase the selection stability but at the cost of some classification performance, as the most relevant features vary across tasks. Our approach is evaluated here in a simulated transfer learning setting and further experimental validations will be conducted.
Different multi-task feature selection methods have been proposed in the
literature (e.g. Group-Lasso with `1/`p regularization [
The present paper answers the growing necessity of considering the selection stability not only as a side-effect of learning accurate predictive models but as an actual goal in a bi-objective framework. It proposes initial approaches to learn Pareto-optimal compromises in such a framework and, hopefully, opens the way to new works and improvements in this area.