ARTICLE IN PRESS
JID: EOR
[m5G;October 28, 2014;10:23]
European Journal of Operational Research 000 (2014) 1–7
Contents lists available at ScienceDirect
European Journal of Operational Research journal homepage: www.elsevier.com/locate/ejor
Cross entropy for multiobjective combinatorial optimization problems with linear relaxations Rafael Caballeroa, Alfredo G. Hernández-Díazb, Manuel Lagunac, Julián Molinaa,∗ a
Department of Applied Economics (Mathematics), University of Malaga, Malaga, Spain Department of Economics, Quantitative Methods and Economic History, Pablo de Olavide University, Seville, Spain c Leeds School of Business, University of Colorado, Boulder, CO 80309-0419, USA b
a r t i c l e
i n f o
Article history: Received 18 February 2013 Accepted 29 July 2014 Available online xxx Keywords: Multiobjective combinatorial optimization Cross entropy EMO Linear relaxation
a b s t r a c t While the cross entropy methodology has been applied to a fair number of combinatorial optimization problems with a single objective, its adaptation to multiobjective optimization has been sporadic. We develop a multiobjective optimization cross entropy (MOCE) procedure for combinatorial optimization problems for which there is a linear relaxation (obtained by ignoring the integrality restrictions) that can be solved in polynomial time. The presence of a relaxation that can be solved with modest computational time is an important characteristic of the problems under consideration because our procedure is designed to exploit relaxed solutions. This is done with a strategy that divides the objective function space into areas and a mechanism that seeds these areas with relaxed solutions. Our main interest is to tackle problems whose solutions are represented by binary variables and whose relaxation is a linear program. Our tests with multiobjective knapsack problems and multiobjective assignment problems show the merit of the proposed procedure. © 2014 Elsevier B.V. All rights reserved.
1. Introduction Cross entropy (CE) has a relatively short history in the realm of optimization methodologies. Rubinstein (1997) developed CE as a method for estimating probabilities of rare event in complex stochastic networks. Two years later, CE was applied for the first time in the context of combinatorial optimization, triggering a number of other operation research applications that are now well documented (de Boer, Kroese, Mannor, & Rubinstein, 2005a). A fairly comprehensive list of CE applications and tutorial materials is found in http://www.cemethod.org. In its most basic form, CE consists of the repeated execution of the following two steps (de Boer, Kroese, Mannor, & Rubinstein, 2005b): 1. Generate a random sample from a pre-specified probability distribution function. 2. Use the sample to modify the parameters of the probability distribution in order to produce a “better” sample in the next iteration.
∗
Corresponding author. E-mail addresses:
[email protected] (R. Caballero),
[email protected] (A.G. Hernández-Díaz),
[email protected] (M. Laguna),
[email protected] (J. Molina).
Performance of CE implementations varies according to the rules (and related parameter values) used to solve the so-called associated stochastic problem (ASP), which is the problem of estimating the probability that the objective function of a vector of random variables that follow a probability distribution function with parameters v exceeds a given value γ . The problem is solved by generating a sequence of (v, γ ) values that are adaptively updated from iteration to iteration. If successful, the search converges to a small neighborhood around or precisely to the optimal values denoted by (v∗ , γ ∗ ). The CE implementations in the literature are variants of this process, where most of the changes consist of the specific rules to obtain (vt+1 , γ t+1 ) in iteration t + 1 from the values (vt , γ t ) in iteration t. The exception is, perhaps, the adaptation introduced by Laguna, Duarte, and Martí (2009), where the method is hybridized with the addition of local search. The main difference between CE and the so-called estimation of distribution algorithms (EDAs) is the assumption of independence of the decision variables. In its simplest form, CE attempts to estimate the probability that at optimality a variable takes on the value of 1, and it does so by treating each variable independently. This is can be a disadvantage if variables are strongly correlated, but, on the other hand, it makes CE very easy to implement and to adjust. Applications to single objective combinatorial optimization problems have shown that the methodology can be quite effective even when not directly addressing the presence of covariance.
http://dx.doi.org/10.1016/j.ejor.2014.07.046 0377-2217/© 2014 Elsevier B.V. All rights reserved.
Please cite this article as: R. Caballero et al., Cross entropy for multiobjective combinatorial optimization problems with linear relaxations, European Journal of Operational Research (2014), http://dx.doi.org/10.1016/j.ejor.2014.07.046
JID: EOR 2
ARTICLE IN PRESS
[m5G;October 28, 2014;10:23]
R. Caballero et al. / European Journal of Operational Research 000 (2014) 1–7
We are aware of only three CE applications in the context of multiobjective optimization. The first one is the work by Ünveren and Acan (2007) that tackles the problem of finding efficient frontiers for problems with multiple multimodal objective functions of continuous variables. They “introduce the notion of clustered nondominated solutions on the Pareto front to adapt the probability distribution parameters” within the CE. In single-objective optimization problem, the updating of the (v, γ ) values is based on the best (elite) solutions found in the random sample of the current iteration. In multiobjective optimization, however, nondominated solutions are found along the Pareto front, rendering the use of a common set of (v, γ ) values impractical. This is why Ünveren and Acan (2007) divide the set of nondominated solutions into clusters and essentially execute one CE procedure in each cluster. The number of clusters is set to be the number of objective functions in the problem plus one. The updating of each set of (v, γ ) values is cluster-dependent and therefore each CE is responsible for finding the best possible set of nondominated solutions in its assigned region of the Pareto front (represented by the cluster). The CE method applied to each cluster developed by Ünveren and Acan (2007) is practically identical to the CE Algorithm for Continuous Optimization described by Kroese, Porotsky and Rubinstein (2006). In particular, the procedure starts with an initial set of values for the mean and variance of the decision variables. A sample is taken from Normal distributions with the appropriate parameters for each variable. The mean and variances are updated using the elite solutions in the sample. Finally, a small (epsilon) value is used to compare the maximum variance and determine whether or not all variable values in the sample have converged to their corresponding mean. The method is compared to eight existing methods from the literature using seven well-known optimization problems with multiple multimodal functions, six with two objective functions and one with three. The key search parameters are set to a sample of 5000 solutions (per CE) and 1000 iterations. This means that for a biobjective problem, the procedure evaluates 5000 × 1000 × 3 = 15 million solutions. A total of 20 million solutions are evaluated in problems with three objective functions. The authors conclude that their procedure “performs better than its competitors for most of the tests cases.” Connor (2008) applies the procedure by Ünveren and Acan (2007) to a real-world problem with two objective functions. In this application, ten clusters instead of the recommended three are employed. The second multiobjective CE method that we have found in the literature is due to Perelman, Ostfeld, and Salomons (2008). Within the field of water resource management, Perelman, Ostfeld and Salomons address the problem of designing a water distribution system. Essentially, the problem consists of choosing where to locate pipes in a water distribution network and also select the size of the pipes. The single-objective version of the problem attempts to find a design that minimizes the overall construction cost while satisfying minimum water-pressure requirements at demandnodes in the network. This is a classical problem in civil engineering that dates back to the work by Schaake and Lai (1969) where a linear programming model was used to find an approximate solution to the New York Tunnels Water Distribution System. Since then, numerous procedures have been applied to the single-objective function problem, including metaheuristics (e.g., Lippai, Heaney, & Laguna, 1999). The multiobjective optimization version consists of making the water-pressure requirements a second objective as opposed to a constraint. This results in two extreme points, the zero-cost solutions, where the demands are not satisfied, and the highest cost solution that corresponds to satisfying all the demands (at a minimum cost). The first point is obtained trivially. The second requires the solution of the single-objective problem where the demand requirements are treated as hard constraints (as explained above). The CE implemen-
tation by Perelman et al. (2008) for the multiobjective optimization problem uses binary variables to select a pipe size in each location. Therefore, corresponding to each location, there is a binary vector of m variables, where m is the number of available pipe diameters (including the “zero” or “do nothing” alternative). A value of 1 indicates the selection of a particular pipe diameter. In order to identify the elite solutions that will serve as the basis for updating the CE parameters, the authors invoke the concept of ranking, introduced by Fonseca and Fleming (1995). The ranking method assigns a rank value of one to all nondominated solutions. Dominated solutions are assigned higher rank values corresponding to the number of solutions that dominate them. In this sense, the rank method incorporates density information corresponding to different regions of the Pareto front. The rank values are used to determine the set of elite solutions in the sample, which in turn become the basis for updating the (v, γ ) values. Bekker and Aldrich (2011) propose a multiobjective CE for continuous problems with characteristics that are similar to the one developed by Perelman et al. (2008) for discrete optimization. The ranking of solutions is employed as a mechanism for selecting the set of elite solutions from the sample to update the (v, γ ) values. In particular, the elite solutions become the basis for constructing empirical probability distribution functions, one for each variable in the problem. The size of the sample is set between 30 and 50 times the number of objective functions in the problem. In the computational experiments with eight standard problems from the literature the number of function evaluations is in the order of 10,000, resulting in high-quality fronts and computational times not exceeding 25 seconds. No comparisons against other methods are provided but convergence and dispersion metrics are used to assess the quality of the solutions. One of the main differences between Bekker and Aldrich (2011) and Perelman et al. (2008) is the selection of the “elite” solutions from the sample. While Perelman et al. (2008) use a fixed size throughout the search, Bekker and Aldrich (2011) use a size that grows over time and that is controlled by the ranking value (e.g., all solutions with ranking less than or equal to two may be selected).
2. Multiobjective combinatorial optimization (MOCO) problems The MOCO problems that we are interested in solving have the following form:
Maximize f1 (x), f2 (x), . . . , fm Subject to x ∈ X where X represents the set of all feasible solutions and may constraint x to take on integer values. As customary in multiobjective optimization, solution x is said to dominate solution y if fi (x) ࣙ fi (y) for all i and fi (x) > fi (y) for at least one i. We focus on two specific MOCO problems in order to test our methodology: the multiobjective knapsack problem and the multiobjective assignment problem. The multiobjective knapsack problem (MOKP) consists of selecting a subset of items from a set of n items in order to maximize the utility (profit) of m different knapsacks without violating their individual capacities. The profit for knapsack i associated with selecting item j is given by pij and the corresponding weight is denoted by wij . The capacity of the ith knapsack is represented by ci . The problem may be formulated as follows:
Maximize fi (x) = Subject to
n
n
pij xj
i = 1, . . . ,m
j=1
wij xj ≤ ci
i = 1, . . . ,m
j=1
x ∈ {0,1}
Please cite this article as: R. Caballero et al., Cross entropy for multiobjective combinatorial optimization problems with linear relaxations, European Journal of Operational Research (2014), http://dx.doi.org/10.1016/j.ejor.2014.07.046
ARTICLE IN PRESS
JID: EOR
[m5G;October 28, 2014;10:23]
R. Caballero et al. / European Journal of Operational Research 000 (2014) 1–7
The MOKP was first approached more than ten years ago by Zitzler and Thiele (1999). They developed the Strength Pareto Evolutionary Algorithm (SPEA) and compared it to four other procedures while introducing nine MOKP instances to the literature. These instances, with number of items ranging from 250 to 750 and number of objective functions ranging from two to four, have become the main MOKP benchmark for subsequent studies. For instance, Jaszkiewicz (2002) utilizes this set of problems to compare his Multiobjective Genetic Local Search (MOGLS) procedure with SPEA and M-PAES (Knowles & Corne, 2000). The most recent heuristic solution method for the MOKP is due to Alves and Almeida (2007) and the most recent exact procedure is due to Bazgan, Hugot, and Vanderpooten (2009). The multiobjective linear assignment problem (MOLAP) consists of assigning n items to n locations in order to minimize the aggregate cost values of m different objectives. The cost—associated with the kth objective—of assigning item i to location j is given by cijk . Each item must be assigned to a location and each location may receive only one item. Mathematically, the problem may be formulated as follows:
Minimize fk (x) = Subject to
n j=1 n
n n
cijk xij
k = 1, . . . , m
i=1 j=1
xij = 1 i = 1, . . . , n xij = 1 j = 1, . . . , n
i=1
x ∈ {0,1} The seminal work by Charnes, Cooper, Niehaus, and Stedry (1969) is the starting point of extensions to the linear assignment problem in order to consider multiple objectives. The work focuses on the assignment of personnel to jobs and the extension from a single criterion for measuring the merit of the assignment to a vector optimization. The models are based on the notion that jobs require a given amount of a number of attributes possessed by the individuals to be assigned to perform the jobs. Then, each assignment of individual to jobs achieves a measure of fulfillment of the required attributes. The approach results in multiobjective problems with a number of objectives given by the product of the number of jobs and attributes. The formulations, however, focus on minimax solutions. Pollatscheck (1976) formulates a bi-objective assignment problem of the form shown above for assigning personnel to jobs, where one of the criteria is the cost of moving an individual to a post and the other is the reliability of the individual in the assigned post. The model seeks to minimize the cost and maximize the reliability. The problem is approached with an interactive procedure that makes one assignment at a time and in which the decision maker chooses first- and secondbest assignments from a given set. White (1984) studied a special MOLAP that is similar to the one proposed by Charnes et al. (1969). Unlike Charnes et al. (1969) and White (1984) is interested in producing efficient solutions instead of concentrating on minimax solutions. One of the main contributions of White’s work consists of demonstrating that “whereas, for many integer problems, the standard scalar weighting factor approach will not produce all the efficient solutions, in [the case of the special MOLAP] it will.” Tuyttens, Teghem, Fortemps, and Van Nieuwenhuyze, (2000) develop an exact procedure and apply MOSA (multiobjective simulated annealing) to the bi-objective linear assignment problem. They show, by example, that the MOLAP has a set of non-supported efficient solutions and therefore it is not possible to find efficient frontiers simply by solving parameterized single-objective problems. For
3
our computational tests, we employ the 10 problem instances with n = 5, 10, . . . , 50 generated by Tuyttens et al. (2000). 3. Multiobjective cross entropy (MOCE) method Before describing our MOCE implementation, it is useful to provide some details of a basic cross entropy procedure for optimization problems with a single objective function S(x) and for which the solutions are represented by a binary vector x of size n. In this situation, the CE methodology prescribes the use of independent Bernoulli distributions to sample the solution space. The parameters associated with the Bernoulli distributions at iteration t are given by vt = (vt (1), . . . , vt (n)), where vt (j) is the probability that xj takes on the value of 1 at iteration t. A sample of size N is generated at each iteration and the elite set consists of the best ρ N solutions in the sample, where ρ is a CE parameter between 0 and 1. This parameter value determines the solution of the ASP, because in the ordered set of elite solutions, γ t is the objective function corresponding to the last solution, that is γt = S(x(ρ N)). Also, the best estimate of γ ∗ is the first solution in the ordered set of elite solutions, that is γ ∗ ࣈ S(x(1) ). A pseudo code of the basic CE procedure for a single-objective binary optimization problem is shown in Fig. 1. The procedure in Fig. 1 starts by setting the probability vt (j) that xj = 1 to 0.5 for all j (line 1). Using these probabilities, the sample of size N is taken in line 2 (inside the while-loop). The sample is ordered (line 3) and the new vt+1 (j) estimates are calculated (line 4). To avoid drastic changes in the estimation of v, the updating is typically done with a smoothing function of the type shown in line 5 of the pseudocode of Fig. 1, where α takes on a value between 0 and 1. The procedure terminates after a pre-specified number of iterations. The procedure operates with the instantiation of 4 input parameters: maxiter, N, ρ and α . We adapt this basic scheme to search for a set of points that estimate the optimal Pareto front instead of searching for a point that optimizes a single objective function. As shown below, the ASP is solved by employing a set of v vectors of size n instead of a single one. The set of v vectors is such that the points associated with solving the APS provide a good approximation of the Pareto front. As indicated above, Ünveren and Acan (2007) solve this problem by dividing the set of nondominated solutions into clusters and executing independent CE procedures in each cluster. This approach, however, has the following drawbacks: 1. There is lack of control regarding the membership of points in the clusters. That is, the same point may belong to different clusters in different iterations, causing the procedure to oscillate arbitrarily. 2. Clusters have different sizes resulting in the overestimation of some areas of the Pareto front while other areas are underestimated. 3. When clusters grow in size, the process may become computationally impractical.
Fig. 1. Basic CE procedure for binary optimization problems with a single objective.
Please cite this article as: R. Caballero et al., Cross entropy for multiobjective combinatorial optimization problems with linear relaxations, European Journal of Operational Research (2014), http://dx.doi.org/10.1016/j.ejor.2014.07.046
ARTICLE IN PRESS
JID: EOR 4
[m5G;October 28, 2014;10:23]
R. Caballero et al. / European Journal of Operational Research 000 (2014) 1–7
4. Due to lack of control mentioned above, convergence becomes a problem. Experiments show that the procedure suggested by Ünveren and Acan (2007) required 15 million objective function evaluations to compete with solutions that state-of-the-art methods produce with 10 thousand evaluations on the same problem instances. (We should point out, however, that in some cases, Ünveren and Acan’s procedure is competitive in terms of computational time with alternative approaches.) In order to address these issues, we propose to divide the Pareto front into uniform areas and then solve independent ASP’s. We create an ɛ-grid on the objective function space by dividing the range of each objective function fi into uniform segments of size ɛi , as illustrated in Fig. 2 in the two-dimensional case. The basic CE procedure is then executed in each area of the ɛgrid. We solve the single-objective linear relaxation of the problem to estimate the range rj of each objective function. We then select a value M to represent the number of segments for each objective and calculate:
εj =
Fig. 3. Fractional and integer solutions within the same area.
rj M
Therefore, a vat vector is associated with each area a and updated at each iteration t. The sample size N is the same for each area, but the updating of the vat vector considers only the nondominated solutions in the area. Because some areas could contain too few points, particularly at the extremes of the function ranges, we impose bounds on the value of any element of the vat vectors:
lb ≤ vat (j) ≤ 1 − lb The bounds avoid premature convergence to a single nondominated point in each area. Even with this bounding mechanism, areas may lack a sufficient number of points to yield a good estimate of the Pareto front. Given that linear relaxations of integer programs are typically easy to solve, we use relaxed solutions to enrich the estimation of vat in each area. Fractional solutions, however, will typically dominate integer solutions and therefore we need a mechanism to allow both to coexist within an area (see Fig. 3). The mechanism consists of expanding the dimensionality of the multiobjective problem by adding an objective m + 1. The form of the additional objective is as follows:
fm+1 (x) =
1 if x is integer 0 otherwise
This additional objective guarantees that, for a maximization problem, an integer solution is not dominated by any fractional solution. In the extended multiobjective problem, both integer and fractional
Fig. 4. MOCE for binary combinatorial optimization problems.
solutions can coexist in the areas formed by the ɛ-grid. The procedure is summarized in Fig. 4. The process starts with the solution of the relaxed single-objective optimization problems which are used to estimate the ranges of these functions and then build the ɛ-grid. The areas in the grid are populated with fractional solutions corresponding to the relaxed multiobjective combinatorial optimization problem (see Step 2). Then, the v vectors are initialized with the average value (for each variable) of the nondominated fractional solutions in each area (Step 3). The procedure then executes Steps 4–8 maxiter times. The iterative part of the procedure consists of sampling each area and then updating the set of nondominated solutions. The v vectors are updated and adjusted if needed in order to fall within the specified bounds. Note that, as pointed out earlier, the set of nondominated solutions in each area includes both fractional and integer solutions and that the entire set is used to update the v values. At the end of the process, the best estimate of the Pareto front consists of the union of the set nondominated integer solutions in each area. 4. Computational experiments Our computational testing includes the two problem classes described in Section 2. We point out that our MOCE implementation does not include a repair mechanism for constraint handling. When applied to MOKP, the procured discards infeasible solutions. By design, no infeasible solutions are generated when tackling MOLAP instances. We start with a comparison of the performance of MOCE and the state-of-the-art in heuristic solution procedures for the MOKP: •
Fig. 2. Sample ε -grid.
SPEA by Zitzler and Thiele (1999). This is the main reference method in the literature for MOKP. Most MOKP publications use SPEA as the standard and the test problems introduced by Zitzler and Thiele (1999) as the benchmark.
Please cite this article as: R. Caballero et al., Cross entropy for multiobjective combinatorial optimization problems with linear relaxations, European Journal of Operational Research (2014), http://dx.doi.org/10.1016/j.ejor.2014.07.046
JID: EOR
ARTICLE IN PRESS
[m5G;October 28, 2014;10:23]
R. Caballero et al. / European Journal of Operational Research 000 (2014) 1–7 •
•
MOGLS by Jaszkiewicz (2002). Experimental testing has shown that MOGLS is able to generate much better approximations to the nondominated sets than SPEA, when both methods are set to employ the same number of function evaluations. MOTGA by Alves and Almeida (2007). This is a Multiobjective Genetic Algorithm based on the Tchebycheff scalarizing function that performs several stages, each one intended for identifying nondominated solutions in a different part of the Pareto front. This is the latest heuristic method in the literature and has demonstrated to outperform both MOGLS and SPEA over the same set of benchmark problems. MOTGA performed more iterations but consumed less computational time than both SPEA and MOGLS to find better approximations of the Pareto front. In fact, Alves and Almeida (2007) do not report the total number of evaluations and it is not possible to calculate them from the information in the article because the number of evaluations performed by the local operators is not reported.
As done in previous studies, we use the nine MOKP instances and the limit on the function evaluations reported in Zitzler and Thiele (1999) and Jaszkiewicz (2002) to compare our work (MOCE) with the outcomes produced by SPEA, MOGLS and MOTGA. The performance measures that we employ are: •
•
•
•
Number of points — This refers to the ability of finding efficient points. One of the strong points of MOGLS consists of finding more potentially nondominated solutions than its competitors (including MOTGA). Hence, we include this measure in order to assess the merit of MOCE as compared to MOGLS in this dimension. SSC — This metric suggested by Zitzler and Thiele (1999) measures the size of the space covered (SSC). In other words, SSC measures the volume of the dominated points. Hence, the larger the SSC value the better. This measure is considered standard in the literature of heuristic approaches to multiobjective optimization problems. k-distance — This density estimation technique, introduced by Zitzler et al. (2001), is based on the kth nearest neighbor method of Silverman (1986). The metric is the distance to the kth nearest efficient point. We use k = 5 and calculate the mean of the k-distance values for all the points in the final approximations. The k-distance value is such that the smaller the better in terms of frontier density. C(A,B) — This is known as the coverage of two sets measure (Zitzler & Thiele, 1999). C(A,B) represents the proportion of points in the estimated efficient frontier B that are dominated by the efficient points in the estimated frontier A. This measure is also standard in the multiobjective optimization literature.
Tables 1–4 show the metrics associated with the four procedures in the comparison set. Table 5 contains the best objective function values found by each of the methods when applied to each of the problems in the test set. The parameter settings for SPEA, MOGLS and MOTGA are the ones suggested in the corresponding publications. The settings for MOCE are M = 10, lb = 0.001, and N = 100. In addition, 100 points are used as seeds in the execution of Step 2 of the procedure. Table 1 shows that MOCE is able to find more nondominated points in each of the test problem than any of the other methods, including MOGLS, which is precisely this method’s strength. Regarding the SSC measure, Table 2 shows that MOCE’s average is better than the competing approaches. The difference with MOTGA is marginal but it is significant when compared to SPEA and MOGLS in all nine problems. In terms of the k-distance measure (reported in Table 3), MOCE is shown to be clearly superior to the rest of methods, achieving the best values for all instances. This means that MOCE is capable of finding approximations of the efficient frontier that have both a
5
Table 1 Number of nondominated solutions. Problem
Objectives
Variables
SPEA
MOGLS
MOTGA
MOCE
1 2 3 4 5 6 7 8 9
2 3 4 2 3 4 2 3 4
250 250 250 500 500 500 750 750 750
60 330 1067 37 410 1054 250 272 1135
99 843 2045 131 1292 2904 162 1292 3393
104 507 1085 166 918 1952 248 1519 3185
145 1674 5710 282 2810 7143 298 2946 7866
Problem
Objectives
Variables
SPEA
MOGLS
MOTGA
MOCE
1 2 3 4 5 6 7 8 9 Average
2 3 4 2 3 4 2 3 4
250 250 250 500 500 500 750 750 750
0.840 0.675 0.581 0.765 0.639 0.547 0.776 0.763 0.538 0.681
0.895 0.753 0.649 0.867 0.748 0.637 0.829 0.510 0.636 0.725
0.909 0.768 0.660 0.889 0.775 0.656 0.852 0.924 0.663 0.788
0.910 0.767 0.666 0.890 0.773 0.652 0.903 0.928 0.650 0.793
Table 2 SSC values.
Table 3 k-Distance measure. Problem
Objectives
Variables
SPEA
MOGLS
MOTGA
MOCE
1 2 3 4 5 6 7 8 9 Average
2 3 4 2 3 4 2 3 4
250 250 250 500 500 500 750 750 750
0.030 0.034 0.043 0.030 0.025 0.036 0.006 0.012 0.027 0.027
0.042 0.041 0.060 0.031 0.033 0.051 0.023 0.011 0.045 0.038
0.035 0.046 0.061 0.023 0.034 0.046 0.014 0.012 0.040 0.034
0.027 0.019 0.020 0.012 0.010 0.011 0.008 0.004 0.008 0.013
Table 4 C(A,B) measure.
SPEA MOGLS MOTGA MOCE
SPEA
MOGLS
MOTGA
MOCE
— 0.779 0.999 0.998
0.127 — 0.789 0.817
0.000 0.009 — 0.340
0.000 0.001 0.068 —
large number of points and a good density, as indicated by the small k-distance values. The C(A,B) measure enables us to make direct comparisons of one approximation of the efficient frontier against another in terms of dominance. Table 4 shows the average C(A,B) values over the entire set of test problems. The values in this table show that the nondominated points generated by MOCE tend to dominate those generated by other methods. That is, C(MOCE,−) > C(−,MOCE). In fact, none of the points generated by SPEA dominates any point generated by MOCE, and only 0.1% of the points generated by MOCE are dominated by points generated by MOGLS, for example. MOTGA—the most competitive method in the experiment—is able to dominate only less than 7% of the points generated by MOCE, while, on the other hand, MOCE dominates 34% of the points generated by MOTGA. Additional insight on the performance of MOCE is gained by examining the best values found for each objective function in the set of test problems, as shown in Table 5. MOCE is able to obtain better
Please cite this article as: R. Caballero et al., Cross entropy for multiobjective combinatorial optimization problems with linear relaxations, European Journal of Operational Research (2014), http://dx.doi.org/10.1016/j.ejor.2014.07.046
ARTICLE IN PRESS
JID: EOR 6
[m5G;October 28, 2014;10:23]
R. Caballero et al. / European Journal of Operational Research 000 (2014) 1–7
Table 5 Best objective function values found by each method in each test problem. Problem
Function
1
f1 f2
2
MOGLS
MOTGA
MOCE
9395 9773
9872 10,059
9880 10,086
9901 10,102
f1 f2 f3
9370 9209 9057
9829 10,045 9728
9829 9962 9677
9960 10,062 9760
3
f1 f2 f3 f4
8958 8978 9091 8298
9684 9855 10,032 8965
9636 9812 9931 8953
9768 9897 10,091 9124
4
f1 f2
19,001 19,383
19,959 20,436
20,066 20,471
20,088 20,480
5
f1 f2 f3
17,923 18,322 18,633
19,435 19,675 20,347
19,488 19,864 20,185
19,696 20,062 20,415
6
f1 f2 f3 f4
17,383 17,425 17,481 17,169
19,660 19,355 19,672 19,436
19,589 19,195 19,588 19,239
19,814 19,551 19,787 19,569
7
f1 f2
28,758 28,621
29,828 29,969
30,060 29,988
30,663 30,208
8
f1 f2 f3
28,002 26,842 27,360
19,435 19,675 20,347
30,957 29,402 30,143
31,180 29,752 30,505
9
f1 f2 f3 f4
25,634 25,297 27,022 25,943
29,789 28,831 30,392 29,718
29,203 28,781 29,785 29,438
29,845 29,324 30,442 30,010
18,086
18,668
19,748
20,002
Average
SPEA
values, when compared to the other three procedures, for each of the 27 objective functions corresponding to the nine problems in the test set. For problems with two objectives, a drawing of the Pareto fronts provides visual evidence of the relative quality of the approximations. Fig. 5 shows, for the bi-objective problems with 750 items, not only the approximations found by each of the four methods but also the Pareto front found by relaxing the integrality constraints (labeled as LP). This front is an upper bound of the exact frontier associated with the integer problem. Note that the MOCE approximation is not only significantly better than the rest of the approximations but it is also almost overlapping the upper bound. The “seeding” step in our MOCE is a key element of the procedure (see Step 2 of Fig. 4). This step consists of initializing the grid with the fractional solutions obtained from solving the relaxed problem. We illustrate the importance of this step with the bi-objective MOKP instance with 500 items. We run MOCE under the following conditions regarding Step 2 of the procedure:
Fig. 5. Approximations of the Pareto front for problems with 2 objectives and 750 items.
Fig. 6. Effect of seed points on the quality of the approximation.
•
•
•
•
0 Points. Step 2 is not performed and therefore no seeds are used to initiate the search. 2 Points. The optimal solutions to the single-objective relaxed problems are used as seeds. These are the two points at the extremes of the upper-bound frontier. 3 Points. In addition to the two extreme points, a third point is chosen randomly from the relaxed solutions. 100 Points. Step 2 is performed with 100 points as the set of seeds to initiate the search.
Fig. 6 shows the results of executing MOCE with the four options shown above. Clearly, the method cannot operate without the seed points. Note that when MOCE is run without executing Step 2, the procedure is only capable of finding one dominated solution, which is near the intersection of the axis in Fig. 6. When two solutions are used, the method is able to find nondominated solutions in the extremes of the frontier. However, it is not able to find any points in the middle part of the frontier. The addition of a third point has a palpable effect in performance by allowing the procedure to create an approximation of the frontier. Note, however that the approximation is denser around the three seed points (i.e., the two extreme points and the middle, around f1 = 17,800 and f2 = 20,000). As expected, the approximation considerably improves when 100 seed points are used. So it is clear how the method strongly depend on the solutions of the linear relaxation, as mentioned. Although it is worth to note that solutions of the linear relaxations are so easy to obtain and most of the most important combinatorial problems have a linear relaxation. We now turn or attention to the multiobjective assignment problem. We are interested in testing how MOCE, without any structural changes or additional parameter tuning, performs on a different problem in the class of binary constrained optimization. One of the strengths of a CE algorithm is low dependence on the parameter setting, so we intend to test the same setting in a completely different problem. In particular, we use MOCE to find an approximation of the efficient frontier of the 10 bi-objective MOLAP instances generated by Tuyttens et al., 2000. In order to apply a similar computational effort to these problems than the one used in the experiments with the MOKP instances, we stop the procedure after a given number of function evaluations. MOKP instances with 750 variables and two objectives consumed 125,000 function evaluations, or equivalently, 166 evaluations per variable. Using a linear relationship, when applied to a MOLAP instance, we allow MOCE to perform 166 × n2 = 203,350 function evaluations (where n is the number of tasks to assign). Since we obtained similar results for all 10 problem instances, instead of producing additional tables with all the metrics shown above, we summarize our findings with observations related to one of the larger instances in the test set, the problem instance with n = 35. Fig. 7 shows the approximation of the Pareto front that MOCE found for the MOLAP instance with n = 35. The exact frontier shown in Fig. 7 is the one found with the exact method of Tuyttens et al., 2000.
Please cite this article as: R. Caballero et al., Cross entropy for multiobjective combinatorial optimization problems with linear relaxations, European Journal of Operational Research (2014), http://dx.doi.org/10.1016/j.ejor.2014.07.046
JID: EOR
ARTICLE IN PRESS
[m5G;October 28, 2014;10:23]
R. Caballero et al. / European Journal of Operational Research 000 (2014) 1–7
7
problems have linear programming relaxations. The application of MOCE to two different problems, MOKP and MOLAP, show that the method is capable of producing high-quality approximations of the efficient frontier without the need to include search strategies designed to exploit problem-specific knowledge or even specific parameter setting. Acknowledgment This work was partially supported by the Government of Andalucía (Project P10-TIC-6618) in Spain. Fig. 7. Approximation of the Pareto front for an instance of MOLAP with 1225 variables and two objectives.
The LP points are the relaxed solutions that are generated as seeds for MOCE. The approximation of the efficient frontier depicted in Fig. 7 representative of the approximations that MOCE produces for all the MOLAP instances in the set. This figure shows how, without any changes, MOCE is capable of finding an excellent approximation of the Pareto front for an instance of a bi-objective assignment problem. In particular, the SSC calculation (i.e., the percentage of the dominated surface of the objective function space) associated with the Pareto fronts in Fig. 7 is such that the exact method has a value of 88.55% while MOCE’s value is 80.37%. We point out that none of the other methods used in the MOKP comparison has been applied to MOLAP instances because they include features specifically related to the knapsack structure of the problems that cannot be generalized to solve assignment problems. However, we did compare MOCE with the simulated annealing procedure (MOSA) developed by Tuyttens et al. (2000). We were able to determine that MOCE outperforms MOSA, even though MOSA is a specialized procedure for MOLAP instances. Tuyttens et al. (2000) measure the performance of MOSA with two distance metrics: D1 (the average distance between the approximation of the efficient frontier and the exact frontier) and D2 (the maximum distance between the two frontiers). For the problem shown in Fig. 7, Tuyttens et al. (2000) report values of D1 = 0.035 and D2 = 0.205 for a MOSA efficient frontier approximation with 47 points. The corresponding values for a MOCE efficient frontier approximation with 49 points are D1 = 0.023 and D2 = 0.099. We observed similar difference in the other 9 test instances in the set. The results obtained with both MOKP and MOLAP instances, without any customization, seem to indicate that MOCE is a flexible and robust procedure for multiobjective combinatorial problems that can be formulated with binary variables. 5. Conclusions We achieved our goal of developing a solution method for multiobjective combinatorial optimization problems based on the cross entropy methodology. In this work, we focused on problems whose solution representation is a vector of binary values and for which the relaxation results in a linear program. The results of our computational experiments are encouraging regarding the potential for expanding MOCE to other combinatorial optimization problems that have “easily” solvable linear relaxations. We illustrated the importance of the presence of relaxed points to initialize the search. However, we do not believe that the reliance on linear relaxation is a major shortcoming of the method because numerous relevant
References Alves, M. J., & Almeida, M. (2007). MOTGA: A multiobjective Tchebycheff based genetic algorithm for the multidimensional knapsack problem. Computers & Operations Research, 34(11), 3458–3470. Bazgan, C., Hugot, H., & Vanderpooten, D. (2009). Solving efficiently the 0–1 multiobjective knapsack problem. Computers and Operations Research, 36(1), 260–279. Bekker, J., & Aldrich, C. (2011). The cross entropy method in multi-objective optimisation: An assessment. European Journal of Operational Research, 211(1), 112–121. Charnes, A., Cooper, W. W., Niehaus, R. J., & Stedry, A. (1969). Static and dynamic assignment models with multiple objectives, and some remarks on organization design. Management Science, 15(8), B-365–B-375. Connor, J. D. (2008). Antenna array synthesis using the cross entropy method. PhD Thesis, Florida State University. de Boer, P. -T., Kroese, D. P., Mannor, S., & Rubinstein, R. Y. (2005a). The cross-entropy method for combinatorial optimization, rare event simulation and neural computation. Annals of Operations Research, 134. de Boer, P. -T., Kroese, D. P., Mannor, S., & Rubinstein, R. Y. (2005b). A tutorial of the cross entropy method. Annals of Operations Research, 134, 19–67. Fonseca, C. M., & Fleming, P. J. (1995). An overview of evolutionary algorithms in multiobjective optimization. Evolutionary Computation, 3(1), 1–16. Jaszkiewicz, A. (2002). On the performance of multiobjective genetic local search on the 0/1 knapsack problem—A comparative experiment. IEEE Transactions on Evolutionary Computation, 6(4), 402–412. Knowles, J. D. & Corne, D. W. (2000). M-PAES: A mimetic algorithm for multiobjective optimization. In Proceedings of the 2000 congress on evolutionary computation (Vol. 1, pp. 325–332). Piscataway, NJ: IEEE Press. Kroese, D. P., Porotsky, S., & Rubinstein, R. Y. (2006). The cross-entropy method for continuous multi-extremal optimization. Methodology and Computing in Applied Probability, 8(3), 383–407. Laguna, M., Duarte, A., & Martí, R. (2009). Hybridizing the cross entropy method: An application to the max-cut problem. Computers and Operations Research, 36(2), 487–498. Lippai, I., Heaney, J. P., & Laguna, M. (1999). Robust water system design with commercial intelligent search optimizers. Computing in Civil Engineering, 13(3), 135–143. Perelman, L., Ostfeld, A., & Salomons, E. (2008). Cross entropy multiobjective optimization for water distribution systems design. Water Resources Research, 44, W09413. doi:10.1029/2007WR006248. Pollatscheck, M. A. (1976). Personnel assignment by multiobjective programming. Mathematical Methods of Operations Research, 20(5), 161–170. Rubinstein, R. Y. (1997). Optimization of computer simulation models with rare events. European Journal of Operations Research, 99, 89–112. Schaake, J. C. & Lai, D. (1969). Linear programming and dynamic programming applications to water distribution network. Report No. 116. Department of Civil Engineering, Massachusetts Institute of Technology, Cambridge, Mass. Silverman, B. W. (1986). Density estimation for statistics and data analysis. London: Chapman and Hall. Tuyttens, D., Teghem, J., Fortemps, Ph., & Van Nieuwenhuyze, K. (2000). Performance of the MOSA method for the bicriteria assignment problem. Journal of Heuristics, 6, 295–310. Ünveren, A., & Acan, A. (2007). Milt-objective optimization with cross entropy method: stochastic learning with clustered Pareto fronts. IEEE Congress on Evolutionary Computation, 3065–3071. White, D. J. (1984). A special multi-objective assignment problem. Journal of the Operational Research Society, 35(8), 759–767. Zitzler, E., Laumanns, M., & Thiele, L. (2001). SPEA2: Improving the strength pareto evolutionary algorithm for multiobjective optimization. In K. C. Giannakoglou, D. T. Tsahalis, J. Periaux, & K. D. Papailiou (Eds.), Evolutionary methods for design, optimisation and control with application to industrial problems (EUROGEN 2001) (pp. 95–100). Zitzler, E., & Thiele, L. (1999). Multiobjective evolutionary algorithms: A comparative case study and the strength Pareto evolutionary algorithm. IEEE Transactions on Evolutionary Computation, 3, 257–271.
Please cite this article as: R. Caballero et al., Cross entropy for multiobjective combinatorial optimization problems with linear relaxations, European Journal of Operational Research (2014), http://dx.doi.org/10.1016/j.ejor.2014.07.046