Elastic SCAD as a novel penalization method for SVM classification tasks in high-dimensional data

Background Classification and variable selection play an important role in knowledge discovery in high-dimensional data. Although Support Vector Machine (SVM) algorithms are among the most powerful classification and prediction methods with a wide range of scientific applications, the SVM does not include automatic feature selection and therefore a number of feature selection procedures have been developed. Regularisation approaches extend SVM to a feature selection method in a flexible way using penalty functions like LASSO, SCAD and Elastic Net. We propose a novel penalty function for SVM classification tasks, Elastic SCAD, a combination of SCAD and ridge penalties which overcomes the limitations of each penalty alone. Since SVM models are extremely sensitive to the choice of tuning parameters, we adopted an interval search algorithm, which in comparison to a fixed grid search finds rapidly and more precisely a global optimal solution. Results Feature selection methods with combined penalties (Elastic Net and Elastic SCAD SVMs) are more robust to a change of the model complexity than methods using single penalties. Our simulation study showed that Elastic SCAD SVM outperformed LASSO (L1) and SCAD SVMs. Moreover, Elastic SCAD SVM provided sparser classifiers in terms of median number of features selected than Elastic Net SVM and often better predicted than Elastic Net in terms of misclassification error. Finally, we applied the penalization methods described above on four publicly available breast cancer data sets. Elastic SCAD SVM was the only method providing robust classifiers in sparse and non-sparse situations. Conclusions The proposed Elastic SCAD SVM algorithm provides the advantages of the SCAD penalty and at the same time avoids sparsity limitations for non-sparse data. We were first to demonstrate that the integration of the interval search algorithm and penalized SVM classification techniques provides fast solutions on the optimization of tuning parameters. The penalized SVM classification algorithms as well as fixed grid and interval search for finding appropriate tuning parameters were implemented in our freely available R package 'penalizedSVM'. We conclude that the Elastic SCAD SVM is a flexible and robust tool for classification and feature selection tasks for high-dimensional data such as microarray data sets.


Background
Classification and prediction methods play important roles in data analysis for a wide range of applications. Frequently, classification is performed on high-dimensional data, where the number of features is much larger compared to the number of samples ('large p small n' problem) [1]. In those cases, classification by Support Vector Machines (SVM), originally developed by Vapnik [2], is one of the most powerful techniques. The SVM classifier aims to separate the samples from different classes by a hyperplane with largest margin.
Often we do not only require a prediction rule but also need to identify relevant components of the classifier. Thus, it would be useful to combine feature selection methods with SVM classification. Feature selection methods aim at finding the features most relevant for prediction. In this context, the objective of feature selection is three-fold: (i) improving the prediction performance of the predictors, (ii) providing faster and more cost-effective predictors, and (iii) gaining a deeper insight into the underlying processes that generated the data.
Three main groups of feature selection methods exist: filter, wrapper and embedded methods [1,[3][4][5][6]. Filter methods simply rank individual features by independently assigning a score to each feature. These methods ignore redundancy and inevitably fail in situations where only a combination of features is predictive. Also, if there is a pre-set limit on the number of features to be chosen (e.g. top 10 features), this limit is arbitrary and may not include all informative features. Because of these drawbacks, the filter methods are not included in this work.
Connecting filtering with a prediction procedure, wrapper methods wrap feature selection around a particular learning algorithm. Thereby, prediction performance of a given learning method assesses only the usefulness of subsets of variables. After a subset with lowest prediction error is estimated, the final model with reduced number of features is built [5]. However, wrapper methods have the drawback of high computational load, making them less applicable when the dimensionality increases. Wrapper methods also share the arbitrariness of filter methods in feature selection.
The third group of feature selection procedures are embedded methods, which perform feature selection within learning classifiers to achieve better computational efficiency and better performance than wrapper methods. The embedded methods are less computationally expensive and less prone to overfitting than the wrappers [7].
Guyon [1] proposed the recursive feature elimination (RFE) method, which belongs to the wrapper methods. RFE iteratively keeps a subset of features which are ranked by their contribution to the classifier. This approach is computationally expensive and selecting features based only on their ranks may not derive acceptable prediction rules.
An alternative to SVM with RFE is to use penalized SVM with appropriate penalty functions. Penalized SVM belongs to embedded methods and provides an automatic feature selection. The investigation of the widely used family of penalization functions such as LASSO, SCAD, Elastic Net [8][9][10] and a novel proposed penalty Elastic SCAD in combination with SVM classification, is the objective of the paper. The ridge penalty [4] corresponds to the ordinary SVM, which does not provide any feature selection, is used as reference with respect to prediction accuracy.
Although feature selection methods can be applied to any high-dimensional data, we illustrate the use of these methods on microarray gene expression data due to their relevance in cancer research. Data from microarray experiments are usually stored as large matrices of expression levels of genes in rows and different experimental conditions in columns. Microarray technology allows to screen thousand of genes simultaneously. Detailed reviews on the technology and statistical methods often used in microarray analyses are presented in [11][12][13].
Since SVM is extremely sensitive to the choice of tuning parameters, the search for optimal parameters becomes an essential part of the classification algorithm [14]. The problem of choosing appropriate tuning parameters is discussed and an interval search technique from Froehlich and Zell [15] is proposed to use for SVM classification.
In this paper, we investigate the behaviour of feature selection SVM classifier techniques including commonly used penalization methods together with a novel penalization method, the Elastic SCAD. We compare them to SVM classification with and without recursive feature elimination (RFE [1]) for situations of 'large p small n' problems.
The RFE SVM is chosen as as a state-of-the-art representative of feature selection methods in applications [16,17].
A simulation study is designed to investigate the behaviour of different penalization approaches. Publicly available microarray data sets are chosen for illustration purposes as applications on real high-dimensional data.

Support Vector Machines
Suppose a training data set with input data vector x i ℝ p and corresponding class labels y i {-1, 1}, i = 1,..., n is given. The SVM finds a maximal margin hyperplane such that it maximises the distance between classes. A linear hyperplane can always perfectly separate n samples in n + 1 dimensions. Since we can assume that high-dimensional data with p ≫ n is generally linear separable [6], increasing complexity by using non-linear kernels is usually not needed. Thus, we use a linear SVM model throughout the paper.
The linear SVM separates classes by a linear boundary where w = (w 1 , w 2 ,..., w p ) is a unique vector of coefficients of the hyperplane with ||w|| 2 = 1 and b denotes the intercept of the hyperplane. We use '·' to denote the inner product operator. The class assignment for a test data vector x test R p is given by y test = sign [f (x test )].

Soft margin SVM
Soft margin SVM allows some data points to be on the wrong side of the margin. To account for erroneous decisions, slack variables ξ i ≥ 0, i = 1,..., n are defined as the distance between a misclassified data point and the corresponding margin. For data points on the correct side of the margin ξ i = 0, for data points inside the margin 0 <ξ i ≤ 1 and for misclassified data points ξ i > 1. The sum of non-zero ξ i is penalized with a cost parameter C and then added to the optimisation function penalty in the minimisation problem: The optimisation problem (2) is called the soft margin SVM. The cost parameter C is a data dependent tuning parameter that controls the balance between minimizing the coefficients of the hyperplane and correct classification of the training data set. C is often chosen by cross validation. Problem (2) can be solved by using convex optimisation techniques, namely by the method of Lagrange multipliers [4]. Convex optimisation techniques provide a unique solution for hyperplane parameters w and b where a i ≥ 0, i = 1,..., n are Lagrange multipliers. The data points with positive a i , are called support vectors (SVs). All data points lying on the correct side of their margin have a i = 0. Thus, they do not have any impact on the hyperplane, and we can rewrite Eq. (3) aŝ where the set of indices of the support vectors S is determined by S := {i : a i > 0}.
The coefficientb can be calculated from y i (ŵ · x i +b) = 1 − ξ i for any i with a i > 0. In praxis, an average of all solutions forb is used for numerical stability.

SVM as a penalization method
Hastie et al. [4] showed that the SVM optimisation problem is equivalent to a penalization problem which has the "loss and penalty" form where the loss term is described by a sum of the hinge loss functions l (y i , f (x i )) = [1 -y i f (x i )] + = max(1 -y i f (x i ), 0) for each sample vector x i , i = 1,..., n. The penalty term is denoted as pen l (w) and can have different forms:

Ridge penalty
The penalty term for ordinary SVM uses the L 2 norm: The L 2 penalty shrinks the coefficients to control their variance. However, the ridge penalty provides no shrinkage of the coefficients to zero and hence no feature selection is performed.

LASSO
The use of a L 1 penalization function is originally proposed by Tibshirani [8] for generalized linear models. The technique for parameter estimation with constraints is called LASSO ('least absolute shrinkage and selection operator'). Later, Bradley [18] adapted the L 1 -regularisation to SVM. Then, the penalty term has the form As a result of singularity of the L 1 penalty function, L 1 SVM automatically selects features by shrinking coefficients of the hyperplane to zero.
However, the L 1 norm penalty has two limitations. First, the number of selected features is bounded by the number of samples. Second, it tends to select only one feature from a group of correlated features and drops the others.
Fung and Mangasarian [19] have published a fast L 1 SVM modification, the Newton Linear Programming Support Vector Machine (NLPSVM), which we use in our analyses.

Smoothly clipped absolute deviation penalty (SCAD)
The SCAD penalty is a non-convex penalty function first proposed by Fan and Li [20]. Later, Zhang et al. [10] combined the SVM technique with the SCAD penalty for feature selection. The SCAD penalty function for a single coefficient w j is defined as where w j , j = 1,..., p are the coefficients defining the hyperplane and a > 2 and λ > 0 are tuning parameters. Fan and Li [21] showed that SCAD prediction is not sensitive to selection of the tuning parameter a. Their suggested value a = 3.7 is therefore used in our analyses.
The penalty term for SCAD SVM has the form The SCAD penalty corresponds to a quadratic spline function with knots λ at and aλ. For small coefficients w j , j = 1,..., p, SCAD yields the same behaviour as L 1 . For large coefficients, however, SCAD applies a constant penalty, in contrast to L 1 . This reduces the estimation bias. Furthermore, the SCAD penalty holds better theoretical properties than the L 1 penalty [21].

Elastic Net
To overcome the limitations of LASSO, Zou and Hastie [9] proposed a linear combination of L 1 and L 2 penalties which they called Elastic Net: The Elastic Net penalty provides automatic feature selection similar to L 1 , but is no longer bounded by the sample size. Moreover, at the same time this penalty manages to select highly correlated features (grouping effect). Increasing l 1 reduces the number of features of the classifier whereas for large l 2 one observes better control of the grouping effect. Wang [22] adapted the Elastic Net penalty to SVM classification problems. Therefore, the Elastic Net SVM optimisation problem can be written as where l 1 , l 2 ≥ 0 are the corresponding tuning parameters. [21] demonstrated the advantages of the SCAD penalty over the L 1 penalty. However, using the SCAD penalty might be too strict in selecting features for non-sparse data. A modification of the SCAD penalty analogously to Elastic Net could keep the advantages of the SCAD penalty, and, at the same time, avoid too restrictive sparsity limitations for non-sparse data.

Fan and Li
We therefore propose a combination of the SCAD and the L 2 penalties. The new penalty term has the form pen λ (w) := p j=1 p SCAD(λ 1 ) (w i ) + λ 2 ||w|| 2 2 , l 1 , l 2 ≥ 0 are the tuning parameters. We expect that the Elastic SCAD will improve the SCAD method for less sparse data. According to the nature of the SCAD and L 2 penalties, the Elastic SCAD should show good prediction accuracy for both, sparse and non-sparse data.
It can be shown that the combined penalty provides sparsity, continuity, and asymptotic normality when the tuning parameter for the ridge penalty converges to zero, i.e. l 2 0. The asymptotic normality and sparsity of Elastic SCAD leads to the oracle property in the sense of Fan and Li [21].
The Elastic SCAD SVM optimisation problem has the form where l 1 , l 2 ≥ 0 are the tuning parameters.

Elastic SCAD SVM: Algorithm
By solving Eq. (9) the same problems as for SCAD SVM occur: the hinge loss function is not differentiable at zero and the SCAD penalty is not convex in w. The Elastic SCAD SVM objective function can be locally approximated by a quadratic function and the minimisation problem can be solved iteratively similar to the SCAD approach [10,21].
For simplicity, we rename the SCAD penalty from p SCAD(λ 1 ) (|w j |) to p λ 1 (|w j |). Accordingly, the first-order derivative of the penalty is denoted by p λ 1 (·). Denote the penalized objective function in Eq. (9) by For each i (with respect to the fact that y 2 i = 1) the loss term can be split according to Given an initial value (b 0 ; w 0 ) close to the minimum of A(b, w), we consider the following local quadratic approximations: When w j0 is close to zero, setŵ j = 0; otherwise use the approximation for the SCAD penalty where due to symmetrical nature of the SCAD penalty |w j | is used instead of w j .
It can be shown that both approximations and their original functions have the same gradient at the point (b 0 , w 0 ). Therefore, the solution of the local quadratic function corresponds approximately to the solution of the original problem.
The local quadratic approximation of A(b, w) has the form By minimisation of A(b, w) with respect to w and b, terms without optimisation parameters w and b can be dropped (due to derivatives of constants): To write the equations in matrix form we define: Moreover, we define the matrix X = [1, x 1 ,..., x p ], where 1 is the vector of 1s with length n and x j is the jth input vector. Set Minimizing A (b, w) is then equivalent to minimizing the quadratic functioñ The solution to Eq. (10) satisfies the linear equation system The Elastic SCAD SVM can be implemented by the following iterative algorithm.
Step 4 Let k = k + 1. Go to step 2 until convergence.
If elements w j are close to zero, for instance, smaller than 10 -4 , then the jth variable is considered to be redundant and in the next step will be removed from the model. The algorithm stops after convergence of (b (k) , w (k) ).

Choosing tuning parameters
All SVM problems with or without feature selection use one or two tuning parameters which balance the tradeoff between data fit and model complexity. Since these parameters are data dependent, finding optimal tuning parameters is part of the classification task.

Fixed grid search
Tuning parameters are usually determined by a grid search. The grid search method calculates a target value, e.g. the misclassification rate, at each point over a fixed grid of parameter values. This method may offer some protection against local minima but it is not very efficient. The density of the grid plays a critical role in finding global optima. For very sparse grids, it is very likely to find local optimal points. By increasing the density of the grid, the computation cost increases rapidly with no guaranty of finding global optima. The major disadvantage of the fixed grid approach lies in the systematic check of the misclassification rates in each point of the grid. There is no possibility to skip redundant points or to add new ones.
When more parameters are included in the model, the computation complexity is increased. Thus, the fixed grid search is only suitable for tuning of very few parameters.

Interval search
Froehlich and Zell [15] suggested an efficient algorithm of finding a global optimum on the tuning parameter space using a method called EPSGO (Efficient Parameter Selection via Global Optimisation).
The main idea of the EPSGO algorithm is to treat the search for an optimal tuning parameter as a global optimisation problem. For that purpose, the Gaussian Process model is learned from the points in the parameter space which have been already visited. Thereby, training and testing of the GP is very efficient in comparison to the calculation of the original SVM models. New points in the parameter space are sampled by using the expected improvement criterion as described in the EGO algorithm [23], which avoids stacking in local minima. The stopping criteria of the EPSGO algorithm is either convergence of the algorithm or no change of the optimum during the last ten iterations.

Stratified cross validation
Using k-fold cross validation, the data set is randomly split into k disjoint parts of roughly equal size, usually k = 5 or k = 10. In addition, the data is often split in a way that each fold contains approximately the same distribution of class labels as the whole data set, denoted by stratified cross validation. For each subset, one fits the model using the other k -1 parts and calculates the prediction error of the selected kth part of the data.
The case k = n is called leave one out cross validation (LOO CV). The choice of k determines a trade-off between bias and variance of the prediction error. Kohavi [24] showed that ten-fold stratified cross validation showed better performance in terms of bias and variance compared to 10 <k <n. Hastie et al. [4] recommended to perform five-or ten-fold cross validation as a good compromise between variance and bias. We used both, five-and ten-fold stratified cross validation for simulation study and real applications, respectively.
In the next two sections the application of penalized SVM classification methods are compared. We used simulated and publicly available data to investigate the behaviour of different feature selection SVMs. For all comparisons the R pack-ages "penalizedSVM" [25] and "e1071" [26] were used which are freely available from the CRAN http://cran.r-project.org/, R version 2.10.1. The R package "e1071" is a wrapper for the well-known LIBSVM software [27]. We used five-and ten-fold stratified cross validation in combination with interval search for tuning parameters as described above.

Results and Discussion
Simulation study Simulation design A comprehensive simulation study evaluating the performance of four feature selection SVM classifiers, L 1 SVM, SCAD SVM, Elastic Net SVM and Elastic SCAD SVM, was performed. We used the ordinary L 2 SVM algorithm with a liner kernel as a reference for prediction accuracy.
Two independent data sets are simulated: a training set for building the classifier and a test set for estimating of the prediction errors of classifiers. First, the training data is generated, and the optimal tuning parameters are found using five-fold stratified cross validation according to the interval search approach [15]. Then, the classification hyperplane is computed using the estimated tuning parameters. Finally, application of the classification rule to the test data provides the prediction characteristics such as misclassification error, sensitivity and specificity.
Training and test input data are represented by a data matrix X = {x i }, i = 1,..., n, where x i ℝ p describes feature patterns for the ith sample. The input data X follows a multivariate normal distribution with mean μ and covariance matrix Σ. The class labels Y = {Y i }, i = 1,..., n are generated according to a logistic regression model where b = {b 1 ,..., b p } is a vector of coefficients of a classifier and u i are realisations of a variable following a U 0[1] distribution.
In our simulations the percentage of relevant features varies between 1% and 20%. Coefficients b j , j = 1,..., p are always defined as with equal numbers of positive and negative coefficients. The intersect b 0 is set to zero.
We also consider to have 'clumps' of correlated features. The clumpy dependency is supposed to describe the most common type of dependency in microarray studies [28]. We define "clumps" of correlated features as blocks of one relevant and four redundant features with a covariance matrix Σ* (k) , where k is the number of the current block. The diagonal elements of Σ* (k) for each block are equal to one and the off-diagonal elements are equal to r = 0.8. In total, we design five blocks of correlated features and therefore the covariance matrix has the form , where r denotes the number of relevant features. Using correlated blocks we investigate the ability of selecting correlated features, the so called grouping effect.
Optimal tuning parameters are found by an interval search in tuning parameter space using five fold cross validation. We select a large tuning parameter interval to be certain not to stick in local optima. The tuning parameter space for L 1 and SCAD SVM is one-dimensional with l 1 [l 1,min , l 1,max ]. Elastic SCAD has two tuning parameters l 1 , l 2 [l l, min , l l, max ], l = 1, 2. Elastic Net applies LARS paths. for fixed l 2 a l 1 path is calculated and the optimal l 1 is identified (for details refer to [17]). Thus, the optimal pair of parameter (λ * 1 , λ * 2 ) for Elastic Net was found in the two-dimensional space ℝ × [l l, min , l l, max ] We set the search interval for both parameters to [l l, min , l l, max ] = [2 -10 , 2 -10 ], l = 1, 2.
The performance of classifiers is characterised by the Youden index. The Youden index describes as equally weighted sum of true positive results ("sensitivity") and false positive results ("1 -specificity"): Youden index = sensitivity + specificity − 1 The maximal Youden index is one, when the true positive rate is one and the false positive rate is zero. For a random classifier the expected Youden index is zero. The sensitivity and specificity have equal weights in this index. Most often the costs and con-sequences of true positives and false positives will differ greatly. Therefore, Gu and Pepe [29] recommend reporting the two measures separately. For our simulated data, we consider the Youden index to be an appropriate indicator for feature selection methods performance, since we weight errors equally.
It is worth to mention, that for discrete classier the Youden index and the area under the curve (AUC) provide the same message due to their linear relationship. According to Greiner et al. [30], if there is only one point in the ROC plot, the ROC curve is estimated by connecting the three points. the point corresponds to the classifier, the (0, 0) and (1, 1) edges of the plot. Then geometrically, the estimated AUC corresponds to the average of estimated sensitivity and specificity. Thus, the Youden index and the AUC have a linear relationship. AUC = (sensitivity + specificity)/2 = (Youden index +1)/2. Optimizing the AUC will lead to the same results as optimizing the Youden index when dealing with discrete classifiers. Nevertheless, for real data application, the AUC values are presented in a separate column due to higher level of familiarity in bioinformatics.
Finally, the misclassification rate, size of the classifiers and frequencies of the selected features within 100 simulation runs are computed.

Simulation results
The performance of the feature selection methods applied to simulated data using p = 1000 features and n = 500 samples for training and testing is presented in the next section. The percentage of relevant features varies between 1% and 20% in four steps, i.e. r = 10, 50, 100, 200. We assume to have correlated blocks of features as described in the design section. The optimal tuning parameters were chosen as described above. Multiple comparisons in performance measures between the proposed prediction methods and the best method (the MCB test) for each simulation step will be done according to Hsu [31] based on 100 simulation runs. We used a non-inferiority margin of a procedure to distinguish methods with similar performance. Table 1 summarises the average misclassification rates depending on the number of relevant features. The numbers in parentheses are the standard errors of the estimates. For very sparse models (10 out of 1000 features are relevant) SCAD showed the lowest misclassification error (18%), followed by Elastic Net and Elastic SCAD (19.4% and 20.8% respectively), where both lie in indifference zone for best methods if the non-inferiority margin was set to Δ = 0.05. For less sparse to nonsparse models (r = 50 and r = 100) Elastic Net showed the best performance. For r = 200 relevant features L 1 and Elastic Net showed nearly the same results (32.9% and 33.1% respectively). The same was observed for SCAD (34.7%) and Elastic SCAD (34.2%). For r ≥ 50 the misclassification rate was indistinguishable for all feature selection methods with exception of the L 1 SVM. The L 2 SVM classifiers showed larger misclassification errors for sparse models (r = 10 and r = 50) than all other feature selection methods. For less sparse models differences in misclassification error levelled out.

Youden index
The average Youden index for very sparse models (r = 10) was considerably high for all feature selection methods: 0.96 for SCAD, 0.95 for Elastic Net, 0.92 for Elastic SCAD, and 0.81 for L 1 SVM (Table 2). By increasing number of informative features, the Elastic Net SVM showed the best Youden index (0.71% -0.27%) among all feature selection methods, closely followed by the Elastic SCAD SVM (0.67% -0.27%), both being indistinguishable.
All methods except the L1 SVM provided significantly comparable Youden indexes at the level a = 0.05 and a relevant difference Δ = 0.10 for r = 10. By increasing model complexity, the Elastic Net SVM showed the best Youden Index among all feature selection methods, closely followed by the Elastic SCAD SVM. Starting from r > 100 the is no significant difference between Elastic Net and Elastic SCAD SVMs. With increasing number of relevant features, the Youden index decreases from 0.9 to 0.27 for 'elastic' methods to 0.14 for the L 1 SVM and to 0.16 for the SCAD SVM. respectively.

Sparsity of the classifier
The SCAD SVM provided the most sparse classifier (in terms of selecting the smallest number of features) for r = 10 and r = 50 out of 1000 features (cf. Table 3). It selected 12 and 61 features, respectively. For less sparse models the Elastic Net and the Elastic SCAD SVMs had similar performance, selecting the smallest number of features.

Selection Frequencies
A frequencies plot for the simulation study is represented in 'Additional file 1 -Frequencies plot'. With increasing number of relevant features (r), a decrease of the proportion of true positives (in red) and an increase of the proportion of false positives (in blue) for all feature selection models was observed, respectively. At the same time we observed an increase of the false positives, which are correlated with the true positives (in green) in classifiers.
The percentage of true positives in the classifiers is shown in Table S1 (Additional file 2 -Tables S1, S2, S3). For r = 10 relevant features the Elastic Net SVM found almost all true positives (99.8%), followed by the Elastic SCAD SVM with 97.6%. For r = 50 the Elastic SCAD SVM achieved the sparsest solution followed by the L 1 SVM. In less sparse models, the L 1 SVM showed highest true positive rates of 84.5% and 86%.

Grouping effect
We further evaluated the ability of feature selection methods to select correlated features of true positives. Although for all scenarios L 1 SVM has found the largest percentage of correlated features, which increases with increasing number of relevant features (23.6 -62.5%), the level of correlated features is comparable to the level of non-relevant features (Table S2).
Comparing Tables S1, S2 and S3 one can observe that the SCAD and the L 1 SVMs failed to find features highly correlated with true positives more often than with independent false positives. The Elastic Net and the Elastic SCAD SVMs managed to discover correlated features (in green) more often than the independent false positives (in blue), at least for sparse models (r = 10 and r = 50).

The false positive rate
For very sparse models, the false positive rate (FPR) was the smallest for the SCAD SVM, followed by the Elastic Net and the Elastic SCAD SVMs (Table S3). For other less sparse models the Elastic Net SVM selected fewer false positives than the remaining methods. The second best method is the Elastic SCAD SVM.

Conclusions
• As expected from theory the SCAD SVM and the L 1 SVM produced classifiers with low prediction error for very sparse situations.
• For less sparse and non-sparse models, the Elastic Net and the Elastic SCAD SVM showed better results than the L 1 and the L 2 SVMs with respect to accuracy, Youden index and sparsity of classifiers.
• The SCAD SVM and the L 1 SVM were not able to find correlated features. The Elastic Net and the Elastic SCAD SVMs found correlated features more frequently than one would expected under random selection. Although the grouping effect strength Training and test data with 1000 features and 500 samples were simulated. The number of relative features (r) were increased from r = 10 to r = 200 in four steps. Each simulation step was based on 100 simulations of training and test data. In bold -the significant best method(s) according to the MCB test at the family-wise significance level a = 0.05 and non-inferiority margin of Δ = 5%. In bold -the significant best method(s) according to the MCB test at the family-wise significance level a = 0.05 and non-inferiority margin of Δ = 0.10. In bold -median number of features that come closest to the true number of relevant features per simulation scenario, (in parentheses -median absolute deviation); underline -the second best.
weakens with increasing number of relevant features, the Elastic Net and Elastic SCAD SVMs still managed the grouping effects.
• In general, the Elastic Net and the Elastic SCAD SVMs showed similar performance. Additionally, the Elastic SCAD SVM provided more sparse classifiers than the Elastic Net SVM.

NKI breast cancer data set
Two studies on breast cancer from the Netherlands Cancer Institute (NKI) were published by the van't Veer group [32], [33]. In the first paper, a set of 78 lymph node negative patients with pre-selected 4919 clones were used to find a predictor for distant metastases. The classifier was trained and validated on patients who developed distant metastases within five years after surgery and patients being metastasis-free for at least the first five years. The resulting predictor was a 70-gene signature also known as MammaPrint(R). We will use the MammaPrint(R) signature as reference in the analysis of the NKI breast cancer data set. The signature was generated based on gene-wise correlations between the gene expression and metastasis occurrence. The data set was taken from http://www.rii. com/publications/2002/vantveer.html. In a subsequent validation study, data from 295 patients (which partially included patients from the first study) were used to validate the signature [33]. Among the patients, 151 were lymph node negative and 144 had lymph node positive disease. The pre-processed data containing 4919 clones is available from http://www.rii. com/publications/2002/nejm.html.
After excluding patients being identical to the training set and 10 patients with no metastasis information, 253 patients remained. Among the 253 patients there are 114 lymph node negative and 139 lymph node positive patients.
In this paper, we combined the 78 lymph node negative sample set from the first publication with 114 lymph node negative patients from the validation study. In total, a data set with 192 lymph node negative samples was used. The estimation of classifier performance was computed by a ten-fold stratified crossvalidation. Table 4 shows the misclassification error, sensitivity, specificity, Youden index and AUC value of four feature selection methods, RFE SVM and standard L 2 SVM based on ten-fold stratified cross validation.

Results on NKI breast cancer data set
RFE SVM was used according to Guyon's approach [1], where at each iteration half of features with lowest ranks are eliminated. To increase the classifier's stability, RFE SVM with five-fold stratified cross validation was repeated five times. According to the average cross validation error the optimal number of features was 2 8 = 256. Optimal tuning parameters for penalized SVM methods were found by the interval search on the tuning parameter space as described in the method section using ten-fold stratified cross validation.
The SCAD SVM reduced the number of features from 4919 to 476, L 1 SVM selected 1573 features, Elastic Net 109 features, and the Elastic SCAD had 459 features in the classifier. For the NKI data set the best predictor with respect to misclassification error was L 1 SVM. Elastic Net and Elastic SCAD SVMs provided similar results, followed by SCAD SVM, which was slightly worse.
The relationship between the true positive rate (TPR, sensitivity) and the false positive rate (FPR, 1-specificity) for each classifier is depicted as a point in the ROC plot ( Figure 1). Isolines with constant Youden index are plotted as dashed lines. Taking the Youden index as an additional criterion, one could prioritise L 1 SVM. RFE SVM and both 'elastic' methods lay clustered in the ROC plot with clear distance to the L 1 classifier. The L 2 was placed in-between L 1 and this cluster, being not far from the cluster.
Interestingly, the MammaPrint(R) signature ("70_sign") neither showed good test accuracy nor a reliable sensitivity or specificity. L 2 SVM and the feature selection methods outperformed the published signature.

Conclusions
For the two data sets from van't Veer group feature selection methods produced signatures with similar prediction accuracy, but being different in size. L 1 SVM with a non-sparse classifier provided the best sensitivity and specificity, followed by more sparse predictors from Elastic Net SVM and Elastic SCAD SVM.

MAQC-II breast cancer data set
This data set is part of the MicroArray Quality Control (MAQC)-II project, which has been designed to investigate numerous data analysis methods and to reach consensus on the "best practices" for development and validation of microarray-based classifiers for clinical and preclinical applications. One biological endpoint is estrogen receptor (ER) status. Out of 230 patients in total, 89 patients have negative ER status and 141 patients positive ER status. A clinical endpoint is pathological complete response (pCR) to preoperative chemotherapy. Among the 230 patients in the data set, 182 patients showed no pCR and 48 had a pCR.
The preprocessed data contains 22283 features and is available from GEO database, accession number GSE20194.

Results on MAQC-II breast cancer data set
The feature selection methods SCAD SVM, L 1 SVM, Elastic Net SVM and Elastic SCAD SVM with internal tenfold stratified cross validation were applied to build classifiers. Additionally, the L 2 SVM and the RFE SVM were used as reference models. To achieve performance measurements ten-fold stratified cross validation was used.

pCR prediction
Based on the minimal average misclassification error, the optimal number of features of RFE SVM classifier was obtained to be 2 11 = 2048 ( Table 5). The penalized SVM methods provided moderately sparse models, Elastic SCAD SVM with 148 features, Elastic Net SVM with 398 features and dense models, L 1 , SCAD and RFE SVMs with more than 1000 features.
The misclassification error rate was similar for all methods with the Elastic SCAD classifier showing the lowest error rate of 15%. With nearly equally high specificity (91-94%), we observed large variations in sensitivity of different feature selection methods as shown in the corresponding ROC plot (Figure 2). The Elastic SCAD outperformed all methods with sensitivity of 52%. Interestingly, the Elastic Net showed the smallest sensitivity of 15% resulting in a small Youden index of 0.06.
Overall, Elastic SCAD showed better classification characteristics than other methods. Moreover, the higher specificity of the Elastic SCAD classifier is of clinical importance. The patients that did not respond to the therapy were recognized with higher probability.

ER status
We also used the MAQC-II data set to predict the ER status. Here, the L 1 SVM failed to derive a sparse solution, whereas SCAD, Elastic Net and Elastic SCAD SVM classifiers were similar (Table 6). Moreover, Elastic SCAD showed the smallest error rate and highest sensitivity over all methods.
All classification methods provided small misclassification errors, high sensitivity and high specificity. The ROC plot in Figure 3 demonstrates this performance of  predictors. As presented in Table 6 the Elastic Net, SCAD and Elastic SCAD SVMs selected small numbers of features, 3, 32 and 59 out of 22283, respectively. The extreme sparseness of the Elastic Net SVM was paid by lower sensitivity and specificity compared to other methods. The misclassification test error was similar for all methods (7-14%). The Elastic SCAD SVM classifier showed the smallest error rate of 7%. For this classification task, the sparse classifier Elastic SCAD and SCAD showed the best characteristics.

Screening on two additional breast cancer data sets
These data sets were recently analysed and published by Johannes et. al. [34]. The first data set, the Mainz cohort, contains of 154 lymph node-negative, relapse free patients and 46 lymph node-negative patients that suffered a relapse (GEO acession number GSE11121). The relapse is defined as appearance of distant metastasis within five years after the treatment. The second data set, the Rotterdam cohort, represents 286 lymph node-negative breast cancer samples including 107 re-lapse events (GSE2034). Both data sets were generated using the Affymetrix HG-U133A platform, normalized with the same methods and relapse as the primary classification endpoint. We trained the feature selection classifiers on the whole cohort, Mainz data or Rotterdam data, and used the other cohort as an independent validation data set, respectively as presented in Tables 7 and 8.
We can see that all feature selection methods had lower misclassification test error than the L 2 SVM containing all features for breast cancer data sets. The classifiers perform different for each data set. The Elastic Net SVM had small error rate for the Rotterdam cohort, but failed to classify the Mainz samples adequately. The L 2 SVM classifier including all features had the second best Youden index for the Mainz set, however for Rotterdam data showed the worst Youden index. Using both, the test error and AUC value as a combined measure of sensitivity and the specificity, one would conclude that the L 1 , SCAD and Elastic SCAD SVMs provide reasonable and robust solutions with respect to the combined analysis of the two breast cancer data sets.
Altogether, Elastic SCAD seems to provide an overall acceptable compromise for sparse and non-sparse data.

Conclusions
In high-dimensional prediction tasks, feature selection plays an important role. In this paper, we proposed a novel feature selection method for SVM classification using a combination of two penalties, SCAD and L 2 . The commonly used penalty functions L 1 , SCAD and Elastic Net were investigated in parallel with the new method on simulated and public data. To address the problem of finding optimal tuning parameters for SVM classification the efficient parameter search algorithm from Froehlich and Zell [15] was implemented.
In almost all cases, the four feature selection classifies outperformed ordinary Support Vector Classification using the L 2 penalty. From the simulation study we  Figure 2 ROC plot for MAQC-II breast data set with pCR as endpoint. The characteristics for the different feature selection methods were derived using ten-fold statrifierd cross validation. TPR and FPR values are presented as points (x axis: 1-specificity = FPR, y axis. sensitivity = TPR). RFE_256 is RFE SVM with 1024 top ranked features, ENet is Elastic Net SVM, ESCAD is Elastic SCAD SVM. Gray dashed lines depict isolines of the Youden index. concluded that for sufficiently large sample sizes, feature selection methods with combined penalties are more robust to changes of the model complexity than using single penalties alone. The SCAD SVM followed by the L 1 SVM, as expected, showed very good performance in terms of pre-diction accuracy for very sparse models, but failed for less sparse models. Combined penalty functions in combination with the SVM algorithm, Elastic Net and Elastic SCAD, performed well for sparse and less sparse models.
Comparisons with commonly used penalty functions in the simulation study illustrated that the Elastic SCAD and the Elastic Net SVMs showed similar performance with respect to prediction accuracy. Both 'elastic' methods were able to consider correlation structures in the input data (grouping effect). However, the Elastic SCAD SVM in general provides more sparse classifiers than the Elastic Net SVM.
Finally, applied to publicly available breast cancer data sets, the Elastic SCAD SVM performed very flexible and robust in sparse and non-sparse situations. Results from the simulation study and real data application render Elastic SCAD SVM with automatic feature selection a promising classification method for high-dimensional applications.  Figure 3 ROC plot for MAQC-II breast data set with ER as endpoint. The characteristics for the different feature selection methods were derived using ten-fold stratified cross validation. TPR and FPR values are presented as points (x axis: 1-specificity = FPR, y axis. sensitivity = TPR). RFE_256 is RFE SVM with 1024 top ranked features, ENet is Elastic Net SVM, ESCAD is Elastic SCAD SVM. Gray dashed lines depict isolines of the Youden index.