New heuristics for the maximum diversity problem

New Heuristics for the Maximum Diversity Problem Geiza C. Silva, Marcos R. Q. de Andrade, Luiz S. Ochi, Simone L. Martins, and Alexandre Plastino Universidade Federal Fluminense, Departamento de Ciência da Computação Rua Passo da Pátria, 156 – Bloco E – 3o andar – Boa Viagem 24210-240, Niterói, RJ, Brazil {gsilva, mandrade, satoru, simone, plastino}@ic.uff.br Abstract. The Maximum Diversity Problem (MDP) consists in identifying, in a population, a subset of elements, characterized by a set of attributes, that present the most diverse characteristics between themselves. The identification of such solution is an NP-hard problem. Some heuristics are available to obtain approximate solutions for this problem. In this paper we propose different GRASP heuristics for the MDP, using distinct construction procedures and path-relinking techniques. Performance comparison between related work and the proposed heuristics is provided. Experimental results show that the new GRASP heuristics are quite robust and are able to find high-quality solutions in a feasible computational time. 1 Introduction Consider P = {p1 , . . . , pn } a set of elements, defined over the index set N = {1, 2, ..., n}, and pik , k ∈ L = {1, . . . , l}, the l attributes associated with each element pi . The maximum diversity problem (MDP) [3, 5, 6] consists in identifying a subset M from the population P , so that the m elements from M present the maximum possible diversity among them. The measure of diversity dij between a pair of elements (i, j) is calculated by a function applied on their attributes. This problem can be formulated as: Maximize z= n−1
n
dij xi xj , (1) i=1 j=i+1 subject to n
xi = m, (2) i=1 where xi is a binary variable indicating if an element i is selected to be a member of the subset M . Many applications [7] can be solved using the resolution of this problem, such as human resource management, measure of biodiversity, and VLSI design. 2 Glover et al. [5] presented mixed integer zero-one formulation for this problem, that can be used to solve small instances by exact methods. They also show that this problem belongs to the class of NP-hard problems. Some heuristics are available to obtain approximate solutions. Weitz and Lakshminarayanan [14] developed five heuristics to find groups of students with the most possible diverse characteristics, such as nationality, age and graduation level. They tested the heuristics using instances based on real data and implemented an exact algorithm for solving them. Constructive and destructive heuristics were presented by Glover et al. [6], who created instances with different sizes of population (maximum value was 30) and showed that the proposed heuristics obtained results close (2%) to the ones obtained by the exact algorithm, but much faster. Ghosh [3] proposed a GRASP (Greedy Randomized Adaptive Search Procedure) that obtained good results for small instances of the problem. In this work, we initially present two new GRASP heuristics for the MDP, using distinct construction procedures. These proposals are developed using the concept of reactive GRASP introduced by Prais and Ribeiro [9]. The proposed GRASP strategies are compared to previous work based on computational experiments for instances created with maximum population of 500 individuals. Thereafter, motivated by the good results achieved, we also propose an improvement on the GRASP heuristic that obtained the best results. We introduce a path-relinking technique into this best GRASP heuristic as an intensification mechanism. This paper is organized as follows. In Section 2 we describe the GRASP proposed by Ghosh [3] and present the new GRASP heuristics developed. Based on computational experiments, these strategies are then compared. In Section 3 we create another GRASP heuristic combining a path-relinking technique with the pure GRASP proposed in the previous section which presented the best behavior. The performance comparison of all these strategies is presented in Section 4. Concluding remarks are presented in Section 5. 2 GRASP Heuristics GRASP [2] is an iterative process, where each iteration consists of two phases: construction and local search. In the construction phase a feasible solution is built, and its neighborhood is explored by a local search. The result is the best solution found over all iterations. The construction phase of GRASP is an iterative process where, at each iteration, the elements c ∈ C that do not belong to the solution are evaluated by a greedy function g : C → ℜ+ , that estimates the gain of including it in the partial solution. They are ordered by their estimated value in a list called restricted candidate list (RCL) and one of them is randomly chosen and included in the solution. The size of the RCL is limited by a parameter α. For a maximization problem, only the elements whose g values are in the range [(1 − α)gmax , gmax ] 3 are placed in RCL. This iterative process stops when a feasible solution is obtained. The solutions generated by the construction phase are not guaranteed to be locally optimal. So a local search is performed to attempt to improve each constructed solution. It works by successively replacing the current solution by a better one from its neighborhood, until no more better solutions are found. In Subsection 2.1, we describe the GRASP proposed by Ghosh [3]. In Subsection 2.2, we present two new GRASP procedures with different construction phases for the maximum diversity problem. 2.1 Ghosh Algorithm Ghosh [3] proposed a GRASP heuristic for the MDP and tested it for small instances. The construction phase consists of m iterations. For each iteration k, Mk−1 is a partial solution with k − 1(1 ≤ k ≤ m) elements. The element i∗ to be inserted in each iteration is selected based on its contribution to the overall diversity z(M ). First, a lower bound ∆zL (i) and an upper bound ∆zU (i) are computed for all i ∈ N \ Mk−1 . Then a random number u is sampled from a uniform distribution U (0, 1) that is used to compute ∆z ′ (i) = (1 − u)∆zL (i) + u∆zU (i). The selected element i∗ is the one that presents the larger ∆z ′ and is included in Mk−1 to obtain Mk . The computations of ∆zL (i) and ∆zU (i) are:   ∆zL (i) = j∈Mk−1 dij + n−m+1≤r≤n−k dri (Qik );   ∆zU (i) = j∈Mk−1 dij + 1≤r≤m−k dri (Qik ); where dri (Qik ) is the rth largest distance in {dij : j ∈ Qik }, Qik = N \ (Mk−1 ∪ {i}). The lower bound ∆zL (i) is computed by adding the distances between i and the elements that are already in the solution to the sum of distances between i and the m− k − 1 elements that are not in solution Mk−1 and which have smaller distances to i. The upper bound ∆zU (i) is computed by adding the same first term of ∆zL (i) to the sum of distances between i and the m − k − 1 elements that are not in solution Mk−1 and which have larger distances to i. After a solution is constructed, a local search phase should be executed in trying to improve the initial solution. The neighborhood of a solution defined by Ghosh [3] is the set of all solutions obtained by replacing an element in the solution by another that does not belong to the set of its elements. The incumbent solution M is initialized with the solution obtained by the construction phase. For each i ∈ M and j ∈ N \ M , the improvement due to exchanging i by j, ∆z(i, j) = u∈M{i} (dju − diu ) is computed. If for all i and j, ∆z(i, j) < 0, the local search is terminated, because no exchange will improve z. Otherwise, after the elements of the pair (i, j) that provides the maximum ∆z(i, j) are interchanged, a new incumbent solution M is created and the local search is performed again. 4 2.2 New Algorithms Usually, the local search phase demands great computational effort and execution time, so the construction phase plays an important role to diminish this effort by supplying good starting solutions for the local search. Next, we propose two construction procedures which are combined with the local search developed by Ghosh, obtaining two new GRASP procedures. These new construction procedures are based on two existing techniques: the widely used filtering of constructed solutions [10] and reactive GRASP proposed by Prais and Ribeiro [9]. Filtering of constructed solutions leads to a more greedy construction. For each GRASP iteration, the construction algorithm is executed X times generating X solutions and only the best solution is selected to be used as the initial solution for the local search phase. Prais and Ribeiro [9] proposed the procedure called reactive GRASP, for which the parameter α used in the construction phase is self adjusted for each iteration. For the first construction iteration, an α value is randomly selected from a discrete set A = {α1 , . . . , αm }. Each element αi has a probability pi associated and, initially, a uniform distribution is applied, thus we have pi = 1/m, i = 1, . . . , m. Periodically the probability distribution pi , i = 1, . . . , m is updated using information collected during the former iterations. The aim is to associate higher probabilities to values of α that lead to better solutions and lower ones to values of α that guide to worse solutions. The two new GRASP algorithms described below use both techniques described before. KLD (K Largest Distances) This algorithm constructs an initial solution by randomly selecting an element from a RCL of size K at each construction iteration. The RCL is created by selecting, for each element i ∈ N , the K elements j ∈ N \{i} that exhibit larger values of dij and computing si by the sum of these K values. Then, a list is created containing all elements i sorted in descending order by their si values and the K first elements are selected to compose the RCL list. The procedure developed to implement the reactive GRASP starts considering m it to be the total number of GRASP iterations. In the first block of iterations B1 = 0.4m it, four different values for K ∈ {K1 , K2 , K3 , K4 } are evaluated by dividing the block into four equal intervals ci , i = 1, . . . , 4. The value Ki is used for all iterations belonging to interval ci . The values of Ki are shown in Table 1, where µ = (n − m)/2. After the execution of the last obtained for each Ki is evaluiteration of block B1 , the quality of the solutions  ated. The average diversity value zmi = 1≤q≤0.1m it z(soliq ) for the solutions soliq , i = 1, . . . 4; q = 1, . . . , 0.1m it; is obtained for each Ki . The values Ki are stored in a list LK ordered by their zmi values. Then the next block of iterations B2 = 0.6m it is divided into four intervals yi , each one with different number of iterations. For each interval a value for K 5 Table 1. K values for block B1 i 1 2 3 4 ci K [1, . . . , 0.1m it] m + µ − 0.2µ (0.1m it, . . . , 0.2m it] m + µ − 0.1µ (0.2m it, . . . , 0.3m it] m + µ + 0.1µ (0.3m it, . . . , 0.4m it] m + µ + 0.2µ is adopted as shown in Table 2. In this way, the values Ki that provide better solutions are used in a larger number of iterations. Table 2. K values for block B2 i 1 2 3 4 yi (0.4m it, . . . , 0.64m it] (0.64m it, . . . , 0.82m it] (0.82m it, . . . , 0.94m it] (0.94m it, . . . , m it] K lk1 lk2 lk3 lk4 At each GRASP iteration, the filter technique is applied by constructing 400 solutions and only the best solution is sent to the local search procedure. The pseudo-code, including the description of the procedure for the construction phase using K largest distances heuristic, is given in Fig. 1. In line 1, the cost of the best solution found in the execution of max sol f ilter iterations is initialized. The value K to be used to build the Restricted Candidate List (RCL) is calculated by the procedure det K in line 2. This procedure defines the value for K, implementing the reactive GRASP described before. In line 3, the RCL is built. From line 4 to line 15, the construction procedure is executed max sol f ilter times and only the best solution is returned to be used as an initial solution by the local search procedure. From line 6 to line 10, a solution is constructed by the random selection of an element from RCL. In lines 11 to 14, the best solution found by the construction procedure is updated. If the GRASP iteration belongs to B1 block, the cost of the solution found using the selected K is stored in line 17. When the first block B1 of iterations ends, the values Ki are evaluated and put in the list LK sorted in descending order, in line 20. The local search strategy employed is the same as the one developed by Ghosh. KLD2 This algorithm is similar to the previously described algorithm. The difference between them is the way that the Restricted Candidate List is built. In the former algorithm, the RCL is computed before the execution of the construction iterations and, for each iteration, the only modification made in the RCL is the removal of the element that is inserted in the solution. 6 procedure constr KLD(it GRASP, m it, numsol, n, m) 1. best cost sol ← 0; 2. K ← det K(it GRASP, m it, LK); 3. RCL ← Build RCL(K); 4. for j = 1, . . . , max sol f ilter do 5. sol ← {}; 6. for k = 1, . . . , m do 7. Randomly select an individual e∗ from RCL; 8. sol ← sol ∪ {e∗ }; 9. RCL ← RCL − {e∗ }; 10. end for; 11. if (z(sol) > best cost sol) then do 12. sol constr ← sol; 13. best cost sol ← z(sol); 14. end if 15. end for; 16. if (it GRASP < 0.4m it) then do 17. U pdate Sol K(K, sol eval, z(sol constr)); 18. end if; 19. if (it GRASP == 0.4m it) then do 20. LK ← Build LK(sol eval); 21. end if; 22. return sol constr. Fig. 1. Construction procedure used to implement the KLD heuristic The first element of the constructed solution is selected in the same way performed by KLD heuristic, which means that an element is randomly selected from the RCL built as described in line 3 of Fig. 1. In KLD2, the following elements are selected from an RCL built using an adaptive procedure. Let Mc be a partial solution with c, 1 ≤ c < m, elements solution Mc+1 . For and i ∈ N \Mc a candidate to be inserted in the next partial  each i, we select the (K − c − 1) elements j ∈ N \(Mc {i}),that present larger values of dij and calculate the sum of the (K − c − 1) values of dij obtaining si . To select the next element to be inserted, an initial candidate list is created based on the greedy function gf (i) shown in (3), where the first term corresponds to the sum of distances from the candidate i to the elements j ∈ Mc , and the second term stands for the sum of distances from element i to the (K − c − 1) elements that are not in the solution Mc and present larger distances to i. The initial candidate list is formed by the elements i, sorted in descending order with respect to gf (i), and the first K elements are selected from this list to build the RCL.
dij + si (3) gf (i) = j∈Mc The Reactive GRASP and the construction filter are implemented in the same way as in KLD. Once this construction algorithm demands much more execution 7 time than KLD algorithm, only 2 solutions, instead of 400, are generated to be filtered. 2.3 Comparing the three strategies The computational experiments were performed on four sets (A, B, C, and D) of test problems with different characteristics [13]. An instance of a set consists of a population P with size n and a diversity matrix DivM at, which contains the diversity dij between elements i and j of P . The instances of a specific set are obtained from the same base diversity matrix, differing among themselves by the population size. Different base matrices are used to generate instances for each set. All sets A, B, and C contain instances with population sizes n = 50, 100, 150, 200, 250. The instances for set D have populations sizes n = 100, 200, 300, 400, 500. For all instances, tests were performed to find subsets M of sizes 10%, 20%, 30%, and 40% of the population size. The algorithms were implemented in C and compiled with gcc 3.3.2. The tests were performed on a 1.4 GHz AMD Athlon with 250 Mbytes of RAM. In Tables 3, 4, 5, and 6, the results of computing 500 iterations for each GRASP heuristic using instances A, B, C, and D are shown. The first and second columns identify two parameters of each instance: the size n of the population and the number m of elements to be selected. Each procedure was executed three times and for each one the average value of the solutions is presented. When different solutions are found, bold values indicate the best result found and underlined values represent the second best result. The three strategies found the same result in 39 from 80 tests. In those 41 tests in which different solutions were reached, both KLD and KLD2 heuristics found better solutions than Ghosh algorithm in 30 cases. The KLD2 algorithm found 31 best solutions, while KLD obtained 25 and Ghosh found 7. Comparing KLD and KLD2, individually to Ghosh heuristic, KLD outperformed Ghosh in 30 tests and KLD2 in 34 tests. Ghosh found better solutions than KLD in 9 tests and than KLD2 in 5 tests. When evaluating KLD2 against KLD, we can observe that, considering the 80 conducted tests, both strategies found the same result in 58 cases, KLD2 outperformed KLD in 16 tests, and KLD found better solutions in 6 cases. Tables 7, 8, 9, and 10 report the CPU times observed for the execution of the same sets of instances. The first and second columns identify the two parameters of each instance. For each GRASP heuristic, the average time for three executions is reported. Among the heuristics, algorithm KLD is the most efficient related to execution time. This strategy executed faster than the other two in 71 from 80 tests. Heuristic KLD2, which achieved better quality solutions, demanded more time than the other two heuristics for all instances. Next, a deeper analysis of the results obtained for the GRASP heuristics is performed. Some tests were executed for two instances from group D: the first one has parameters n = 200 and m = 40, and the second one, n = 300 and m = 90. Each GRASP heuristic was executed until a solution was found with a greater or equal cost compared to a target value. Two target values were used for each 8 Table 3. Comparison of Ghosh, KLD, and KLD2 for instances from set A n 50 50 50 50 100 100 100 100 150 150 150 150 200 200 200 200 250 250 250 250 m 5 10 15 20 10 20 30 40 15 30 45 60 20 40 60 80 25 50 75 100 Ghosh 86733.0 334976.0 692704.0 1171416.0 352286.0 1267277.0 2674152.0 4544642.0 753041.0 2758381.0 5825682.0 9959068.0 1321476.6 4788086.0 10211323.0 17544448.0 2041049.0 7388210.5 15788177.0 27148842.0 KLD 86733.0 334976.0 692704.0 1171416.0 353730.0 1267277.0 2674152.0 4544642.0 772982.0 2758381.0 5825149.5 9960461.0 1330354.6 4788086.0 10211323.0 17544448.0 2047828.0 7382340.5 15795352.0 27161422.0 KLD2 86733.0 334976.0 692704.0 1171416.0 353730.0 1267277.0 2674152.0 4544642.0 772982.0 2758381.0 5825682.0 9960461.0 1330146.6 4788078.0 10207669.0 17544448.0 2046995.4 7385769.5 15801103.0 27159518.0 Table 4. Comparison of Ghosh, KLD, and KLD2 for instances from set B n 50 50 50 50 100 100 100 100 150 150 150 150 200 200 200 200 250 250 250 250 m 5 10 15 20 10 20 30 40 15 30 45 60 20 40 60 80 25 50 75 100 Ghosh 80120.0 316409.0 659250.0 1094343.0 337420.0 1207522.0 2509045.0 4219476.0 737342.0 2613286.0 5519104.0 9374611.0 1294817.0 4630545.0 9811296.0 16759895.0 1983723.0 7178042.5 15303499.0 26047022.0 KLD 84517.0 316409.0 659250.0 1094343.0 337420.0 1207522.0 2509045.0 4219476.0 737342.0 2613286.0 5519104.0 9374611.0 1294817.0 4630545.0 9811296.0 16759895.0 1983723.0 7178042.5 15303499.0 26047022.0 KLD2 84517.0 316409.0 659250.0 1094343.0 337420.0 1207522.0 2509045.0 4219476.0 737342.0 2613286.0 5519104.0 9374611.0 1294817.0 4630545.0 9811296.0 16759895.0 1983723.0 7178042.5 15303499.0 26047022.0 9 Table 5. Comparison of Ghosh, KLD, and KLD2 for instances from set C n 50 50 50 50 100 100 100 100 150 150 150 150 200 200 200 200 250 250 250 250 m 5 10 15 20 10 20 30 40 15 30 45 60 20 40 60 80 25 50 75 100 Ghosh 93007.0 380750.3 851353.0 1502908.0 399449.0 1570800.0 3475608.0 6067776.0 896529.3 3498361.0 7748205.0 13609447.6 1590285.3 6199891.3 13744341.6 24124858.0 2483263.0 9668626.6 21440150.3 37718428.0 KLD 93007.0 381379.0 851353.0 1502908.0 399449.0 1570800.0 3475608.0 6067776.0 898585.0 3502567.0 7748205.0 13611261.0 1595768.0 6207580.0 13749403.0 24133320.0 2488888.0 9685273.0 21464360.0 37753120.0 KLD2 93007.0 381379.0 851353.0 1502908.0 399449.0 1570800.0 3475608.0 6067776.0 898585.0 3502567.0 7748205.0 13611261.0 1595768.0 6207580.0 13749403.0 24133320.0 2488888.0 9685430.0 21464360.0 37753120.0 Table 6. Comparison of Ghosh, KLD, and KLD2 for instances from set D n 100 100 100 100 200 200 200 200 300 300 300 300 400 400 400 400 500 500 500 500 m 10 20 30 40 20 40 60 80 30 60 90 120 40 80 120 160 50 100 150 200 Ghosh 318.0 1178.0 2457.0 4142.0 1233.0 4443.0 9437.0 16225.0 2666.0 9652.0 20725.0 35880.0 4615.0 16900.5 36298.7 62442.3 7079.3 26219.0 56571.0 97255.0 KLD 333.0 1195.0 2457.0 4142.0 1247.0 4446.3 9437.0 16182.3 2686.0 9667.7 20640.0 35871.0 4635.3 16916.3 36175.0 62313.0 7124.0 26197.0 56572.0 97213.0 KLD2 333.0 1195.0 2457.0 4142.0 1247.0 4447.0 9437.0 16225.0 2688.5 9682.5 20726.0 35878.0 4647.0 16930.0 36272.5 62478.0 7112.5 26229.0 56572.0 97316.5 10 Table 7. Computational time for Ghosh, KLD, and KLD2 for instances from set A n 50 50 50 50 100 100 100 100 150 150 150 150 200 200 200 200 250 250 250 250 m 5 10 15 20 10 20 30 40 15 30 45 60 20 40 60 80 25 50 75 100 Ghosh 2.4 4.6 6.7 8.7 15.5 38.7 64.4 93.5 56.2 146.4 288.1 449.3 146.6 425.4 820.5 1108.8 322.5 982.7 1854.2 2696.0 KLD 2.6 4.5 6.9 9.8 7.9 21.4 40.7 62.6 21.9 76.5 164.1 267.7 51.6 214.8 490.0 796.4 108.7 505.3 1167.2 1885.5 KLD2 3.9 7.2 10.7 14.0 28.7 62.1 100.9 139.4 102.7 240.8 415.4 607.0 263.1 644.6 1149.8 1700.1 536.5 1384.3 2568.6 3876.5 Table 8. Computational time for Ghosh, KLD, and KLD2 for instances from set B n 50 50 50 50 100 100 100 100 150 150 150 150 200 200 200 200 250 250 250 250 m 5 10 15 20 10 20 30 40 15 30 45 60 20 40 60 80 25 50 75 100 Ghosh 2.5 4.5 6.0 6.9 16.5 36.8 52.5 67.1 57.7 133.9 210.4 226.3 141.8 399.3 637.9 874.7 312.7 919.6 1701.9 2472.7 KLD 2.7 4.8 7.3 9.8 8.4 21.8 40.2 64.1 21.9 75.8 159.3 262.6 52.6 218.5 468.9 711.4 111.0 511.7 1158.9 1735.4 KLD2 3.8 7.4 11.2 15.0 28.5 62.8 103.1 159.3 103.0 246.2 432.5 652.3 263.5 676.6 1220.2 1788.4 559.3 1497.2 2860.7 4111.5 11 Table 9. Computational time for Ghosh, KLD, and KLD2 for instances from set C n 50 50 50 50 100 100 100 100 150 150 150 150 200 200 200 200 250 250 250 250 m 5 10 15 20 10 20 30 40 15 30 45 60 20 40 60 80 25 50 75 100 Ghosh 2.3 4.1 5.7 7.5 14.3 34.4 56.9 84.0 51.0 124.8 235.7 329.1 132.8 372.4 715.4 1159.4 301.6 979.0 2000.2 3000.2 KLD 2.3 4.0 6.1 8.9 7.6 20.6 39.1 58.7 19.6 69.1 149.3 241.9 47.8 201.8 468.5 716.5 143.4 605.2 1291.3 1864.1 KLD2 3.6 6.7 9.8 12.8 26.3 58.7 96.9 133.4 94.7 223.2 388.5 564.1 244.9 623.0 1137.5 1626.5 570.9 1578.5 2925.0 4273.1 Table 10. Computational time for Ghosh, KLD, and KLD2 for instances from set D n 100 100 100 100 200 200 200 200 300 300 300 300 400 400 400 400 500 500 500 500 m 10 20 30 40 20 40 60 80 30 60 90 120 40 80 120 160 50 100 150 200 Ghosh 18.0 47.4 88.3 130.9 147.5 561.9 1075.3 1647.0 712.7 2477.2 5603.2 9740.2 2469.3 9651.9 20210.3 32217.8 5122.4 26128.9 54392.6 80252.1 KLD 10.4 28.2 57.4 89.2 71.7 318.3 626.6 804.3 412.7 1690.5 3260.6 4453.0 1208.7 5266.5 10097.5 15750.3 2934.9 11921.7 29051.4 39043.3 KLD2 28.0 64.8 111.7 165.6 263.6 784.7 1958.8 3051.8 1159.4 3450.4 6829.5 11458.5 3202.7 10464.3 22038.4 32843.6 6988.6 26513.3 55898.1 84678.4 12 instance: the worst value obtained by these heuristics and an average of the values generated by them. Empirical probability distributions for the time to achieve a target value are plotted in Fig(s). 2 and 3. To plot the empirical distribution for each variant, each GRASP heuristic was executed 100 times using 100 different random seeds. In each execution, the time to achieve a solution whose cost was greater or equal to the target cost was measured. The execution times were sorted in ascending order and a probability pi = (i − 0.5)/100 was associated for each time ti and the points zi = (ti , pi ) were plotted for i = 1, . . . , 100. Each point zi in the graph represents the probability pi of finding the target value in time ti . This kind of analysis was first used by Aiex et al. [1]. We can observe that KLD behaves better than Ghosh and KLD2. KLD algorithm presented a higher probability of finding a target solution in a given amount of time than the others. Fig. 2 (target=4442) shows that the probability of finding the target value in 10 seconds is approximately 4% for both Ghosh and KLD2, while for KLD is 60%. Considering 100 seconds, we have the probabilities 26%, 50%, and 98% respectively. KLD2 outperformed Ghosh algorithm in these experiments, except in one test (Fig. 3, target=20640). 3 Path-relinking Path-relinking is a technique proposed by Glover [4] to explore possible trajectories connecting high quality solutions, obtained by heuristics like tabu search and scatter search. The pure GRASP metaheuristic is a memoryless method, because all iterations are independent and no information about the solutions is passed from one iteration to another. The objective of introducing path-relinking to a pure GRASP algorithm is to retain previous good solutions and use them as guides in the search of new good solutions. Laguna and Marti [8] were the first to use path-relinking within a GRASP strategy. Several extensions, improvements and successful applications of this technique can be found in the literature [11]. From the results presented in Subsection 2.3, we can see that the algorithm KLD2 presents some better results than KLD, but demands much more computational time. Besides, KLD presents a higher probability of reaching good solutions in a given amount of time. So we decided to introduce a path-relinking strategy to KLD in order to improve the trade-off between solution quality and computation time. Path-relinking is applied to a pair of solutions {si , sg }. First, the set of moves which should be executed to reach sg (the guiding solution) starting from si (initial solution) is computed. At each step, all moves needed are evaluated, the one which generates the best cost solution is executed and the set of available moves is updated. This procedure terminates when sg is reached. 13 1 0.8 n=200 m=40 Target: 4442 Probability 0.6 0.4 0.2 Ghosh KLD KLD2 0 1 10 100 1000 10000 Time (sec) 1 0.8 n=200 m=40 Target: 4443 Probability 0.6 0.4 0.2 Ghosh KLD KLD2 0 1 10 100 1000 10000 Time (sec) Fig. 2. Comparison of GRASP heuristics for the instance n = 200, m = 40 with targets values 4442 and 4443 14 1 n=300 m=90 Target: 20640 0.8 Probability 0.6 0.4 0.2 Ghosh KLD KLD2 0 1 10 100 1000 Time (sec) 1 n=300 m=90 Target: 20693 0.8 Probability 0.6 0.4 0.2 Ghosh KLD KLD2 0 1 10 100 1000 10000 Time (sec) Fig. 3. Comparison of GRASP heuristics for the instance n = 300, m = 90 with targets values 20640 and 20693 15 To use path-relinking within a GRASP procedure an elite set E is maintained, in which good solutions found in GRASP iterations are stored. Two basic strategies for introducing path-relinking into GRASP may be used [11]: – it is performed after each GRASP iteration using a solution from the elite set and a local optimum obtained after the GRASP local search. – it is applied to all pairs of elite solutions, either periodically or after all GRASP iterations terminate. The first strategy seems to be more effective than the second one [12], so this strategy was chosen. The elite set E is maintained, in which good solutions found in GRASP iterations are stored to be combined with solutions created by other GRASP iterations. The path-relinking is activated only after the elite set achieves its size. For each solution Sol generated in a GRASP iteration, one of the solutions e ∈ E is randomly selected and a path-relinking is applied between them. A path-relinking is performed by starting from an initial solution si and gradually incorporating attributes from a guide solution sg to it, until si becomes equal to sg . Several ways to select si and sg are used [11]: – Backward relinking: the best cost solution between Sol and e is set as the initial solution si and the other one as sg . – Forward relinking: the worst cost solution between Sol and e is set to si and the best one to sg . – Backward and Forward relinking: two different paths are explored, from si to sg and from sg to si . – Mixed relinking: both trajectories are explored in alternate way. At each iteration of path relinking, the roles of initial and guide solutions are inverted. Ribeiro et al. [12] observed that using the backward and forward strategy demands almost twice the time used to explore only one of them, and presents very small improvement in solution quality. They have also observed that best solutions are found when the relinking strategy uses the best cost solution as the initial solution si . The reason for this behavior is that the neighborhood of the initial solution is much more carefully explored than that of the guiding one, so there is a better chance to investigate in more detail the neighborhood of the most promising solution. Therefore, we implemented a backward path relinking strategy. Figure 4 shows the GRASP with path-relinking developed in this work. In line 1, the instance parameters are obtained and the elite set is initialized in line 3. For each GRASP iteration a solution is constructed in line 5 and, in line 6, a local search is executed using the KLD heuristic. If the elite set is full, the path relinking should be performed. In line 8, the element of the elite set is randomly selected according to uniform probability function and, in line 9, the path-relinking is performed. In line 11, the best solution between that generated by the GRASP iteration and the one by path-relinking is evaluated to verify if it should be inserted in the elite set. If the elite set is not full and the solution is 16 not already in the elite set, then it is inserted into it. If the elite set is full and the solution cost is better than the worst-cost elite set solution, then this new solution replaces the worst one. In lines 13 and 14 the best solution and best cost are updated. procedure GRASP PR 01. Read Input Data(n, m, DivM at, M axIter); 02. z best ← −∞; 03. E ← {}; 04. for i = 1, . . . , M axIter do 05. Sol ← Build Random Greedy Solution KLD(n, m, DivM at); 06. Sol ← Local Search KLD(Sol); 07. if E is full then do 08. SolElite ← Select Sol Elite(E); 09. Sol ← P ath Relinking(Sol, SolElite); 10. end if; 11. U pdate Elite(E, Sol); 12. if z(Sol) > z best do 13. Best Sol ← Sol; 14. z best ← z(Sol); 15. end if; 16. end for; 17. return Best Sol; Fig. 4. Algorithm for GRASP with path-relinking We show in Fig. 5 the details on our implementation of path-relinking for the MDP. In line 1, the solution generated by the GRASP iteration Sol and the solution from the elite size SolElite are compared to determine which has the better cost. BS will contain the best solution and W S the worst. In lines 3 and 4, the initial and guide solutions are set. The initial solution is set as the best-cost solution between Sol and SolElite. From line 5 to 12, the steps of path-relinking are performed until the initial solution reaches the guide solution. In line 6, the elements that are in the guide solution and are not in the initial solution are found. In line 7, a step of the path-relinking is performed. Intermediate solutions are obtained by replacing an element that belongs to the initial solution and do not belong to the guide solution with all other elements that belong to the guide solution and do not belong to the initial solution. The solution which presents the best cost among all these intermediate solutions is selected (IntSol) as the result of an iteration of path-relinking. This solution then is more similar to the guide solution because one element from the initial solution was replaced by another one from the guide solution. In line 9, if this new intermediate solution has a better cost than the current best intermediate solution (BestSolP R), then the latter is updated by this new one. In line 11, the intermediate solution is set as the initial solution for the next step of the path-relinking. 17 procedure Path Relinking(Sol, SolElite) 01. F ind Best W orst(Sol, SolElite, BS, W S); 02. BestSolP R ← BS; 03. Initial ← BS; 04. Guide ← W S; 05. while Initial not equal to Guide do 06. Dif El ← Obtain Dif f erent Elements(Initial, Guide); 07. IntSol ← Obtain Int Sol(Initial, Dif El); 08. if z(IntSol) > z(BestSolP R) then do 09. BestSolP R ← IntSol; 10. end if; 11. Initial ← IntSol; 12. end while; 13. return BestSolP R; Fig. 5. Algorithm for path-relinking The computational experiments implemented to evaluate the performance of this heuristic are presented in next section. 4 Computational Results The performance of the new developed GRASP with path-relinking was evaluated using all instances from set D. These instances are larger and, according to the results presented in Subsection 2.3, they are more difficult to solve. The new algorithm was also implemented in C and compiled with gcc 3.3.2. The tests were performed on a 1.4 GHz AMD Athlon with 250 Mbytes of RAM. In Table 11, we show the results of computing 500 iterations for the pure GRASP (KLD) and for KLD with path-relinking (KLD+PR) using instances from set D. The elite set size was 5 and the path-relinking strategy starts to be performed only after this set set is full, i.e., during the first five iterations, pathrelinking is not applied. The first and second columns identify the parameters of each instance: the size n of the population and the number m of elements to be selected. For all instances, each procedure was executed three times. The third and fourth columns show the average values obtained by KLD and KLD+PR strategies. For these average values, when two different results are reached, a bold value indicates the best average found. In order to compare the results obtained by KLD+PR to those reached by other proposed strategies (Ghosh, KLD, KLD2, or any other presented in the related literature [13]), the last column shows the best known average value for each test. The two algorithms obtained the same results for population size 100. Considering the other 16 tests, with population size larger or equal to 200, different average diversity values were reached in 12 tests. KLD+PR heuristic found better results than KLD algorithm in all these 12 cases. We can also observe that the average diversity solutions reached by KLD+PR are equal to the best known average values in almost all cases: in 16 from 20 tests. 18 Table 11. Comparing KLD and KLD+PR average results for instances from set D n 100 100 100 100 200 200 200 200 300 300 300 300 400 400 400 400 500 500 500 500 m 10 20 30 40 20 40 60 80 30 60 90 120 40 80 120 160 50 100 150 200 KLD 333.0 1195.0 2457.0 4142.0 1247.0 4446.3 9437.0 16182.3 2686.0 9667.7 20640.0 35871.0 4635.3 16916.3 36175.0 62313.0 7124.0 26197.0 56572.0 97213.0 KLD+PR 333.0 1195.0 2457.0 4142.0 1247.0 4450.0 9437.0 16225.0 2691.0 9684.0 20740.3 35879.0 4655.7 16935.0 36302.0 62467.3 7124.0 26246.3 56572.0 97333.3 Best Average 333.00 ([13],KLD,KLD2,KLD+PR) 1195.0 ([13],KLD,KLD2,KLD+PR) 2457.0 ([13],Ghosh,KLD,KLD2,KLD+PR) 4142.0 ([13],Ghosh,KLD,KLD2,KLD+PR) 1247.0 ([13],KLD,KLD2,KLD+PR) 4450.0 (KLD+PR) 9437.0 ([13],Ghosh,KLD,KLD2,KLD+PR) 16225.0 ([13],Ghosh,KLD2,KLD+PR) 2692.0 ([13]) 9684.0 (KLD+PR) 20740.3 (KLD+PR) 35881.0 ([13]) 4655.7 (KLD+PR) 16946.7 ([13]) 36302.0 (KLD+PR) 62478.0 (KLD2) 7124.0 (KLD,KLD+PR) 26246.3 (KLD+PR) 56572.0 (KLD,KLD2,KLD+PR) 97333.3 (KLD+PR) Ghosh heuristic obtained only four best known average values, KLD reached 8 values, and KLD2 9 values. Distinct algorithms developed in [13] found 10 best known average values. For seven instances, KLD+PR achieved a new best average which were not previously reached by any other strategy. In Table 12, we show the best diversity values achieved by both KLD and KLD+PR strategies in the conducted results. The third and fourth columns show the best values obtained by them. For these best values, when two different results are reached, a bold value indicates the best one. In order to compare the results obtained by KLD+PR to those reached by the other proposed strategies, the last column shows the best known value for each test. Again, the two algorithms obtained the same results for population size 100. Considering the other 16 tests, different best diversity values were reached in 10 tests. In all of them, KLD+PR heuristic found better results than KLD algorithm. The KLD+PR heuristic was able to find the best known solution in 15 from 20 tests. Ghosh heuristic obtained only six best known values, KLD reached 8 values, and KLD2 9 values. The algorithms developed in [13] together found 17 best known average values. However, we observe that, from these strategies, the one with best behavior achieved 12 best results. 19 Table 12. Comparing KLD and KLD+PR best results for instances from set D n 100 100 100 100 200 200 200 200 300 300 300 300 400 400 400 400 500 500 500 500 m 10 20 30 40 20 40 60 80 30 60 90 120 40 80 120 160 50 100 150 200 KLD 333.0 1195.0 2457.0 4142.0 1247.0 4448.0 9437.0 16207.0 2691.0 9689.0 20640.0 35871.0 4653.0 16925.0 36175.0 62313.0 7130.0 26201.0 56572.0 97213.0 KLD+PR 333.0 1195.0 2457.0 4142.0 1247.0 4450.0 9437.0 16225.0 2691.0 9689.0 20743.0 35879.0 4658.0 16945.0 36306.0 62487.0 7130.0 26254.0 56572.0 97344.0 Best Result 333.00 ([13],KLD,KLD2,KLD+PR) 1195.0 ([13],KLD,KLD2,KLD+PR) 2457.0 ([13],Ghosh,KLD,KLD2,KLD+PR) 4142.0 ([13],Ghosh,KLD,KLD2,KLD+PR) 1247.0 ([13],KLD,KLD2,KLD+PR) 4450.0 ([13],KLD+PR) 9437.0 ([13],Ghosh,KLD,KLD2,KLD+PR) 16225.0 ([13],Ghosh,KLD2,KLD+PR) 2694.0 ([13]) 9689.0 (KLD,KLD+PR) 20743.0 (KLD+PR) 35881.0 ([13],Ghosh,KLD2) 4658.0 ([13],KLD+PR) 16956.0 ([13]) 36315.0 ([13]) 62487.0 (KLD+PR) 7131.0 ([13]) 26254.0 ([13],KLD+PR) 56572.0 ([13],Ghosh,KLD,KLD2,KLD+PR) 97344.0 ([13],KLD+PR) In two cases (n = 300, m = 90 and n = 400, m = 160), KLD+PR achieved a new best solution which were not previously found by any other strategy. Table 13 report the CPU times observed for the execution of the same instances. For each GRASP heuristic, the average time for the three executions is reported. Heuristic GRASP with path-relinking, which achieved better quality solutions, demanded more time than the pure GRASP heuristic for all instances. Next, the same analysis performed in Subsection 2.3 is conducted here to evaluate the empirical probability distributions for the time to achieve a target value. We choose three instances from set D with parameters: (n = 200,m = 20), (n = 200,m = 80), and (n = 300,m = 120). Each GRASP heuristic was executed until a solution was found with a greater or equal cost compared to a target value. Two target values were used for each instance. Empirical probability distributions for the time to achieve a target value are plotted in Fig(s). 6, 7, and 8. We can observe that KLD+PR behaves better than KLD. KLD+PR algorithm presented a higher probability of finding a target solution in a given amount of time than KLD in all tests. Fig. 6 (target=1244) shows that the probability of finding the target value in 1 second is approximately 6% for KLD, while for KLD+PR is 75%. Considering 10 seconds, we have the probabilities 50% and 100% respectively. Also in Fig. 8 (target=35873), we observe that the 20 Table 13. Computational time for KLD and KLD+PR for instances from set D n 100 100 100 100 200 200 200 200 300 300 300 300 400 400 400 400 500 500 500 500 m 10 20 30 40 20 40 60 80 30 60 90 120 40 80 120 160 50 100 150 200 KLD 10.4 28.2 57.4 89.2 71.7 318.3 626.6 804.3 412.7 1690.5 3260.6 4453.0 1208.7 5266.5 10097.5 15750.3 2934.9 11921.7 29051.4 39043.3 KLD+PR 11.0 28.4 57.7 89.7 77.4 321.3 743.6 1122.8 413.7 1842.7 4183.1 7013.3 1300.6 5957.1 13850.7 20456.7 3072.5 15857.2 34661.4 52497.4 probability of finding the target value in 9 seconds is approximately 60% for KLD and for KLD+PR is 80%. 5 Concluding Remarks This paper presented some versions of GRASP heuristics to solve the maximum diversity problem (MDP). The main goal of this work was to analyse the influence of inserting the path-relinking technique to a pure GRASP. Experimental results show that the version which uses path-relinking technique (KLD+PR) significantly improves the average and the best quality of solutions generated by an important pure GRASP approach proposed in the literature (KLD). Our experiments also show that path-relinking speeds up convergence to target values. References 1. Aiex, R. M., Resende, M. G. C., Ribeiro, C. C.: Probability distribution of solution time in GRASP: An experimental investigation, Journal of Heuristics 8 (2002), 343–373. 2. Feo, T. A., Resende, M. G. C.: Greedy randomized adaptive search procedures, Journal of Global Optimization 6 (1995), 109–133. 21 1 n=200 m=20 Target: 1240 0.8 Probability 0.6 0.4 0.2 KLD KLD+PR 0 0.1 1 10 100 Time (sec) 1 0.8 n=200 m=20 Target: 1244 Probability 0.6 0.4 0.2 KLD KLD+PR 0 0.01 0.1 1 10 100 Time (sec) Fig. 6. Comparison of GRASP heuristics for the instance n = 200, m = 20 with targets values 1240 and 1244 22 1 n=200 m=80 Target: 16171 0.8 Probability 0.6 0.4 0.2 KLD KLD+PR 0 0.1 1 10 100 Time (sec) 1 n=200 m=80 Target: 16210 0.8 Probability 0.6 0.4 0.2 KLD KLD+PR 0 1 10 100 Time (sec) Fig. 7. Comparison of GRASP heuristics for the instance n = 200, m = 80 with targets values 16171 and 16210 23 1 n=300 m=120 Target: 35855 0.8 Probability 0.6 0.4 0.2 KLD KLD+PR 0 1 10 100 Time (sec) 1 n=300 m=120 Target: 35873 0.8 Probability 0.6 0.4 0.2 KLD KLD+PR 0 1 10 100 Time (sec) Fig. 8. Comparison of GRASP heuristics for the instance n = 300, m = 90 with targets values 35855 and 35873 24 3. Ghosh, J. B.: Computational aspects of the maximum diversity problem, Operations Research Letters 19 (1996), 175–181. 4. Glover, F., Laguna, M., Marti, R.: Fundamentals of scatter search and pathrelinking, Control and Cybernetics 19 (1977), 653–684. 5. Glover, F., Hersh, G., McMillan, C.: Selecting subsets of maximum diversity, MS/IS Report No. 77-9, University of Colorado at Boulder, (1977). 6. Glover, F., Kuo, C-C., Dhir, K. S.: Integer programming and heuristic approaches to the minimum diversity problem, Journal of Business and Management 4 (1996), 93–111. 7. Kochenberger, G., Glover, F.: Diversity data mining, Working Paper, The University of Mississipi, (1999). 8. Laguna, M., Mart, R.: GRASP and path relinking for 2-layer straight line crossing minimization, INFORMS Journal on Computing 11 (1999), 44–52. 9. Prais, M., Ribeiro, C. C.: Reactive GRASP: An application to a matrix decomposition problem in TDMA traffic assignment, INFORMS Journal on Computing 12 (2000), 164–176. 10. Resende, M. G. C., Ribeiro, C. C.: Greedy randomized adaptive search procedures, Handbook of Metaheuristics (F. Glover and G. Kochenberger editors), (2003), 219– 249. 11. Resende, M. G. C., Ribeiro, C. C.: GRASP with path-relinking: Recent advances and applications, Metaheuristics: Progress as Real Problem Solvers (T. Ibaraki et al. editors), (2005), 29–63. 12. Ribeiro, C. C., Uchoa, E., Werneck, R. F.: A hybrid GRASP with perturbations for the Steiner problem in graphs, INFORMS Journal on Computing 14 (2002), 228-246. 13. Silva, G. C., Ochi, L. S., Martins, S. L.: Experimental comparison of greedy randomized adaptive search procedures for the maximum diversity problem, Lecture Notes on Computer Science 3059 (2004), 498–512. 14. Weitz, R., Lakshminarayanan, S.: An empirical comparison of heuristic methods for creating maximally diverse groups, Journal of the Operational Research Society 49 (1998), 635–646.