A Hybrid Metaheuristic for Biclustering Based on Scatter Search and Genetic Algorithms

A Hybrid Metaheuristic for Biclustering Based on Scatter Search and Genetic Algorithms Juan A. Nepomuceno1, Alicia Troncoso2, and Jesús S. Aguilar–Ruiz2 1Department of Computer Science, University of Sevilla, Spain [email protected] 2Area of Computer Science, Pablo de Olavide University of Sevilla, Spain {ali,aguilar}@upo.es Abstract. In this paper a hybrid metaheuristic for biclustering based on Scatter Search and Genetic Algorithms is presented. A general scheme of Scatter Search has been used to obtain high–quality biclusters, but a way of generating the initial population and a method of combination based on Genetic Algorithms have been chosen. Experimental results from yeast cell cycle and human B-cell lymphoma are reported. Finally, the performance of the proposed hybrid algorithm is compared with a genetic algorithm recently published. Keywords: Biclustering, Gene Expression Data, Scatter Search, Evolutionary Computation. 1 Introduction Recently, data mining techniques are being applied to microarray data analysis in order to extract useful information [1]. Clustering techniques find groups of genes with similar patterns from a microarray. However, genes are not necessary related to every condition. Thus, the goal of the biclustering is to identify genes with the same behavior only under a specific group of conditions. In the context of microarray analysis, many approaches have been proposed for biclustering [2]. Biclustering techniques have two important aspects: the search algorithm and the measure to evaluate the quality of biclusters. Most of proposed approaches in the literature are focussed on different search methods. Thus, in [3] an iterative hierarchical clustering is applied to each dimension separately and biclusters are built by the combination of the obtained results for each dimension. In [4] an iterative search method which built biclusters adding or removing genes or conditions in order to improve the measure of quality called Mean Squared Residue (MSR) was presented. An exhaustive biclusters enumeration by means a bipartite graph-based model that nodes were added or removed in order to find maximum weight subgraphs was generated in [5]. The FLOC algorithm [6] improved the method presented in [4] obtaining a set of biclusters simultaneously and adding missing values techniques. In [7], a simple linear model for gene expression was used assuming normally distributed expression level for each gene or condition. Also, geometrical characterizations such as hyperplanes in a high dimensional data space have been used to find biclusters [8]. In the last few years, global optimization techniques, such as Simulated Annealing [9] or Evolutionary Computation [10,11], have been applied to obtain biclusters due to their good performance in several environments. Recently, several papers were focussed on the measure proposed to evaluate the quality of biclusters. In [12] an analysis of the MSR was made, showing that this measure is good to find biclusters with shifting patterns but not scaling patterns. A new measure based on unconstrained optimization techniques was proposed in [13] as alternative to the MSR in order to find biclusters with certain patterns. In this paper a hybrid metaheuristic for biclustering based on Scatter Search and Genetic Algorithms (SS&GA) is presented. A general scheme of Scatter Search has been used to obtain high–quality biclusters, but a way of generating the initial population and a method of combination based on Genetic Algorithms have been chosen. Finally, the performance of the proposed hybrid algorithm is compared with a genetic algorithm recently published [10]. A Scatter Search has been selected due to the recent success obtained to solve different hard optimization problems and to references about the application of Scatter Search for biclustering have not been found in the literature. This paper is organized as follows. Section 2 presents basic concepts about Scatter Search. The description of the proposed metaheuristic is described in Section 3. Some experimental results from two real datasets and a comparison between the proposed method and a genetic algorithm are reported in Section 4. Finally, Section 5 outlines the main conclusions of the paper and future works. 2 Scatter Search Scatter Search [14] is an optimization algorithm based on populations introduced in the seventies. Recently, Scatter Search algorithms have been applied to many nonlinear and combinatorial optimization problems providing remarkable outcomes due to its flexibility to adopt different search strategies mainly. Basically, a standard Scatter Search can be summarized by the following steps: 1. Generate an initial population in a deterministic manner to assure the diversity of the population regarding a distance. 2. A set, called set of reference, is built with the best individuals from this population. The best individuals is not limited to a measure of individuals provided by a fitness function but a individual that improves the diversity can be added to the reference set. 3. New individuals are created by the deterministic combination of individuals of the reference set and all individuals of the reference set are selected to be combined. 4. The reference set is updated using the new individuals and the combination is repeated until the reference set does not change. 5. The reference set is rebuilt and if the maximum number of iterations is not reachedgotostep3. Therefore, the search strategies of a Scatter Search depend on a diversification method to generate the initial population, a method to built the reference set, a method to combine individuals and a method to rebuilt the reference set. The main differences between a Genetic Algorithm and a Scatter Search are the way of generating the initial population, as the initial population is generated randomly and deterministic, respectively; the selection of individual to create offspring, as a probabilistic procedure is applied to select parents in Genetic Algorithms and all individuals of the reference set are combined in Scatter Search; the evolution of the population, based on the survival of the best depending on the fitness function in Genetic Algorithms and the rebuilding method of reference set used in Scatter Search. Finally, the size of the population in Genetic Algorithms is bigger than that of the reference set in Scatter Search. A typical size in Genetic Algorithms is 100 and 10 in Scatter Search, due to that the combination method in Scatter Search takes into account all pairs of individuals to create new individuals. In short, the underlying idea of Scatter Search is to emphasize systematic processes against random procedures to generate populations, to create new individuals and to inject diversity to the population. 3 Description of the Algorithm In this section the proposed SS&GA algorithm to obtain biclusters is described, detailing the steps aforementioned in the previous section such as combination, generation, updating and rebuilding methods. The pseudocode of the proposed SS&GA algorithm is presented in algorithm 1. 3.1 Biclusters Codification and Generation Formally, a microarray is a real matrix composed by Ngenes and Mconditions. The element (i, j)of the matrix means the level of expression of gene iunder the condition j. A bicluster is a submatrix of the matrix Mcomposed by n≤N rows or genes and m≤Mcolumns or conditions. Biclusters are encoded by binary strings of length N+M[10]. Each of the first N bits of the binary string is related to the genes and the remaining M bits to the conditions from microarray M. For instance, the bicluster shown in Fig. 1 is encoding by the following string, 0010110000|01100 (1) Thus, this string codifies the bicluster composed by genes number 3,5and 6 and conditions 2and 3from a microarray comprising 10 genes and 5conditions. The initial population of biclusters is strictly randomly generated (typical in Genetic algorithms) without taking into account the diversity (typical in Scatter Search). Random strings composed by 0and 1are generated until nB biclusters are built, where nB is the size of the starting population, i.e. the number of biclusters. Algorithm 1. SS&GA for Biclustering INPUT Microarray M, penalization factors M1and M2, size of population nB,size of reference set, S, and maximum number of iterations, MaxIter. OUTPUT The reference set, RefSet begin Initialize Prandomly with nB biclusters //Building reference set R1←S/2best biclusters from P(according to their fitness function) R2←S/2most scatter biclusters, regarding R1,fromPR1(according to a distance). RefSet =R1∪R2 P=PRefSet //Initialization stable ←FALSE iter =0 while (iter ≤MaxIter)do //Updating reference set while (NOT stable) do A←RefSet B←CombinationMethod(RefSet) RefSet ←Sbest biclusters from RefSet ∪B if (A=RefSet)then stable ←TRUE end if end while //Rebuilding reference set R1←S/2best biclusters from RefSet (according to their fitness function) R2←S/2most scatter biclusters, regarding R1,fromPR1. RefSet =R1∪R2 P←PRefSet iter =iter +1 end while end 3.2 Building Reference Set The reference set comprises the best Sbiclusters of the initial population, P, where Sis the number of biclusters that belong to this set. The reference set is built taken into account both quality and scattering of biclusters. The quality of biclusters is measured evaluating the fitness function considered in the evolutive process. Thus, a bicluster is better than another if the fitness function value is lower than that of the second one. On the other hand, a distance must be defined in order to show how the scattering is introduced in the search space. In the proposed SS&GA approach the distance used is the Hamming distance. The Hamming distance between two binary strings is defined by the number of positions for which their corresponding 0/1 values are different. For example, the Hamming distance for strings 001001001|001 and 001011001|101 is 2. 1C 3C 2C 5C 4C 3.0 6.2 - 3.5 6.3 2.2 1G 2 - 2 1.2 - 1.3 - 5.1 3.1 2G 3G 7.4 0.1 0.1 4.0 9.7 4.1 1.3 2.2 3.0 - 8.3 - 4G 5G 5.7 0.1 0.1 3.2 - 1.2 6G 4.0 0.1 0.1 3.0 4.0 1.3 5.2 - 5.2 - 3.8 2.3 7G 1.0 3.0 1.4 1.3 5.2 8G 2.0 2.9 9.6 4.0 1.3 9G 1.0 - 3.0 3.0 5.0 3.0 01G 1.0 1.0 1.0 1.0 1.0 1.0 bicluster 00101100000|1100 codification microarray Fig. 1. Microarray and bicluster along with its codification Therefore, the reference set is formed by the S/2best biclusters from P(set R1) according to their fitness function and the S/2biclusters from PR1 (set R2) with the highest distances to the set R1according to the Hamming distance. 3.3 Combination Method and Updating Reference Set Combination method is the mechanism to create new biclusters in Scatter Search. All pairs of biclusters are combined generating S∗(S−1)/2new biclusters. In the SS&GA algorithm the typical uniform crossover operator used in Genetic Algorithms is the proposed combination method. This crossover operator is shown in Fig. 2. A binary mask is randomly generated and a child is composed by values from the first parent when there is a 1 in the mask, and from the second parent when there is a 0. The reference set is updated with the Sbest biclusters from the reference set and the new biclusters generated by the combination method according to the fitness function. This process is repeated iteratively until the reference set does not change. parent 1 1 0 1101 1 0 1 1100 0 1 1 1100 1child parent 2 mask 1 0 0100 1 Fig. 2. Uniform crossover operator of Genetic Algorithms 3.4 Rebuilding Reference Set After getting the stability of reference set in the updating process, this set is rebuilt to introduce diversity in the search process. This task is made by mutation operators in Genetic Algorithms. Thus, the reference set is composed by the S/2 best biclusters from the updated reference set (set R1) according to the fitness function and the S/2most distant from PR1according to the Hamming distance. 3.5 Biclusters Evaluation The fitness function is fundamental in order to evaluate the quality of biclusters. Cheng and Church proposed the MSR which measures the correlation of a bicluster. Given a bicluster comprising the subset of genes Iand the subset of conditions J, the MSR is defined as follows, MSR(I,J)= 1 |I||J| i∈I,j∈J R(i, j)2(2) where R(i, j)=eij −eIj −eiJ +eIJ (3) eIj =1 |I| i∈I eij (4) eiJ =1 |J| j∈J eij (5) eIJ =1 |IJ| i∈I,j∈J eij (6) In this work, biclusters with low residue and high volume are preferred. Therefore, the fitness function is defined by: f(B)=MSR(B)+M11 G+M21 C(7) where MSR(B)is the MSR of the bicluster B,M1and M2are penalization factors to control the volume of the bicluster B,andGand Care the number of genes and conditions of the bicluster B, respectively. The use of MSR in the fitness function considered in the proposed SS&GA algorithm allows to establish a comparison with a previous evolutionary-based biclustering method and the Cheng and Church algorithm. 4 Experimental Results Two well known datasets [4] have been used to shows the performance of the proposed SS&GA algorithm. The first dataset is the yeast Saccharomyces cerevisiae cell cycle expression and the second one is the human B-cells expression data originated from [15] and [16], respectively. Original data were preprocessed in [4] replacing missing values with random numbers. The Yeast dataset contains 2884 genes and 17 experimental conditions and the Human dataset consists of 4026 genes and 96 conditions. The main parameters of the proposed SS&GA algorithm are as follows: 200 for the initial population; 10 for the reference set and 20 for the maximum number of iterations. The penalization factor for the number of conditions has been chosen of one order of magnitude larger than the range in which the fitness function varies for both datasets. However, that of the number of genes has been chosen of same order of magnitude to the range of values of the fitness function for both datasets. The main goal of this choice is to test the influence of the penalization factors on the volume of the biclusters. 4.1 Yeast Data Set Table 1 presents several biclusters obtained by the application of the SS& GA approach from Yeast dataset. For each bicluster is shown an identifier of the bicluster, the value of its MSR, the number of genes and the number of conditions. It can be observed that high–quality biclusters have been obtained as the values of the MSR are lower than 220. Moreover, the volume of the obtained biclusters is satisfactory showing that the SS& GA approach find non–trivial biclusters. Concretely, biclusters and no clusters are obtained since the number of conditions is less than 17 always. Biclusters presented in Table 1 are shown in Fig. 3. Although biclusters are good taking into account the values of their MSR, in this figure their tends cannot be observed easily. This is due to the overlapping among biclusters as the same gene can be found in different biclusters. Fig. 4 shows the evolution of the average MSR, fitness function values and volume of the reference set throughout the evolutionary process for the Yeast dataset. The values of the MSR and the volume are represented in the axis on the left and that of the fitness function in the axis on the right. It can be noticed that the initial reference set improves the average MSR throughout the iterations Table 1. Results obtained by SS&GA algorithm for Yeast dataset Bicluster MSR Genes Conditions bi.1 74.72 10 13 bi.2 106.25 13 13 bi.3 125.9 22 13 bi.4 216.16 25 14 bi.5 97.04 26 11 bi.6 117.25 14 14 bi.7 136.67 25 13 bi.8 159.44 39 13 bi.9 121.89 26 11 2 4 6 8 10 12 0 100 200 300 400 Conditions Expression Value bi1 24681012 100 200 300 400 Conditions Expression Value bi2 24681012 150 200 250 300 350 400 Conditions Expression Value bi3 2 4 6 8 10 12 14 0 100 200 300 400 Conditions Expression Value bi4 246810 0 100 200 300 400 Conditions Expression Value bi5 2 4 6 8 10 12 14 150 200 250 300 350 400 Conditions Expression Value bi6 2 4 6 8 10 12 150 200 250 300 350 Conditions Expression Value bi7 24681012 100 200 300 400 500 Conditions Expression Value bi8 246810 0 100 200 300 400 Conditions Expression Value bi9 Fig. 3. Biclusters from Yeast dataset and the SS&GA algorithm converges in 8 iterations approximately. The average volume of the reference set decreases versus the number of iterations due to the non too large penalization factors have been chosen. 4.2 Lymphoma Data Set Table 2 presents information about several biclusters found by the SS&GA approach for Human dataset. The values of the MSR are considerably low since all are lower than 1100. Thus, it can be stated that obtained biclusters have a remarkable quality. Moreover, in general the obtained biclusters have a large number of genes, specially the bicluster number 1, 2 and 4. These biclusters are also represented in Fig. 5. Figure 6 presents the performance of the proposed algorithm for Human dataset. The evolution of the average MSR, fitness function values and volume 100 150 200 250 300 350 1 2 3 4 5 6 7 8 9 1011121314151617181920 Number of iterations Average residue 0 100 200 300 400 500 600 Average fitness function Volume MSR Fitness Function Fig. 4. Performance of the proposed SS&GA algorithm for Yeast dataset Table 2. Results obtained by SS&GA for Human dataset Bicluster MSR Genes Conditions bi.1 855.17 109 13 bi.2 813.70 127 12 bi.3 642.13 85 11 bi.4 815.74 122 10 bi.5 771.69 48 12 bi.6 595.69 44 9 bi.7 1074.10 56 13 bi.8 507.17 67 8 bi.9 794.07 70 11 Table 3. Comparison of the results obtained by SS&GA, SEBI and CC algorithms Algorithm-Dataset Avg. Residue Avg. gene num. Avg. cond. num. SS&GA–Yeast 128.37 (40.71) 22.23 (8.86) 12.78 (1.09) SS&GA–Human 763.27 (165.73) 80.89 (36.61) 11 (1.73) SEBI–Yeast 205.18 (4.49) 13.61 (10.38) 15.25 (1.37) SEBI–Human 1028.84 (29.19) 14.07 (5.39) 43.57 (6.20) CC–Yeast 204.29 (42.78) 166.71 (226.37) 12.09 (4.39) CC–Human 850.04 (153.91) 269.22 (204.71) 24.5 (20.92) for the reference set is shown. A good performance of the SS&GA technique and a fast convergence can be appreciated. The values of the fitness function decrease quickly and only ten iterations approximately are needed to find high–quality biclusters. In this case, the choice of penalization parameters to keep under control the volume of the biclusters provides a nearly constant volume in the last iterations. Finally, a comparison between the results obtained with the SS&GA algorithm and two representative techniques reported in the literature is provided. Concretely, the SS&GA algorithm is compared to SEBI [10] and Cheng and