The Artificial Neural Networks Based on Scalarization Method for a Class of Bilevel Biobjective Programming Problem.

A two-stage artificial neural network (ANN) based on scalarization method is proposed for bilevel biobjective programming problem (BLBOP). The induced set of the BLBOP is firstly expressed as the set of minimal solutions of a biobjective optimization problem by using scalar approach, and then the whole efficient set of the BLBOP is derived by the proposed two-stage ANN for exploring the induced set. In order to illustrate the proposed method, seven numerical examples are tested and compared with results in the classical literature. Finally, a practical problem is solved by the proposed algorithm.

1. Introduction

The bilevel programming problem (BLP) is a nested optimizations problem with two levels in a hierarchy: the upper and lower level decision-makers. The upper level maker makes his decision firstly, followed by the lower level decision-maker. The objective function and constraint of the upper level problem not only rely on their own decision variables but also depend on the optimal solution of the lower level problem. The decision-maker at the lower level has to optimize his own objective function under the given parameters from the upper level decision-maker. Since many practical problems, such as engineering design, management, economic policy, and traffic problems, can be formulated as hierarchical problems, BLP has been studied and received increasing attention in the literatures. During the past decades, some surveys and bibliographic reviews were given by several authors [1-4]. Reference books on bilevel programming and related issues have emerged [5-8].

The bilevel programming problem is a nonconvex problem, which is extremely difficult to solve. As we know, BLP is a NP-Hard problem [9-11]. Vicente et al. [12] also showed that even the search for the local optima to the bilevel linear programming is NP-Hard. Even so, many researchers are devoted to develop the algorithms for solving BLP and propose many efficient algorithms. To date a few algorithms exist to solve BLP; it can be classified into four types: Karus-Kuhn-Tucker approach (KKT) [13-16], Branch-and-bound method [17], penalty function approach [18-21], and descent approach [22, 23].

Unfortunately, the bilevel programming problem is nonconvex and the properties such as differentiation and continuity are necessary when proposing the traditional algorithms. Thus, many researchers tend to propose the heuristic algorithms for solving BLP because of their key characteristics of minimal problem restrictions such as differentiation. Mathieu et al. [24] firstly developed a genetic algorithm (GA) for bilevel linear programming problem because of its good characteristics such as simplicity, minimal problem restrictions, global perspective, and implicit parallelism. Motivated by the same reason, other kinds of genetic algorithm for solving bilevel programming were also proposed in [25-28]. Because of the prominent advantage that neural computing can converge to the equilibrium point (optimal solution) rapidly, the neural network approach was used to solve bilevel programming problem in [29-31]. In additional, since McCulloch and Pitts [32] and Pyne [33] utilized logical calculus to emulate nervous activities, there have been various types of analogue neural networks proposed for computation. Sheng et al. [34] firstly proposed a neural network approach based on Frank-Wolfe method for a class of BLP problems. Shih et al. [35] and Lan et al. [29] presented neural network for solving linear BLP problem, respectively. Recently, Lv et al. [30, 36] investigated the nonlinear bilevel programming and the convex quadratic bilevel programming by neural network. However, there are few results on application of neural networks to the BLBMPP in the literature.

Particle swarm optimization (PSO) is a relatively novel heuristic algorithm inspired by the choreography of a bird flock. Due to its high speed of convergence and relative simplicity, the PSO algorithm has been employed for solving BLP problems. For example, Li et al. [37] proposed a hierarchical PSO for solving BLP problem. Kuo and Huang [38] applied the PSO algorithm for solving bilevel linear programming problem. Jiang et al. [39] presented the PSO based on CHKS smoothing function for solving nonlinear bilevel programming problem. Gao et al. [40] presented a method to solve bilevel pricing problems in supply chains using PSO. Zhang et al. [41] presented a new strategic bidding optimization technique which applies bilevel programming and swarm intelligence. The hybrid algorithms based on PSO are also proposed to solve the bilevel programming problems [42-44]. Besides, Tabu search [45-47], simulated annealing [48], ant colony optimization [49], and λ-cut and goal-programming-based algorithm [50] are also typical intelligent algorithms for solving bilevel programming problem.

However, the algorithms mentioned above are only for the simple single objective bilevel programming problems. In fact, the multiobjective characteristics widely existing in the BLPP and the bilevel multiobjective programming problem (BLMPP) have attracted many researchers' interesting. For example, Shi and Xia [51, 52], Abo-Sinna and Baky [53], and Nishizaki and Sakawa [54] presented an interactive algorithm for BLMPP. Eichfelder [55] developed a numerical method for solving nonlinear nonconvex bilevel multiobjective optimization problems. In recent years, the metaheuristic has attracted considerable attention as an alternative method for BLMPP. For example, Deb and Sinha [56-58] as well as Sinha and Deb [59] discussed BLMPP based on evolutionary multiobjective optimization principles. Based on those studies, Deb and Sinha [60] proposed a viable and hybrid evolutionary-local-search based algorithm and presented challenging test problems. Sinha [61] presented a progressively interactive evolutionary multiobjective optimization method for BLMPP. Lately, Zhang et al. [62] proposed an improved PSO for BLMPP and a framework of PSO for solving BLMPP is established. Subsequently, Zhang et al. [63] proposed an elite quantum behaved PSO for relatively complex BLMPP. In 2013, Zhang et al. [64] proposed a hybrid particle swarm optimization algorithm with crossover operator to solve high-dimensional BLMPP. Almost all the research object of the BLMPPs is the BLBOP, so we mainly consider the BLBOPs in this paper.

As we known, the authenticity of the lower level Pareto optimal solution is very important for the BLBOP. If the obtained optimal Pareto solutions possess the fraudulence, it can lead to the failure to solve the whole problem. In this paper, the induced set of the BLBOP is firstly expressed as the set of minimal solutions of a biobjective optimization problem by a scalar approach which can greatly improve the accuracy of the lower level Pareto optimal solutions. Based on the efficient set of the BLBOP, a two-stage ANN is presented for solving whole problem which can reduce the computation burden.

The remaining of this paper is organized as follows. In Section 2, we give the formulation of the model and related definitions. In Section 3, we will introduce the scalar approach for the induced set of the BLBOP and the two-stage ANN algorithm for the whole problem. In Section 4, some numerical examples and practical problem are given to demonstrate the feasibility and efficiency of the proposed algorithm, while the conclusion is reached in Section 5.

2. Problem Formulation and Main Theory Results

Let X be a nonempty subset of R, Y be a nonempty subset of R, and F : R × R → R, f : R × R → R, G : R × R → R, and g : R × R → R be vector-valued mappings. We consider the following bilevel multiobjective programming problem (BLMPP): where F(x, y) and f(x, y) are the upper level and the lower level objective functions, respectively. x ∈ I⊆X and g(x, y) denote the upper level and the lower level constraints, respectively. Let I = [c, d] ⊂ X and G = {(x, y)∣g(x, y) ≥ 0, x ∈ I}. For the fixed x ∈ I, let denote the Pareto optimal solutions to the lower level problem; the induced set of problem (1) is denoted as . Note that the constraint x ∈ I is uncoupled from the lower level variable y ∈ Y. Particularly, if m1 = m2 = 2, we also call the BLMPP as bilevel biobjective programming problem (BLBOP). In the following, we shall focus on the BLBOP and assume that m1 = m2 = 2 and n1 = 1.

Definition 1 .

Let K be a closed pointed convex cone. A point is called a Pareto optimal solution of the lower level problem with respect to K if .

Definition 2 .

For a fixed x ∈ I, if y is a Pareto optimal solution to the lower level problem, then (x, y) is a feasible solution to problem (1).

Definition 3 .

If (x, y) is a feasible solution to problem (1) and there are no (x, y) ∈ Ω, such that F(x, y)≺F(x, y), then (x, y) is a Pareto optimal solution to problem (1), where “≺” denotes Pareto preference.

Definition 4 .

If (x, y) is the optimistic solution for problem (1), then (x, y) is given by .

Remark 5 .

The optimistic solution to the BLBOP is the one that optimizes the leader's objective function over the set of efficient solutions to the follower, assuming that the follower has no preferences among the efficient solutions obtained for each leader's decision x or that the follower will choose the solution that most benefits the leader. In this paper, we only consider the optimistic BLBOP.

For problem (1), it is noted that a solution (x, y) is feasible for the upper level problem if and only if y is an optimal solution for the lower level problem with x = x. In practice, we often make the approximate optimal solutions of the lower level problem as the optimal response feedback to the upper level problem, and this point of view is accepted usually. On the other hand, the authenticity of the lower level Pareto optimal solution is very important for the BLMP. If the obtained optimal Pareto solutions possess the fraudulence, it can lead to the failure to solve the whole problem. In this paper, we proposed the scalar approach for the lower level Pareto optimal solutions in order to improve the accuracy of the lower level Pareto optimal solutions.

3. Algorithm

3.1. The Scalar Approach for the Induced Set of BLBOP

In bilevel optimization, the constraint set of the upper level problem is given by the solution set of the lower level optimization problem. According to the Theorem 4.1 of the literature [55], the induced set Ω of problem (1) is equivalent to the Pareto optimal set of the multiobjective optimization problem

Thus, to solve the induced set of BLBOP is transformed to solve the Pareto optimal solution set of problem (2). Inspired by the scalar approach adopted in [55], the approximation for the induced set of problem (1) can be solved by the following algorithm.

Algorithm 1 .

Step 1. Discretize the upper level variable.

Step 1.1. For one-dimension upper level decision variable, choose θ > 0 and discretize I = [c, d] by x1 = c, x = x1 + (k − 1)θ  (k = 2,3,…, n − 1), x = d.

Step 1.2. For n-dimension (n ≥ 2) upper level decision variable, choose θ > 0 and discretize x ∈ I = [c1, d1] × [c2, d2] × ⋯×[c, d] by {x ∈ R∣x = c + nθ, x ≤ d}.

Step 2. Initialize the loop variables k = 1 and predefine an accuracy measure α > 0.

Step 3. Execute the following steps for x = x  (k = 1,2,…, n) and let .

Step 3.1. Solve and . Then, determine and . Denote the approximate Pareto optimal solution set by . Let with a small δ > 0 and t = 2.

Step 3.2. If a1 ≤ a1, solve problem (2) using the scalarization method in [56] with , r = (0,1) for the minimal solution y. Set A0, = A0, ∪ {(x, y)}.

Step 3.3. Calculate the M, N according to Theorem 4.3 and determine according to (5.2) in [56]. Set and let t = t + 1; go to Step 4.

Step 4. Output the approximation Pareto optimal solution set A0; that is, A0 = ⋃A0,.

3.2. The Two-Stage ANN

If θ > 0 in Step 2 is small, then set A0 will consist of many points. Then, the determination of all nondominated points for the upper level problem can be very expensive. In this subsection, the two-stage ANN is presented for determining the Pareto optimal solution of the BLBOP on the induced set A0. The first stage is to map the vectors from R to R and the second stage is to determine the Pareto optimal solutions of problem (1). The detail of each stage is described as follows.

The first stage of the ANN is a feed-forward artificial neural network (FFANN) which is composed by two subnetworks A and B with the same structure and neuron output function. For each subnetwork, the nodes are organized into two layers and the weighted arcs only link nodes in lower layers to nodes in higher layers. The first stage is to compute the objective function value for the upper lever.

The second stage is a quasineural artificial network; namely, the network has no connectivity weight and the output value can be computed directly by software. The input layer of the quasineural artificial network is the output layer of the first stage network. For the hidden layer, the input of and the output are defined by (3) and (4), respectively.

For the output layer, the input of and the output are defined by (5) and (6), respectively.

Based on the set A0 obtained by Algorithm 1, the first approximation Pareto optimal solutions of problem (1) can be achieved by the following algorithm.

Algorithm 2 .

Step 1. Discrete induced set A0, denote the discrete set by A0, and let s = 0.

Step 2. Divide the set A0 into two subsets A0 and A0 randomly; that is A0 = A0 ∪ A0.

Step 3. Input x ∈ A0 and x ∈ A0 into subnetwork A and subnetwork B, respectively.

Step 4. Select Pareto optimal solutions.

Step 4.1. If O2(S) = 1, then A0 = A0∖x. Let i = i + 1, go to Step 4.

Step 4.2. If O2(S) = −1, then A0 = A0∖x. Let j = j + 1, go to Step 4.

Step 4.3. If O2(S) = 0, then x, x ∈ A0. Go to Step 4.

Step 5. If A0 = ϕ, A0 ≠ ϕ, then (∀x ∈ A0) ⊂ A0.

Step 6. If A0 ≠ ϕ, A0 = ϕ, then (∀x ∈ A0) ⊂ A0.

Step 7. If the stopping criterion is met, then stop. Otherwise, go to Step 2.

In Step 2, the two subsets A0 and A0 have approximately equal number of feasible solutions. The output of Algorithm 2 is the first approximation Pareto optimal solution set of problem (1); we denote the set by P(A0).

3.3. The Algorithm for BLBOP

Based on Algorithms 1 and 2, as well as the refinement strategy used by the literature [9] which is employed in this paper, the algorithm for solving the BLBOP can be described as follows.

Algorithm 3 .

Step 1. Based on the induced set A0 obtained by Algorithm 1, we determine the approximation Pareto optimal solution set P(A0) of problem (1) by Algorithm 2.

Step 2. Let i = 0 and choose the distance γ > 0.

Step 3. For any point (x, y) ∈ P(A0), determine the refinement induced set around this point by the refinement strategy according to (5.4) in [55].

Step 4. Set ; update the approximation Pareto optimal solution set A by Algorithm 2.

Step 5. If approximation of solution set of problem (1) is sufficient, then stop. Otherwise set i = i + 1, choose γ > 0, and go to Step 3.

4. Results

In this section, we considered seven numerical examples and a practical problem to illustrate the feasibility of the proposed algorithm for problem (1). In order to evaluate the closeness between the obtained Pareto optimal front and the theoretical Pareto optimal front, as well as the diversity of the obtained Pareto optimal solutions along the theoretical Pareto optimal front, we adopted the following evaluation metrics.

4.1. Performance Evaluation Metrics

4.1.1. Generational Distance (GD)

This metric used by Deb [65] is employed in this paper as a way of evaluating the closeness between the obtained Pareto optimal front and the theoretical Pareto optimal front. The GD metric denotes the average distance between the obtained Pareto optimal front and the theoretical Pareto optimal front:where n is the number of the obtained Pareto optimal solutions by the proposed algorithm and d is the Euclidean distance between each obtained Pareto optimal solution and the nearest member of the theoretical Pareto optimal set.

4.1.2. Spacing (SP)

This metric is used to evaluate the diversity of the obtained Pareto optimal solutions by comparing the uniform distribution and the deviation of solutions as described by Deb [65]:where d = min(|F1(x, y) − F1(x, y)| + |F2(x, y) − F2(x, y)|), i, j = 1,2,…, n, is the mean of all d, d is the Euclidean distance between the extreme solutions in obtained Pareto optimal solution set and the theoretical Pareto optimal solution set on the mth objective, M is the number of the upper level objective function, and n is the number of the obtained solutions by the proposed algorithm.

All results presented in this paper have been obtained on a personal computer (CPU: AMD 2.80 GHz; RAM: 3.25 GB) using a c# implementation of the proposed algorithm.

4.2. Numerical Experiment

In this section, we will present seven BLBOPS to illustrate the proposed algorithm for the bilevel biobjective programming. Problem 1 and problem 2 are low-dimensional problems. Problems 3–6 are high-dimensional problems. For problem 7, the theoretical optimal front is unknown. In this paper, we refined every upper problem's feasible solutions three times and the obtained results are compared with the classical literature.

Example 1 .

Example 1 is taken from [57]. Here x ∈ R1, y ∈ R2. In this example, the parameters γ0, γ1, γ2 of Algorithm 3 are 0.23, 0.11, and 0.03.

Figure 1(a) shows the obtained Pareto front by the proposed algorithm and the method in [57]. Table 1 shows the comparison results between the two algorithms considering the metrics previously described. It can be seen that the performance of the proposed algorithm is better with respect to the generational distance, although it places slightly below the method in [57] with respect to spacing. By looking at the obtained Pareto fronts, some nondominated vectors produced by the method in [57] are not part of the true Pareto front of the problem; however, the proposed algorithm is able to cover the full Pareto front. Figure 2(b) shows that all the obtained solutions by the proposed algorithm, which follow the relationship; that is, y1 = −1 − y2, , and . However, some solutions obtained by the method in [57] do not meet the relationship.

Figure 1

The obtained Pareto front and solutions of Example 1.

Table 1

Results of the Generation Distance (GD) and Spacing (SP) metrics for Examples 1 and 2.

Prob.	GD		SP
	The method in [57]	The proposed method	The method in [57]	The proposed method
1	0.00216	0.00097	0.01135	0.01201
2	0.01013	0.00312	0.00203	0.00969

Figure 2

The obtained Pareto front and solutions of Example 2.

Example 2 .

Example 2 is taken from [57]. Here x ∈ R1, y ∈ R2. In this example, the parameters γ0, γ1, γ2 of Algorithm 3 are 0.3, 0.21, and 0.12.

Figure 2(a) shows the obtained Pareto front by the proposed algorithm and the method in [57]. Obviously, both the proposed algorithm and the method in [57] almost have the same spacing. However, there are some areas of the Pareto optimal solution obtained by the method in [57] which is sparse. In Figures 2(a) and 2(b), it can be seen that the obtained solutions by the proposed algorithm almost follow the relationship; that is, x = y1,  y1 ∈ [0.5,1], y2 = 0. However, some areas of the solutions obtained by the method in [57] are sparse and some solutions do not meet the relationship.

Example 3 .

Example 3 is taken from [60]. Here x ∈ R10, y ∈ R10. In this example, the parameters γ0, γ1, γ2 of Algorithm 3 are 0.5, 0.36, and 0.12.

This problem is more difficult compared to the previous problems (Examples 1 and 2) because the lower level problem of this example has multimodalities, thereby making the lower level problem difficult in finding the upper level Pareto optimal front. Figure 3 shows the graphical results produced by the method in [60] and the proposed algorithm in this paper. Tables 2 and 3 show the comparison of results between the two algorithms considering the metrics previously described. It can be seen that the performance of the proposed algorithm is the best with respect to the generational distance. By looking at the Pareto fronts of this test problem, some nondominated vectors produced by the method in [60] are not part of the true Pareto front and there are some areas are sparse. However, the proposed algorithm is able to cover the full Pareto optimal front. Furthermore, two obtained lower level Pareto optimal fronts are given when x1 = 2 and x1 = 2.5.

Figure 3

The obtained Pareto front of Example 3.

Table 2

Results of the Generation Distance (GD) metrics for Examples 3, 4, 5, and 6.

Prob.	The proposed algorithm	The method in [60]
3	0.00019	0.00768
4	0.00015	0.06391
5	0.00038	0.01677
6	0.00021	0.00652

Table 3

Results of the Spacing (SP) metrics for Examples 3, 4, 5, and 6.

Prob.	The proposed algorithm	The method in [60]
3	0.00076	0.00197
4	0.00273	0.00269
5	0.00299	0.01737
6	0.00130	0.00127

Example 4 .

Example 4 is taken from [60]. Here x ∈ R10, y ∈ R10. In this example, the parameters γ0, γ1, γ2 of Algorithm 3 are 0.6, 0.38, and 0.16.where

For Example 4, the upper level problem has multimodalities, thereby causing an algorithm difficulty in finding the upper level Pareto optimal front. Figure 4 shows the obtained Pareto front by the method in [60] and the proposed algorithm. Tables 2 and 3 show the comparison of results among the three algorithms considering the metrics previously described. It can be seen that the performance of the proposed algorithm is better than the method in [60] with respect to the generational distance though they almost have the same performance with respect to the spacing.

Figure 4

The obtained Pareto front of Example 4.

Example 5 .

Example 5 is taken from [60]. Here x ∈ R10, y ∈ R10. In this example, the parameters γ0, γ1, γ2 of Algorithm 3 are 0.7, 0.49, and 0.23.

Figure 5 shows the obtained Pareto front of Example 5 by the method in [60] and the proposed algorithm. Tables 2 and 3 show the comparison of results between the two algorithms considering the metrics previously described. For this example, the graphical results again indicate that the method in [60] does not cover the full Pareto front. Since the nondominated vectors found by the method in [60] are clustered together, the spacing metric provides very good results. Graphically, it can be seen that the proposed algorithm is able to cover the entire Pareto front. It is also interesting to note that the Pareto optimal fronts for the upper level lie on constraint boundaries and every lower level Pareto optimal front has an unequal contribution to the upper level Pareto optimal front.

Figure 5

The obtained Pareto front of Example 5.

Example 6 .

Example 6 is taken from [60]. Here x ∈ R1, y ∈ R9. In this example, the parameters γ0, γ1, γ2 of Algorithm 3 are 0.5, 0.23, and 0.17.

Figure 6 shows the obtained Pareto front of Example 6 by the method in [60] and the proposed algorithm. Tables 2 and 3 show the comparison of results between the two algorithms considering the metrics previously described. It can be seen that our algorithm and the method in [60] almost have the same performance of the spacing, but some nondominated vectors produced by the method in [60] are slightly off the true Pareto front. Moreover, three obtained lower level Pareto optimal fronts are given when y1 = 1, y1 = 1.5, and y1 = 2. It can be seen that only one Pareto optimal point from each participating lower level problem qualifies to be on the upper level Pareto optimal front.

Figure 6

The obtained Pareto front of Example 6.

Example 7 .

Example 7 is taken from [55]. Here x ∈ R1, y ∈ R2. In this example, the parameters γ0, γ1, γ2 of Algorithm 3 are 0.3, 0.16, and 0.07, respectively.

Figure 7 shows the final archive solutions by the proposed algorithm. For this problem, the exact Pareto optimal front is not known, but the obtained Pareto optimal front after four approximations by the proposed algorithm is similar to that reported in the previous study [55].

Figure 7

The obtained Pareto optimal solutions with A3.

4.3. Application of the Algorithm for a Practical Problem

In a company, the CEO's goal is usually to maximize net profits and quality of products, whereas a branch head's goal is to maximize its own profit and worker satisfaction. The problem involves uncertainty and is bilevel in nature, as a CEO's decision must take into account optimal decisions of branch heads. We present a deterministic version of the case study from [66] in (4). Figure 8 shows the obtained Pareto optimal front of this practical problem by the proposed algorithm. Note that Zhang et al. [66] only obtained a single optimal solution x = (146.2955,28.9394) and y = (0,67.9318,0) which lies on the maximum of F2 using weighted sum method. In contrast, a set of Pareto optimal solutions is obtained by the proposed algorithm. However, the fact that the single optimal solution in [66] is included in the obtained Pareto optimal solutions illustrates the feasibility of proposed algorithm. In this problem, the parameters γ0, γ1, γ2 of Algorithm 3 are 0.2, 0.11, and 0.02.

Figure 8

The obtained Pareto front of the practical problem.

5. Conclusions

In this paper, a two-stage ANN based on scalarization method is presented for solving BLBOPs. Seven numerical examples and a practical problem are used to state the feasibility and efficiency of the proposed algorithm. The experimental results indicate that the obtained Pareto front by the proposed algorithm is very close to the theoretical Pareto optimal front, and the solutions are also distributed uniformly on entire range of the theoretical Pareto optimal front. The proposed algorithm is easy to implement, which provides another appealing method for further study on the general BLMPP.