Incremental wrapper-based gene selection from microarray data for cancer classification

Incremental wrapper-based gene selection from microarray data for cancer classification Roberto Ruiza,∗, José C. Riquelmea, Jesús S. Aguilar-Ruizb aDepartment of Computer Science, University of Seville, Avda. Reina Mercedes s/n. 41012 Seville, Spain bPolytechnic, Pablo de Olavide University, Ctra. Utrera, km 1, 41013 Seville, Spain Abstract Gene expression microarray is a rapidly maturing technology that provides the opportunity to assay the expression levels of thousands or tens of thousands of genes in a single experiment. We present a new heuristic to select relevant gene subsets in order to further use them for the classification task. Our method is based on the statistical significance of adding a gene from a ranked-list to the final subset. The efficiency and effectiveness of our technique is demonstrated through extensive comparisons with other representative heuristics. Our approach shows an excellent performance, not only at identifying relevant genes, but also with respect to the computational cost. Keywords: Microarray; Gene selection; Classification; Feature selection 1. Introduction All cells have a nucleus, and inside nucleus there is DNA, which encodes the “program” for making future organisms. DNA has coding and non-coding segments, and coding segments, called “genes”, specify the structure of proteins, which are large molecules, like hemoglobin, that do the essential work in every organism. Practically all cells in the same organism have the same genes, but these genes can be expressed differently at different times and under different conditions. Genes make proteins in two steps. First, DNA is transcribed into messenger RNA or mRNA, which in turn is translated into proteins (see [1] for a review). In recent years there has been an explosion in the rate of acquisition of biomedical data. Advances in molecular genetics technologies, such as DNA microarrays allow us for the first time to obtain a “global” view of the cell. Analysis of microarray data presents unprecedented opportunities and challenges for data mining in areas such as gene clustering [2], sample clustering and class discovery [2,3], sample classification [3,4] and gene selection [5–9]. ∗Corresponding author. Tel.: +34 954553867; fax: +34 954557139. E-mail address: [email protected] (R. Ruiz). In this work, we address the gene selection issue under a classification framework. The aim is to build a classifier that accurately predicts the classes (diseases or phenotypes) of new unlabeled samples. Three well-known machine learning classifiers (naïve Bayes, instance-based and decision trees), with completely different approaches to learning, are applied to perform the class prediction. A typical data set may contain thousands of genes but only a small number of samples (often less than 200). Theoretically, having more genes should give us more discriminating power. However, this can cause several problems: increase computational complexity and cost; too many redundant or irrelevant genes; and estimation degradation in the classification error. In addition to reducing noise and improving the accuracy of classification, the selected subsets of genes may have important biological interpretation and may be used for drug target discovery or identifying future possible research directions. The problem of feature selection received a thorough treatment in pattern recognition and machine learning. Most of the feature selection algorithms approach the task as a search problem, where each state in the search specifies a distinct subset of the possible attributes [10]. The search procedure is combined with a criterion in order to evaluate the merit of each candidate subset of attributes. There are a lot of possible combinations between each procedure search and each attribute measure [11]. However, search methods can be prohibitively expensive in high-dimensional data sets, especially when a data mining algorithm is applied as evaluation function. There are various ways in which feature selection algorithms can be grouped according to the attribute evaluation measure: depending on the type (filter or wrapper techniques) or on the way the features are evaluated (individual or subset evaluation). The filter model relies on general characteristics of the data to evaluate and select gene subsets without involving any mining algorithm. The wrapper model requires one predetermined mining algorithm and uses its performance as the evaluation criterion. It searches for features better suited to the mining algorithm, aiming to improve mining performance, but it also is more computationally expensive [12,13] than filter models. Feature ranking (FR), also called feature weighting [10,14], assesses individual features and assigns them weights according to their degrees of relevance, while the feature subset selection (FSS) evaluates the goodness of each found feature subset. (Unusually, some search strategies in combination with subset evaluation can provide a ranked list). In order to compare the effectiveness of gene selection, gene sets chosen by each technique are tested with three well-known learning algorithms: a probabilistic learner (naïve Bayes), an instance-based learner (IB1) and a decision tree learner (C4.5). These three algorithms have been chosen because they represent three quite different approaches to learning, and their long standing tradition in classification studies. The comparison is performed in four DNA microarray data sets involved in diagnosis of cancer (colon, leukemia, lymphoma and global cancer map (GCM)). The paper is organized as follows. In the next two sections, we will review previous work, and notions of feature relevance and redundancy, respectively. In Section 4, we will present our proposed measures of gene relevance and redundancy using a wrapper approach, and our algorithm is described. Experimental results are shown in Section 5, and the most interesting conclusions are summarized in Section 6. 2. Related work Traditional gene selection methods often select the topranked genes according to their individual discriminative power [3]. This approach is efficient for high-dimensional data due to its linear time complexity in terms of dimensionality. They can only capture the relevance of genes to the target concept, but cannot discover redundancy and basic interactions among genes. In the FSS algorithms category, candidate gene subsets are generated based on a certain search strategy. Different algorithms address these issues distinctively. In [11], a great number of selection methods are categorized. We found different search strategies, namely exhaustive,heuristic and random search, combined with several types of measures to form different algorithms. The time complexity is exponential in terms of data dimensionality for exhaustive search and quadratic for heuristic search. The complexity can be linear to the number of iterations in a random search, but experiments show that in order to find the best feature subset, the number of iterations required is usually at least quadratic to the number of features [15]. The most popular search methods in pattern recognition and machine learning cannot be applied to these gene expression data sets due to the large number of genes (sometimes tens of thousands). One of the few used search techniques in these domains is sequential forward (SF, also called hill-climbing or greedy search). Different subset evaluation measures in combination with SF search engine can be found. We are specially interested in wrapper approach [6,16]. A key issue of wrapper method is how to search into the space of subsets of genes. Although several heuristic search strategies such as greedy sequential search, best-first search, and genetic algorithm exist, most of them are still computationally expensive O(N2)(Nbeing the number of genes of the original data set), which prevents them from scaling well to data sets containing thousands of genes. A rough estimate of the time required by most of these techniques is in the order of thousands of hours, assuming that the method does not get caught in a local minima first and stops prematurely. For example, if we have chosen 50 genes from 20 000 (0.0025% of the whole set) through a greedy search, the subset evaluator would be run approximately one million times (Ntimes to find the best single gene, then tries each of the remaining genes in conjunction with the best to find the most suited pair of genes N−1 times, and so on, more or less 20 000 ×50). Assuming 4 s on average by each evaluation, the results would take more than 1000 h. The limitations of both approaches, FR and FSS, clearly suggest that we should pursue a hybrid model. Recently, a new framework of feature (gene) selection has been used, where several above-mentioned approaches are combined. Yu and Liu [17] proposed a fast correlation-based filter algorithm (FCBF) which used correlation measure to obtain relevant genes and to remove redundancy. There are other methods based on relevance and redundancy concepts. Recursive feature elimination (RFE) is a proposed feature selection algorithm described by Guyon et al. [8]. The method, given that one wishes to find only rdimensions in the final subset, works by trying to choose the rfeatures which lead to the largest margin of class separation, using an SVM classifier. This combinatorial problem is solved in a greedy fashion at each iteration of training by removing the input dimension that decreases the margin the least until only r input dimensions remain (this is known as backward selection). Ding and Peng [5] have used mutual information for gene selection that has maximum relevance with minimal redundancy by solving a simple two-objective optimization. Xing et al. [9] proposed a hybrid of filter and wrapper approaches to feature selection. In [18], a rank search method is proposed to compare feature selection algorithms. Rank search techniques rank all genes, and subsets of increasing size are evaluated from the ranked list (i.e., the first attribute, the two first ones, etc.). The best attribute set is reported. The authors apply the wrapper approach to data sets up to 300 attributes and state that for the ADS data set (1500 attributes) the estimated time to only generate the ranking in a machine with 1.4 GHz processor would be about 140 days and to evaluate the ranked list of attributes would take about 40 days. In contrast, our method can be tested on data sets with 20 000 genes in a similar machine in a few hours. This paper presents a new gene selection method, based on the hybrid model, and attempts to take advantage of all of the different approaches by exploiting their best performances in two steps: first, a filter or wrapper approach provides a ranked list of genes, and second, ordered genes are added using a wrapper subset evaluation ensuring good performance (the search algorithm is valid for any gene rankedlist). This approach provides the possibility of efficiently applying wrapper model in high-dimensional domains, obtaining better results than the filter model. The final subset is obviously not the optimum, but it is unfeasible to search for every possible subset of genes through the search space. The main goal of our research is to obtain a few features with high predictive power. 3. Preliminary concepts 3.1. Relevance The purpose of a feature subset algorithm is to identify relevant features according to a definition of relevance. However, the notion of relevance in machine learning has not yet been rigorously defined in common agreement [19].Kohavi and John [13] include three disjointed categories of feature relevance: strong relevance, weak relevance and irrelevance. These groups are important to decide what features should be conserved and which ones can be eliminated. The strongly relevant features are, in theory, important to maintain a structure in the domain, and they should be conserved by any feature selection algorithm in order to avoid the addition of ambiguity to the sample. Weakly relevant features could be important or not, depending on the other features already selected and on the evaluation measure that has been chosen (accuracy, simplicity, consistency, etc.). Irrelevant attributes are not necessary at all. Bell and Wang [19] make use of information theory concepts to define the entropic or variable relevance of a feature with respect to the class. Blum and Langley [10] collect several relevance definitions. The above notions of relevance are independent of the specific learning algorithm being used. There is no guarantee that just because a feature is relevant, it will necessarily be useful to an algorithm (or vice versa). The definition of incremental relevance by Caruana and Freitag [20] makes it explicit, being considered especially suited to obtain a predictive feature subset. Definition 1 (Incremental usefulness). Given a sample of data D, a learning algorithm L, and a feature subset F, the feature xiis incrementally useful to Lwith respect to Fif the accuracy of the hypothesis that Lproduces using the group of features {xi}∪Fis better than the accuracy achieved using just the subset of features F. We consider this definition to be especially suited to obtain a predictive feature subset. In the next section, concepts can be applied to avoid a subset which contains attributes with the same information. 3.2. Redundancy Notions of feature redundancy are normally in terms of feature correlation. It is widely accepted that two features are redundant to each other if their values are completely correlated. There are two widely used types of measures for the correlation between two variables: linear and nonlinear. In the first case, the Pearson correlation coefficient is used, and in the second one, many measures are based on the concept of entropy, or measure of the uncertainty of a random variable. Symmetrical uncertainty is frequently used, defined as SU(X, Y ) =2IG(X|Y) H(X)+H(Y), where H(X)=− iP(x i)log2(P (xi)) is the entropy of a variable Xand IG(X|Y)=H(X)−H(X|Y)is the information gain from Xprovided by Y. The above-mentioned definitions are between pairs of variables. However, it may not be as straightforward in determining feature redundancy when one is correlated with a set of features. Koller and Sahami [21] apply a technique based on cross-entropy, named Markov blanket filtering, to eliminate redundant features. This idea is formalized in the following definition. Definition 2 (Markov blanket). Given a feature xi∈F(a set of attributes) and the class C, the subset M⊆F(x i/∈M) is a Markov blanket of xiif, given M,xiis conditionally independent of F−M−{xi}and C. Two attributes (or set of attributes) X, Y are said to be conditionally independent given a third attribute Z(or set) if, given Zmakes Xand Yindependent, i.e., the distribution of X, known Yand Z, is equal to the distribution Xknown Z, therefore, Ydoes not have influence on X(P(X|Y, Z) = P(X|Z)). Theoretically, it can be shown that once we find a Markov blanket Mof feature xiin a feature set F, we can safely remove xifrom Fwithout increasing the divergence from the original distribution. Furthermore, in a sequential filtering process, in which unnecessary features are removed one by one, a feature tagged as unnecessary based on the existence of a Markov blanket Mremains unnecessary in later stages when more features have been removed. The Markov blanket condition requires that Massumes not only the information that xihas about C, but also about all the other features. In [21] it is stated that the cardinality of set Mmust be small and fixed. Xing et al. [9] and Yu and Liu [17] are among the most cited works at present following the above-mentioned framework (FR +FSS). Both of them are based on this concept of Markov blanket. In the first one, the number of attributes of Mare not provided, but it is a fixed number among the highly correlated features. In the second one, a FCBF is implemented, where Mis formed by only one attribute, and gradually eliminates redundant attributes with respect to M from the first to the final attribute of an ordered list. Other methods based on relevance and redundancy concepts can be found in [5,8]. 4. Incremental performance of wrapper-search over ranking In this section, we first introduce our ideas of relevance and redundancy taking into account the aim of applying a wrapper model to gene expression data sets, and then our approach is described. As previously indicated, the wrapper model makes use of the algorithm that will build the final classifier to select a gene subset. Thus, given a classifier L, and given a set of genes G, a wrapper method searches in the space of G, using cross-validation to compare the performance of the trained classifier Lon each tested subset. While the wrapper model is more computationally expensive than the filter model, it also tends to find gene sets better suited to the inductive biases of the learning algorithm and therefore provides superior performance. In this work, we propose a fast search over a minimal part of the gene space. Beginning with the first gene from the list ordered by some evaluation criterion, genes are added one by one to the subset of selected genes only if such inclusion improves the classifier accuracy. Then, the learning algorithm of the wrapper approach is always run N(number of genes) times, usually with a few genes. A gene ranking algorithm makes use of a scoring function computed from the values of each gene and the class label. By convention, we assume that a high score is indicative of a valuable gene and that we sort genes in decreasing order of this score. We consider ranking criteria defined for individual genes, independently of the context of others. When a ranking of genes is provided from a highdimensional data set, a large number of genes with similar scores is generated, and a common criticism is that it leads to the selection of redundant subsets. However, according to Guyon and Elisseeff [14], noise reduction and consequently, better class separation may be obtained by adding variables that are presumably redundant. Moreover, a very high attribute correlation (in absolute value) does not mean absence of attribute complementarity. Therefore, our idea of redundancy is not based on correlation measures, but on the learning algorithm target (wrapper approach), in the sense that a gene is chosen if additional information is gained by adding it to the selected subset of genes. 4.1. Incremental ranked usefulness In gene subset selection, it is a fact that two types of genes are generally perceived as being unnecessary: genes that are irrelevant to the target concept, and genes that are redundant given other genes. Our approach is based on the concept of Markov blanket, which is described in [21]. This idea was formalized using the notion of conditionally independent attribute, which can be defined by several approaches [9,17]. We set this concept by a wrapper model, defining incremental ranked usefulness in order to devise an approach to explicitly identify relevant genes and do not take into account redundant genes. Let Dbe a sample of labeled data; Gbe a subset of genes of D; and Lbe a learning algorithm; correct rate (or accuracy) (D/G, L) is the ratio between the number of instances correctly classified by Land the total number of instances of Dconsidering only the subset G. Let R={gi},i=1...N be a ranking of all the genes in Dsorted in descending order, and Gbe named the subset of the ifirst genes of R. Definition 3 (Incremental ranked usefulness). The gene gi+1in Ris incrementally useful to Lif it is not conditionally independent of the class Cgiven G, therefore the correct rate of the hypothesis that Lproduces using the group of genes {gi+1}∪Gis significantly better (denoted by 5) than the correct rate achieved using just the subset of genes G. Therefore, if (D/G ∪{gi+1},L)5/(D/G, L), then gi+1is conditionally independent of class Cgiven the subset G, then we should be able to omit gi+1without compromising the accuracy of class prediction. A fundamental question in the previous definition is how the significant improvement is analyzed. A five-fold crossvalidation is used to estimate if the accuracy of the learning scheme for a set of genes is significantly better (5)than the accuracy obtained for another set. We conducted a Student’s paired two-tailed t-test in order to evaluate the statistical significance (at 0.1 level) of the difference between the previous best subset and the candidate subset. This last definition allows us to select genes from the ranking, but only those that increase the classification rate significantly. Although the size of the sample is small (five-folds), our search Fig. 1. BIRS Algorithm. method uses a t-test. We want to obtain a heuristic, not to do an accurate population study. However, on the one hand it must be noticed that it is a heuristic based on an objective criterion, to determine the statistical significance degree of difference between the accuracies of each subset. On the other hand, the confidence level has been relaxed from 0.05 to 0.1 due to the small size of the sample. Statistically significant differences at the p<0.05 significance level would not allow us to add more features, because it would be difficult for the test to obtain significant differences between the accuracy of each subset. Obviously, if the confidence level is increased, more genes can be selected, and vice versa. 4.2. Algorithm There are two phases in the algorithm, named BIRS (best incremental ranked subset), shown in Fig. 1: firstly, the genes are ranked according to some evaluation measure (lines 1–4). In the second phase, we deal with the list of genes once, crossing the ranking from the beginning to the last ranked gene (lines 5–12). We obtain the classification accuracy with the first gene in the list (line 9) and it is marked as selected (lines 10–12). We obtain the classification rate again with the first and second genes. The second will be marked as selected depending on whether the accuracy obtained is significantly better (line 10). Repeat the process until the last gene on the ranked list is reached. Finally, the algorithm returns the best subset found, and we can state that it will not contain irrelevant or redundant genes. The first part of the above algorithm is efficient since it requires only the computation of Nscores and sorting them, while in the second part, time complexity depends on the learning algorithm chosen. It is worth to note that the Table 1 Example of gene selection process by BIRS Rank g5g7g4g3g1g8g6g2g9 Subset Evaluated Acc. p-val Acc.-best sub 1g580 80 g5 2g5,g 782 3g5,g 481 4g5,g 383 5g5,g 184 <0.1 84 g5,g1 6g5,g1,g 884 7g5,g1,g 686 8g5,g1,g 289 <0.1 89 g5,g1,g2 9g5,g1,g2,g 987 learning algorithm is run N(number of genes) times with a small number of genes, only the selected ones. Therefore, the running time of the ranking procedure can be considered as negligible regarding the global process of selection. In fact, the results obtained from a random order of genes (without previous ranking) showed the following drawbacks: (1) the solution was not deterministic; (2) greater number of genes were selected; (3) the computational cost was higher because the classifier used in the evaluation contains more genes since the first iterations. Consider the situation depicted in Table 1: an example of the gene selection process done by BIRS. The first line shows the genes ranked according to some evaluation measure. We obtain the classification accuracy with the first gene in the list (g5:80%). In the second step, we run the classifier with the first two genes of the ranking (g5,g 7:82%), and a paired t-test is performed to determine the statistical significance degree of the differences. Since it is greater than 0.1, g7is not selected. The same happens with the next two subsets (g5,g 4:81%, g5,g 3:83%). Later, the gene g1is added, because the accuracy obtained is significantly better than that with only g5(g5,g 1:84%), and so on. In short, the classifier is run nine times to select, or not, the ranked genes (g5,g 1,g 2:89%): once with only one gene, four times with two genes, three with three genes and once with four genes. Most of the time, the learning algorithm is run with few genes. In short, this wrapper-based approach needs much less time than others with a broad search engine. As we can see in the algorithm, the first gene is always selected. This does not mean a great shortcoming in highdimensional databases, because usually several different sets of genes share similar information. The main disadvantage of sequential forward generation is that it is not possible to consider certain basic interactions among genes, i.e., genes that are useless by themselves can be useful together. Backward generation remedies some problems although there still will be many hidden interactions (in the sense of being unobtainable), but it demands more computational resources than the forward approach. The computer-load necessities of the forward search might become very inefficient in high-dimensional domains, as it starts with the original set of attributes and removes genes increasingly. 5. Experimental results The aim of this section is to evaluate our approach in terms of classification accuracy, degree of dimensionality and speed in selecting genes, in order to see how good BIRS is in situations where there is a large number of genes. Four public microarray data sets were used to assess the performance of the algorithm. Colon cancer data set: This data set is a collection of expression measurements from colon biopsy samples reported by Alon et al. [2]. The data set consists of 62 samples of colon epithelial cells. These samples were collected from colon cancer patients. The “tumor” biopsies were collected from tumors, and the “normal” biopsies were collected from healthy parts of the colons of the same patients. The final assignments of the status of biopsy samples were made by pathological examination. Gene expression levels in these 62 samples were measured using high-density oligonucleotide arrays. Of the ≈6000 genes represented in these arrays, 2000 genes were selected based on the confidence of the measured expression levels. Leukemia data set: This data set is a collection of expression measurements reported by Golub et al. [3]. The data set contains 72 samples. These samples are divided into two variants of leukemia: 25 samples of acute myeloid leukemia (AML) and 47 samples of acute lymphoblastic leukemia (ALL). The source of the gene expression measurements was taken from 63 bone marrow samples and 9 peripheral blood samples. Gene expression levels in these 72 samples were measured using high-density oligonucleotide microarrays. The expression levels of 7129 genes are reported. Lymphoma data set: This data set is a collection of expression measurements reported by Alizadeh et al. [4]. The data comprise 96 samples described by the expression values of 4026 genes. There are nine different subtypes of lymphoma. GCM data set: This data set is a collection of expression measurements reported by Ramaswamy et al. [22]. The data set contains 190 samples. These samples are divided to 14 variants of tumor. The expression levels of 16 063 genes are reported. Colon and leukemia are broadly studied, but not lymphoma and GCM, perhaps due to their difficulties at classification tasks. We chose three different learning algorithms, C4.5, IB1 and Naı¨ve Bayes, to evaluate the accuracy of selected genes for each gene selection algorithm. As already mentioned, the proposed search was done over a ranking of genes, and any evaluation measure could be used for it. In the experiments, we used two criteria: one belongs to the wrapper model, and other to the filter model. In the wrapper approach, denoted by BIRSW, we order genes according to their individual predictive power, using as criterion the performance of the target classifier built with a single gene. In the filter approach, a ranking is provided using a non-linear correlation measure. We chose symmetrical uncertainty (denoted by BIRSF), based on entropy and information gain concepts. Due to the high dimensionality of data, we limited our comparison to SF techniques and FCBF algorithm [17]. We chose three representative subset evaluation measures in combination with SF search engine. One, denoted by SFW, uses a target learning algorithm to estimate the worth of gene subsets; the other two are subset search algorithms which exploit SF search and use correlation measures (variation of CFS, correlation-based feature selection algorithm [23]) or consistency measure (variation of FOCUS [24]) to guide the search, denoted by CFSSF and FOCUSSF , respectively (both of them used in [17]). The experiments were conducted using the WEKA’s implementation of all these existing algorithms and our algorithm is also implemented in the WEKA environment [25]. We must take into account that the proper way to conduct a cross-validation for feature selection is to avoid using a fixed set of features selected with the whole training data set, because this induces a bias in the results. Instead, one should withhold a pattern, select features, and assess the performance of the classifier with the selected features using the left out examples. The results reported in this section were obtained with a 10-fold cross-validation over each data set, i.e. a feature subset was selected using the 90% of the instances, then, the accuracy of this subset was estimated over the unseen 10% of the data. This was performed 10 times, each time proposing a possible different feature subset. In this way, estimated accuracies, selected attribute numbers and time needed were the result of a mean over 10 cross-validation samples. Ambroise and McLachlan [26] recommend to use 10-fold rather than leave-oneout cross-validation, because the last one can be highly variable. Standard methods have been used for the experimental section (SF; Naïve Bayes, IB1 and C4.5 classifiers; and the t-Student statistical test). There exist other methods following the wrapper approach to extract relevant genes, which involve the selection process into the learning process (neural networks, Bayesian networks, support vector machines), although the source code of these methods is not freely available and therefore the experiments cannot be reproduced. In fact, some of them are designed for specific tasks, so the parameter setting is quite different for the learning algorithm. Table 2 shows the results obtained with the two BIRS approaches. Tables 3–5, report accuracy and number of genes selected by Naı¨ve Bayes, IB1 and C4.5, respectively, by each gene selection algorithm and the original set. We conducted a Student’s paired two-tailed t-test in order to evaluate the statistical significance of the difference between two averaged accuracy values: one resulted from BIRSWand the other resulted from one of BIRSF,SFW,CFSSF ,FOCUSSF , FCBF and the original set. The symbols “+” and “−”, Table 2 BIRS accuracy of Naïve Bayes (NB), IB1 (IB) and C4.5 (C4) on selected genes Data Colon Leukemia Lymphoma GCM Classif. Rank Acc. #g Acc. #g Acc. #g Acc. #g NB BIRSW85.48 3.50 93.04 2.50 82.14 10.30 67.37 44.00 BIRSF85.48 7.40 93.04 2.80 76.11 9.80 65.26 40.90 IB BIRSW79.76 6.30 93.04 3.30 85.56 16.40 58.95 37.00 BIRSF79.05 7.60 88.75 4.20 83.44 19.60 60.00 33.10 C4 BIRSW83.81 2.90 88.57 1.20 80.00 8.80 46.84 9.80 BIRSF83.81 2.90 84.64 1.10 81.33 6.60 44.21 24.80 Acc. records 1 ×10CV classification rate (%) and #g records the number of genes selected by each algorithm. Rank indicates if the classifier is using a wrapper (W) or a filter ranking (F). Table 3 Accuracy of Naïve Bayes on selected genes Data Colon Leukemia Lymphoma GCM Algorithm Acc. #g Acc. #g Acc. #g Acc. #g BIRSW85.48 3.50 93.04 2.50 82.14 10.30 67.37 44.00 SFW84.05 5.90 87.32 3.20 83.56 7.10 N/A CFSSF 82.62 22.10 91.43 40.30 75.11 153.22 N/A FOCUSSF 77.14 4.60 84.82−2.40 70.07 3.90 56.84 12.20 FCBF 77.62 14.60 95.89 45.80 78.22 290.90 68.95 60.90 Original 53.33−98.57 75.11 65.79 Acc. records 1 ×10CV accuracy rate (%) and #g records the number of genes selected by each algorithm. N/A—not available. Table 4 Accuracy of IB1 on selected genes Alg. Colon Leukemia Lymphoma GCM Acc. #g Acc. #g Acc. #g Acc. #g BIRSW79.76 6.30 93.04 3.30 85.56 16.40 58.95 37.00 SFW66.67−4.80 88.93 2.30 80.11 8.40 N/A CFSSF 80.71 22.10 90.18 40.30 92.78 153.22 N/A FOCUSSF 69.29 4.60 81.96−2.40 61.22−3.90 46.84−12.20 FCBF 80.71 14.60 94.46 45.80 91.89 290.90 61.05 60.90 Original 77.62 86.25 84.33 57.37 Acc. records 1 ×10CV accuracy rate (%) and #g records the number of genes selected by each algorithm. N/A—not available. Table 5 Accuracy of C4.5 on selected genes Alg. Colon Leukemia Lymphoma GCM Acc. #g Acc. #g Acc. #g Acc. #g BIRSW83.81 2.90 88.57 1.20 80.00 8.80 46.84 9.80 SFW80.71 3.30 87.32 1.60 73.00 8.20 N/A CFSSF 86.90 22.10 84.82 40.30 79.22 153.22 N/A FOCUSSF 79.05 4.60 88.93 2.40 62.44−3.90 49.47 12.20 FCBF 88.33 14.60 83.21 45.80 78.22 290.90 52.63 60.90 Original 82.14 82.14 81.44 60.00 Acc. records 1 ×10CV accuracy rate (%) and #g records the number of genes selected by each algorithm. N/A—not available. Table 6 Running time (s) for each gene selection algorithm Data Classif. BIRSWBIRSFSFWCFSSF FOCUSSF FCBF Colon NB 27.60 19.89 127.09 13.73 1.40 0.42 IB 42.90 32.98 120.06 13.71 1.39 0.42 C4 56.27 55.07 156.92 13.71 1.40 0.42 Leukemia NB 105.47 71.30 226.12 498.33 6.65 1.96 IB 145.86 99.26 228.30 498.46 6.63 1.96 C4 165.63 118.78 189.74 498.87 6.66 1.98 Lymphoma NB 309.44 304.58 1357.65 1571.52 5.79 8.71 IB 235.75 202.38 1032.81 1571.80 5.78 8.60 C4 579.95 428.43 4540.74 1570.68 5.79 8.68 Subtotal 1668.87 1332.68 7979.42 6250.81 41.50 33.14 GCM NB 11060.86 10626.81 N/A 177.33 26.28 IB 8669.74 7888.31 N/A 177.39 26.24 C4 32328.50 25624.29 N/A 177.19 26.28 Total 53727.97 45472.08 N/A 573.41 111.95 N/A—not available. respectively, identify statistic significance, at 0.05 level, wins or losses over BIRSW. Before comparing our technique to the others, note the similarity in Table 2 between the results obtained with the two approaches of our algorithm, one based on a rankingwrapper (BIRSW)and the other on a ranking-filter (BIRSF). As we can see in Table 2, in all the cases these accuracy differences are not statistically significant. The number of genes selected are similar too, although BIRSFis a little bit faster than BIRSWbecause of the time needed to build the ranking for the wrapper-ranking approach (see Table 6). Apart from the previous comparison, we studied the behavior of BIRSWin three ways in Tables 3–5: with respect to a whole set of genes (last row, original); with respect to another wrapper approach (SFW); and to three filter approaches (CFSSF ,FOCUSSF and FCBF). Table 6 reports the running time for each gene selection algorithm, showing three different results for each wrapper approach, depending on the learning algorithm chosen. Obviously, time needed for filter approaches are approximately the same in the three cases, because filters do not depend on the classifier used. Classification accuracies obtained with our wrapper approach are not statistically different than those obtained with the original set of genes, except for colon data set with NB, where BIRSWwins (see Table 3). We noticed that the number of selected genes was drastically low with regard to the original set, retaining 0.0018% of the genes on average for the three classifier. 5.1. BIRSWversus SFW No significant statistical differences are shown between the accuracy of our wrapper approach and the accuracy of the SF wrapper procedure (SFW), except for colon data set and IB classifier, where BIRSWwins (Table 4). On the other hand, the advantage of BIRSWwith respect to the SFWfor NB, IB1 and C4.5 is clear. We can observe (see Table 6) that BIRSWis consistently faster than SFW, because the wrapper subset evaluation is run less times. For example, for lymphoma data set and C4.5 classifier, BIRSWand SFWretain 8.80 and 8.20 genes, respectively, on average. To obtain these subsets the first one evaluated 4026 genes individually (to generate the ranking) and 4026 subsets, while the second one evaluated 32 180 subsets (4026 genes+ 4025 pairs of genes +···+4019 sets of eight genes). The time savings of BIRSWbecame more obvious when the computer-load necessities of the mining algorithm increased. In many cases the time savings were in degrees of magnitude (leukemia and lymphoma), and we must take into account that SFWdid not report any results on GCM data set after 1 week running. These results verify the computational efficiency of incremental search applied by BIRSWover greedy sequential search used by SFW, with similar number of genes selected and without significant statistical differences on accuracy. 5.2. BIRSWversus FOCUSSF With respect to the sequential search combined with consistency measure as subset evaluator (FOCUSSF ), classification accuracies obtained with this filter are lower than those obtained with BIRSWin five cases for the three classifier. In this sense, FOCUSSF retains very few genes, and in all cases, except for two, accuracies obtained with the original set are greater than those obtained with FOCUSSF . We noticed that the computer-load necessities of this filter procedure can be considered as negligible regarding wrapper models (Table 6). 5.3. BIRSWversus CFSSF With respect to CFSSF algorithm, where the sequential search is carried out together with a subset evaluation based on correlation measure, the error rate produced by our approach is not significantly different. However, gene subsets and time needed are not similar. Firstly, for the last data set (GCM) no results were provided by CFSSF because the program ran out of memory after a long period of time due to its quadratic space complexity. Secondly, considering that gene subsets provided by filters are not dependent on the classifier used later and those by wrappers are, CFSSF obtained on average 22.10,40.30 and 153.22 genes for colon, leukemia and lymphoma data sets, respectively, while BIRSWretained 3.5, 2.5 and 10.3 genes for NB classifier, 6.3, 3.3 and 16.4 genes for IB, and 2.9, 1.2 and 8.8 genes for C4. The time needed to reduce leukemia and lymphoma data sets by CFSSF is almost four times the time used by BIRSW. It is certainly true that CFSSF has problems when data sets have high dimensionality. 5.4. BIRSWversus FCBF Accuracies obtained with FCBF algorithm are very similar to those obtained with CFSSF , but with higher number of genes and much less time needed. Differences in accuracy between BIRSWand FCBF are not statistically significant. However, the number of genes selected by FCBF are much higher than that of our algorithm. Gene subsets provided by FCBF are 12 times greater than those provided by BIRS, i.e. FCBF retains 0.0224% of the genes on average for the four data sets, while BIRSWretain only 0.0018% of the genes on average for all data sets and the three classifiers. The computational cost of FCBF is very low (see Table 6). We used WEKA implementation of the FCBF algorithm with default values. However, if the threshold by which genes can be discarded is modified the results obtained might vary. In Tables 3–5, accuracy obtained on average with FCBF for all data sets and all classifiers is 79.27% retaining 0.0224% of the original set of genes. If the threshold is set to 0.25, results obtained are 77.35% of accuracy with 0.0156% of genes, and if 0.50 is fixed as a threshold, results are 59.81% and 0.0011%, respectively. This percentage of genes retained is similar to the obtained with our algorithm, 0.0018%, although BIRSWprovides a higher averaged accuracy, 78.74%. 5.5. Biological interpretation To complete our study, we reduced each data set by running our gene selection method on the original data sets. BIRS procedure always choose the top gene of each ranking, but generally the rest of the genes are not located at consecutive positions. For colon data set and NB classifier, BIRS chooses the genes M63391 (human desmin gene), H25136 (Inositol 1,4,5-trisphosphate-binding protein type 2 receptor), M64231 (human spermidine synthase gene) and R80427 (C4-dicarboxylate transport sensor protein DCTB) ranked at positions 1, 127, 159 and 160, respectively. For IB, H77597 (H. sapiens mRNA for metallothionein), X12671 (human gene for heterogeneous nuclear ribonucleoprotein (hnRNP) core protein A1) and H24956 (protooncogene tyrosine-protein kinase receptor rate precursor) at positions 1, 6 and 215, respectively. And J02854 (myosin regulatory light chain 2), D00860 (ribose-phosphate pyrophosphokinase I) and M26383 (human monocyte-derived neutrophil-activating protein mRNA) at positions 1, 6 and 1033, respectively, for the C4 classifier. The first gene of each subset appears in the list of relevant genes detected by previous studies over these data sets [6,27,28]. Similar behavior appears in the rest of data sets. For leukemia and NB classifier, BIRS chooses the genes M84526, M27891, M31523 and M23197 among the top20 genes of the ranking and M36652. For IB, M23197, M27891, M31523 and M11722, all of them among the top20. And only one gene, the first (M27891) for C4 classifier. The first gene of each subset appears in the list of relevant genes detected by previous studies over this data sets [6,27,28]. For lymphoma and NB classifier, BIRS chooses three genes among the top-20 genes and five further. For IB, seven and three respectively. And two genes and six for C4 classifier (see [29]). Finally, for GCM, BIRS chooses three, three and five genes for NB, IB and C4, respectively, among the top-20 genes of the ranking [30,31]. 6. Conclusions The success of many learning schemes, in their attempts to construct data models, hinges on the reliable identification of a small set of highly predictive attributes. Traditional gene selection methods often select the top-ranked genes according to their individual discriminative power. However, the inclusion of irrelevant, redundant and noisy genes in the model building process phase can result in poor predictive performance and increased computation. The most popular search methods in machine learning cannot be applied to microarray expression data sets due to the very high dimensionality, especially when a wrapper approach is used as evaluation function. We use the incremental ranked usefulness definition to decide at the same time whether or not a gene is relevant and non-redundant. The technique extracts the best non-consecutive genes from the ranking, trying to statistically avoid the influence of unnecessary genes in further classifications. Our approach, named BIRS, uses a very fast search through the attribute space and any classifier can be embedded into it as evaluator. Very high-dimensional datasets take a lot of computational resources when wrappers are chosen. BIRS reduces the search space complexity as it works directly on the ranking, transforming the combina-