NatureInspired Algorithms and Applications: Selected Papers from CIS2013
View this Special IssueResearch Article  Open Access
A Differential Evolution with Two Mutation Strategies and a Selection Based on an Improved ConstraintHandling Technique for Bilevel Programming Problems
Abstract
Two mutation operators are used in the differential evolution algorithm to improve the diversity of population. An improved constrainthandling technique based on a comparison mechanism is presented, and then it is combined with the selection operator in the differential evolution algorithm to fulfill constraint handling and selection simultaneously. A differential evolution with two mutation strategies and a selection based on this improved constrainthandling technique is developed to solve bilevel programming problems. The simulation results on some linear and nonlinear bilevel programming problems show the effectiveness and efficiency of the proposed algorithm.
1. Introduction
Bilevel programming problem (BLPP) deals with mathematical programming whose feasible set is implicitly determined by another embedded optimization problem [1]. As a hierarchical system, BLPP is widely used to lots of fields such as economy, control, engineering, and management (see [2â€“4]). In such a problem, there is a leader that makes the decision first and a follower that makes its decision based on the leaderâ€™s. More and more practical applications can be formulated as bilevel programming models; so it is important to design some effective algorithms to solve different types of BLPPs in order to satisfy the needs of realworld problems.
Various approaches have been developed in different branches of numerical optimization to solve BLPPs under various assumptions, including extreme point methods, branch and bound algorithms, cutting plane techniques, mixed integer programming, parametric complementary pivoting methods, penalty function methods, trust region methods, and heuristic algorithms. A detailed list of references can be found in [2â€“4].
As it is known, BLPPs being intrinsically difficult, it is not surprising that most algorithmic research to date has focused on the simplest cases of BLPPs, that is, problems having nice properties such as linear, quadratic or convex objective, and/or constraint functions. In particular, the most studied instance of BLPPs has been for a long time the linear version [5]. Linear BLPP was shown to be a strongly NPhard problem [6]. Furthermore, it is a NPhard problem to find the local optimal solution of BLPP [7]. Linear BLPP has the important property that at least one optimal solution is attained at an extreme point of the constraint region. Based on this property, the thbest approach was used successfully to compute global solutions of linear BLPPs by enumerating the extreme points of the constraint region [1]. A heuristic algorithm was proposed, in which classical extreme point enumeration techniques were combined with genetic algorithms by associating chromosomes with extreme points of the polyhedron [8]. Another customary way to solve linear BLPPs is to use the KarushKuhnTucker conditions to reduce the problem to a onelevel complementary slackness problem. To solve this nonlinear slackness problem, many distinct approaches were proposed, for example, the branch and bound algorithms [6, 9, 10], the hybrid neural network approach [11], the penalty function method [10], the branch and cut algorithm [12], genetic algorithms [13, 14], and so forth.
Due to their inherent nonconvexity and nondifferentiability, the general nonlinear BLPPs are more intractable than linear version. Moreover, there exist some types of BLPPs which have nondifferentiable upperlevel objective functions and can not be solved by traditional optimization algorithms based on gradient methods. In this case, heuristic algorithms such as evolutionary algorithm [20], differential evolution [21], and particle swarm optimization [22] are alternative methods. For solving nonlinear BLPPs, a variety of solution techniques based on heuristic algorithms can been found in the literature. For example, Liu [23] presented a genetic algorithm for StackelbergNash equilibrium. Oduguwa and Roy [24] proposed a genetic algorithm for the bilevel optimization problems. Wang et al. [25] proposed an evolutionary algorithm for nonlinear BLPPs based on a new constrainthanding scheme. Li et al. [26] proposed a hierarchical particle swarm optimization for solving general BLPPs. The BLPP is transformed to solve the upperlevel and lowerlevel problems iteratively by two variants of PSO. Wang et al. [27] developed a new evolutionary algorithm for solving a class of nonlinear BLPPs. Angelo et al. [17] proposed an algorithm which uses differential evolution to solve both the upper and lowerlevel optimization problems. Sinha et al. [28] presented a standard test suite of 12 problems and have solved the proposed test problems using a nested bilevel evolutionary algorithm.
In our previous work [29], a linear BLPP is transformed into a linear singlelevel implicit programming problem. A differential evolution (DE) algorithm with two mutation strategies and a selection based on an improved constrainthandling technique is developed to solve the linear implicit programming problem. Preliminary results were obtained for linear BLPPs, but this algorithm is considered to have a general applicability. In this paper, this algorithm is extended from linear version of BLPPs to nonlinear bilevel programming problems. Two types of bilevel problems including linear and nonlinear versions are used to test the effectiveness of the proposed algorithm. The results obtained are compared with other available algorithms in the literature.
This paper is organized as follows. The general BLPP is described in the next section. A DE algorithm with two mutation strategies and a selection based on an improved constrainthandling method is presented in Section 3. Computational results and comparison are provided in Section 4. Section 5 concludes this paper.
2. Statement of General BLPP
When two independent decision makers, the leader and follower, have conflicting objectives, this bilevel programming model arises. Without loss of generality, we consider the minimization problems in our discussion. The general BLPP under discussion is expressed as where and are box constraint sets of leaderâ€™s variable and followerâ€™s variable , respectively. , , , and .
The constraint region of problem (1) is denoted by : Denote the projection of onto the leaderâ€™s decision space by The followerâ€™s rational reaction set for is defined as Denote the inducible region (feasible set) by
The main feature of BLPP is that it involves two mathematical programming, the solution of one being part of the constraints of the other. That is to say, the followerâ€™s solution depends implicitly on the leaderâ€™s variable; in return, the leaderâ€™s objective optimal value is restricted by the followerâ€™s solution. So it is quite evident that problem (1) can be converted into the following equivalent singlelevel implicit programming problem:
It is also rewritten as follows:
We suppose that and may be continuous, yet nonconvex, even nondifferentiable, and and are differentiable and convex in for fixed. Under these assumptions, the constraint region of problem (1) is bounded and closed set.
If the lowerlevel problem is linear programming, we can use available linear programming solver to obtain for the given ; if the lowerlevel problem is quadratic programming, we can use available quadratic programming solver to obtain for the given ; for the generalized lowerlevel problem, we must turn to some available effective algorithms for solving this constrained optimization problem, such as trust region algorithm, interior algorithm, and activeset algorithm. In this paper, we adopt the optimization toolbox of MATLAB to solve lowerlevel problems. Certainly, the lowerlevel problem can be solved by evolutionary algorithms, but it is over timeconsuming.
3. The Proposed Differential Evolution Algorithm
Differential evolution (DE), introduced by [21], is one of the most prominent evolutionary algorithms (EAs) for solving realvalued optimization problems. DE has exclusive advantages, such as a simple and easytounderstand concept, ease of implementation, powerful search capability, and robustness. Using only a few control parameters, DE has gradually become more popular and has been used in many practical cases [33].
3.1. Individual Coding
For solving implicit programming problem (3), a vector is used to express an individual, which represents a decision variable of leader, where , .
If the bounds of each () of upperlevel variable are not given, the following method is used to establish the bounds of upperlevel variable ; that is, , :
3.2. Fitness Function
The fitness function is expressed as follows:
3.3. A ConstraintHandling Technique Based on a Comparison Mechanism
A constrainthandling technique is developed by improving Debâ€™s feasibilitybased comparison method [34], which is based on the following criteria.(i)Between two feasible solutions, the one with the lowest objective function value wins for minimization problem.(ii)If one solution is feasible and the other one is infeasible, the feasible solution wins.(iii)If both solutions are infeasible, the one with the lowest sum of constraint violations is preferred.
There are two changes in our constrainthandling method. Sum of constraint violations is replaced by maximum of constraint violations. The infeasible individuals in certain range (close to the boundary of the feasible region) with the good values of the objective function are selected for the next generation.
This improved constrainthandling technique is presented detailedly as follows.
Definition 1. An individual is called a feasible individual, if and only if .
Definition 2. For a given , obtain . Let for , for . is called constraint violation of an individual .
Definition 3. An individual is better than another individual , if one of the following criteria is satisfied.(1)If and are two feasible individuals, that is, and , as well as .(2)If both individuals and are infeasible, that is, and , as well as .(3)If one individual is feasible and another one is infeasible, the following cases are considered.(a)If and , as well as ;(b)if and , as well as ;(c)if and ,
where is a given small positive parameter related to the evolutionary generation . In this paper, for a given maximum evolutionary generation , the values of are listed as follows:
Definition 4. The best individual in the current generation is determined by the constraint violations and fitness values of all individuals in the current population according to the following rules. (1)If the constraint violations of all individuals are less than 0.0001, one individual with the lowest fitness value is the best individual in the current generation.(2)If the constraint violations of all individuals are greater than 0.0001, one individual with the lowest constraint violation is the best individual in the current generation.(3)If there is only part of all individuals whose constraint violations are less than 0.0001, one with the lowest fitness value in this part is the best individual in the current generation.
3.4. Hybrid Mutation Operation
Two mutation operators are adopted to generate two subpopulations, each of whom comprises individuals ( denotes population size and is even); then combine them together to form a population with individuals.
According to Definition 4, determine the best individual in the th generation denoted by .
The first subpopulation with individuals is generated by the following mutation operator: where , , , and are mutually distinct indexes in not equal to . denotes scaling factor.
The second subpopulation with individuals is generated by the following mutation operator: where , , , and are mutually distinct indexes in not equal to . is scaling factor.
If a mutant individual , , is not in box constraints, that is, , , then randomly select a real number in to substitute it.
3.5. Crossover Operation
For each pair of individuals constituted by the parent individual and mutation individual , the crossover individual is generated by the following strategy: where is a uniform random number, is the randomly selected index chosen once for each , and is a realvalued crossover rate constant in the range .
3.6. Selection Operation Based on the ConstraintHandling Technique
The abovementioned constrainthandling technique is introduced to the selection operator to decide the population for the next generation.
For , the trial vector and the parent vector are compared. The individual of the next generation is decided according to Definition 3:
3.7. The Proposed Algorithm
Now we present the differential evolution algorithm with two mutation strategies and a selection based on an improved constrainthandling technique as follows.
Step 1 (parameter setting). Set population size , scaling factor , crossover rate , and maximum number of generations .
Step 2 (initialization). Let . Generate randomly an initial population . For each individual , obtain lowerlevel optimal solution by available optimization solver. Evaluate the constraint violations and fitness values of all individuals and determine the best individual according to Definition 4.
Step 3 (mutation operation). For to , select randomly and . The first subpopulation is generated according to (6).
For to , select randomly and . The second subpopulation is generated according to (7).
Step 4 (crossover operation). Crossover individuals are produced by (8). For each crossover individual , obtain lowerlevel optimal solution by using available optimization solver. Evaluate the constraint violations and fitness values of all crossover individuals and determine the best individual according to Definition 4.
Step 5 (selection operation). According to selection rule, form the next population.
Step 6 (stopping criterion). If the stopping criterion is met, then stop and output the result. Otherwise, set and go to Step 3.
4. Computational Results and Comparison
To evaluate the performance of the proposed algorithm, we test two types of benchmark problems: one type is linear BLPP and another type is nonlinear BLPP, which includes linearquadratic programming, quadraticquadratic programming, linearlinear fractional programming, general nonlinear BLPP, and complex nonlinear BLPP with nondifferentiable leaderâ€™s objective function. The results of the proposed algorithm are compared with other algorithms in the literature.
During the simulations, we adopt the following parameter suite:(i)population size: ;(ii)scaling factor: , where is a random number in the range ;(iii)crossover rate ;(iv)maximum number of generations .
All experiments are performed on a notebook computer with 2.53â€‰GHz Dualcore Processor and 1.94â€‰GB of RAM in MATLAB software. Each test problem implemented 50 independent runs. For solving the lowerlevel optimization problems, we use corresponding function calls in the optimization toolbox of MATLAB to solve linear programming, quadratic programming, or other nonlinear programming.
Performance of the proposed algorithm is evaluated by the following criteria:(i)success rate (SR), that is, the percentage of convergence to the optimal solution in 50 independent runs, which can be used to evaluate the reliability of the proposed algorithm;(ii)mean number of fitness evaluations (MNFE) when the optimal solution is reached first in 50 independent runs, which is assessed the computational cost;(iii)standard deviation of fitness values (SD), which is used to measure the quality of the solution obtained;(iv)best values of upperlevel objective functions and lowerlevel objective functions , as well as best solution obtained .
In order to demonstrate the effectiveness of the improved constrainthandling technique, we compare it with Debâ€™s feasibilitybased method [34] in the same framework of algorithm and same parameters. Also, we compare the obtained results with those presented in the corresponding references.
4.1. Computational Results for Linear BLPPs
To examine the effectiveness of the proposed algorithm for linear BLPPs, we first test 9 linear bilevel programming problems and compare the results of the proposed algorithm with other algorithms in the literature. The linear test problems are provided in Appendix A.
The comparative results obtained by two constrainthandling methods are shown in Table 1. As shown in Table 1, the improved constrainthandling technique is slightly better than Debâ€™s method in terms of SR, MNFE, and SD. Furthermore, the proposed algorithm has a good performance according to stability, efficiency, and quality of solution obtained. Therefore it can be concluded that the improved constrainthandling technique is effective.

Table 2 presents the best results obtained by the proposed algorithm for linear BLPPs, including best solution, upperlevel and lowerlevel objective function values.

In order to demonstrate the improvement in algorithm performance by the use of hybrid mutation operator, we compare hybrid mutation operator with single mutation operator in the same framework of algorithm and same parameters. The comparative results are shown in Table 3. From this table, it can be seen that the reliability of algorithm is enhanced by using hybrid mutation operator in terms of SR. Although hybrid mutation operator requires a larger computational cost than single mutation operator (6) except for Problems A.4 and A.7, the stability of single mutation operator (6) is worse than that of hybrid mutation operator according to SR and SD. The computational cost of single mutation operator (7) is larger than that of hybrid mutation operator, but the quality of solution obtained is better than hybrid mutation operator for Problem A.7 in terms of SR and SD. Therefore the algorithmic performance is improved by the use of those two mutation operators.

The proposed algorithm is compared with the algorithms in the literature according to the best values of upperlevel and lowerlevel objective functions. The comparative results are shown in Table 4. From this table, the proposed algorithm exhibits the better results than other algorithms. For Problems A.1, A.8, and A.9, the best results obtained by the proposed algorithm are better than those of other algorithms in the literature. For remaining problems, they have same optimization results.

4.2. Computational Results for Nonlinear BLPPs
To validate the effectiveness of the proposed algorithm for nonlinear BLPPs, we test two sets of nonlinear bilevel programming problems. The first set has 10 examples selected from different sources of the literature, and the second set consists in solving the test collection (SMD1 to SMD12) proposed in [28]. Finally, we compare the results of the proposed algorithm with the other algorithms in the literature. The first set of nonlinear BLPPs is provided in Appendix B, and the second set of nonlinear BLPPs can be obtained in [28].
With same algorithmâ€™s framework and parameter setting, the improved constrainthandling method is compared with Debâ€™s feasibilitybased method [34] for first set of nonlinear BLPPs to illuminate the effectiveness of the improved constrainthandling method again. Comparative results of two constrainthandling methods are shown in Table 5. From this table, it can be seen that the improved constrainthandling method is better than the Debâ€™s method in terms of success rate, mean number of fitness evaluations, and standard deviation of fitness values. Moreover, this table shows that the proposed algorithm with the improved constrainthandling method is effective and can find the global optimal solutions with a 100 success rate, except for Problem B.4 with a 76 success rate.

In order to display the computational results obtained by the proposed algorithm for first set of nonlinear BLPPs, the best solution , the upperlevel objective value , and the lowerlevel objective value for all problems are provided in Table 6. This table shows that the best solutions obtained by the proposed algorithm are same as the global optimal solutions of all problems.

Table 7 illustrates the comparative results of the proposed algorithm with available algorithms in the literature for the first set of nonlinear BLPPs. As shown in Table 7, the proposed algorithm is able to find the global optimal solutions for the first set of nonlinear BLPPs. Compared against the other algorithms in the literature, the proposed algorithm can obtain equal optimal results for all problems.

In order to investigate further the performance of the proposed algorithm, the second set of nonlinear BLPPs is tested, which is a complex test collection with different dimensions in [28]. We performed 50 runs for each of the test problems with 5 and 10 dimensions and 11 runs for each of the test problems with 20 dimensions. In the case with five dimensions, , , and are used in SMD1 to SMD5 and SMD7 to SMD12, and , , , and are used in SMD6. In the case with ten dimensions, , , and are used in SMD1 to SMD5 and SMD7 to SMD12, and , , , and are used in SMD6. In the case with twenty dimensions, , , and are used in SMD1 to SMD5 and SMD7 to SMD12, and , , , and are used in SMD6, where , , , and are parameters of test problems in [28]. Abovementioned parameters of the proposed algorithm continue to be used here. When the accuracy of best upperlevel function value becomes less than , the algorithm terminates.
The statistical results for 5, 10, and 20dimensional test problems of second set are provided in Table 8, including success rate (SR), mean number of fitness evaluations (MNFE), and standard deviation of fitness values (SD). This table shows that the proposed algorithm can locate the optimal solution for 5 and 10dimensional test problems. Except for SMD10â€“ SMD12 with 20 dimensions, it can also obtain the optimal solution for other 20dimensional test problems. The proposed algorithm fails to find the optimal solution within a given accuracy (less than ) for SMD10SMD12 with 20 dimensions.

Table 9 presents the best results obtained by the proposed algorithm for 5 and 10dimensional test problems of second set, including the best solution , the best upperlevel objective value , and the best lowerlevel objective value .

The comparison of median accuracy achieved by the proposed algorithm and those reported in [28] is provided in Table 10. For 20dimensional test problems, the accuracy for the upper level (UL) and lower level (LL) is not given in [28]. From this table, it is obvious that the proposed algorithm is better than the nested bilevel evolutionary algorithm in [28] according to median value of the accuracy at the upper level and lower level for 5 and 10dimensional test problems, except for SMD5. Furthermore, the proposed algorithm can obtain the good accuracy for 20dimensional test problems except for SMD10â€“SMD12. The nested bilevel evolutionary algorithm in [28] cannot find feasible solutions for 10dimensional SMD9â€“SMD12, which are constrained problems, and 20dimensional SMD7â€“SMD12. This also shows that the improved constrainthandling method in this paper performs better than the constraint handling in [28], which belongs to Debâ€™s method [34].
