PBI function based evolutionary algorithm with precise penalty parameter for unconstrained many-objective optimization

Abstract

Fixed or experiential penalty parameter of the penalty-based boundary intersection (PBI) function method cannot simultaneously ensure the convergence and diversity for all shape of Pareto front (PF). Too large penalty parameter may lead to bad convergence while too small parameter can not ensure the diversity. Specially, if the penalty parameter is too small, some reference weight vectors may have no solution on it. This error is hard to be rectified. In this paper, we prove that the lower bound of the penalty parameter is determined by three factors. The first one is the shape of the PF. The second one is the cosine distance between two adjacent reference vectors. The third one is the number of objectives. We deduce the lower bound of the penalty parameter. Once the penalty parameter was calculated, an individual with minimal PBI function is attached to the corresponding reference vector. The minimal-PBI-function-first principle is used in the environmental selection to guarantee the wideness and uniformity of the solution set. The time complexity is low. The proposed method is compared with other three state-of-the-art many-objective evolutionary algorithms on the unconstrained test problems MaOP, DTLZ and WFG with up to fifteen objectives. The experimental results show the competitiveness and effectiveness of the proposed algorithm in both time efficiency and accuracy.

Introduction

Multi-Objective Problems (MOPs) with more than three objectives are informally known as Many Objective Optimization Problems (MaOPs) [1,2]. Generally, a MaOP with only box constraints can be stated as follows.

$$\text{minimize}\mathit{F}\left(\mathit{x}\right)={\left(f_1\left(\mathit{x}\right),\cdots ,f_m\left(\mathit{x}\right)\right)}^T\text{subject to x}\in \Omega$$

$\text{Ω}\subseteq {\mathbb{R}}^n$ is the decision (variable) space, $\mathbf{x}={\left(x_1,\cdots ,x_n\right)}^T\in \text{Ω}$ is a candidate solution. Here n is the number of decision variables, and m is the number of objectives. $F:\text{Ω}\to {\mathbb{R}}^m$ constitutes $m$ objective functions, and ${\mathbb{R}}^m$ is the objective space. The attainable objective set is defined as $\text{Θ}=\left\{F\left(\mathbf{x}\right)|\mathbf{x}\in \text{Ω}\right\}$. The ${\mathit{x}}^{\left(1\right)}$ is said to dominate ${\mathit{x}}^{\left(2\right)}$ (denoted as ${\mathit{x}}^{\left(1\right)}\preccurlyeq {\mathit{x}}^{\left(2\right)}$) iff for $\forall$ $i\in \left\{1,\cdots ,m\right\}{f}_i\left({\mathit{x}}^{\left(1\right)}\right)\le f_i\left({\mathit{x}}^{\left(2\right)}\right)$ and $\exists j\in \{1,\cdots ,m\}$ $f_j\left({\mathit{x}}^{\left(1\right)}\right)<f_j\left({\mathit{x}}^{\left(2\right)}\right)$. A solution ${\mathit{x}}^{\text{*}}$ is a Pareto-optimal solution to formula (1) if there is no other solution $\text{x}\in \text{Ω}$ such that $\mathit{x}\preccurlyeq {\mathit{x}}^{\ast }$. $\mathit{F}\left({\mathit{x}}^{\text{*}}\right)$ is then referred to as Pareto-optimal (objective) vector. The set of all Pareto-optimal solutions is referred to as the Pareto-optimal set (PS). Accordingly, the set of all Pareto-optimal vectors, $\text{ }PF=\left\{F\left(\mathit{x}\right)\in R^m|\mathit{x}\in \text{PS}\right\}$ is referred to as the Pareto front (PF) [11].

Recently, MaOPs are extensively researched within the Evolutionary Multi-objective Optimization (EMO) community. A number of Many-Objective Evolutionary Algorithms (MOEAs) have been proposed to deal with MaOPs [[4], [5], [6], [7], [8], [9], [10], [11]]. Based on the Pareto optimization and dominance, lots of MOEAs are proposed. The improved Strength Pareto Evolutionary Algorithm (SPEA2) [12] is proposed by Zitzler et al. The Nondominated Sorting Genetic Algorithm II (NSGA-II) is proposed in Ref. [13]. Besides, the Pareto Archived Evolution Strategy (PAES) [14] and Pareto envelope-based selection algorithm II (PESAII) [15] belong to such category. Due to the obstacles caused by the Dominance Resistance (DR) and the Active Diversity Promotion (ADP) phenomena [4], the Pareto based methods encounter difficulties when handling MaOPs. The DR phenomenon refers to the incomparability of solutions in terms of the Pareto dominance relation. The main reason is that the proportion of non-dominated solutions in a population rise rapidly as the number of objectives increases [16]. The probability of any two solutions being comparable is η = 1/2^m−1 in an m-dimensional objective space [17]. When m = 10, η is as low as 0.002. As a result of DR phenomenon, the dominance-based primary selection criterion fails to distinguish between solutions, and then the diversity based secondary selection criterion is activated to determine the survival of solutions in the environmental selection. This is the so-called ADP problem which may be harmful to the convergence of the approximated Pareto fronts. Furthermore, for the 10-objective problem, the diversity maintenance mechanism even makes the population gradually move away from the true Pareto front. To address the DR phenomenon, a number of modified or relaxed Pareto dominance relations were proposed to enhance the selection pressure toward the true Pareto front. The typical methods include the $\mathrm{\epsilon }$-dominance [18], α-dominance [19], fuzzy Pareto dominance [20], L-dominance [21], θ-dominance [6], dominance area control [22], grid-dominance [23] and preference order ranking [24].

Another category of MOEA is the decomposition-based method, which decomposes a MOP into a set of subproblems and optimizes them in a collaborative way. The MOEA/D proposed in Ref. [25] is a representative of this category. In this algorithm, a predefined set of weight vectors are used to specify search directions towards different parts of the PF. In MOEA/D, each solution is associated with a subproblem, and each subproblem is optimized by using information from its neighborhoods. Since the weight vectors of subproblems are widely distributed, the obtained solutions are expected to have a wide spread over the PF. The original study on MOEA/D investigated three decomposition methods namely—the WS, the weighted TCH and the penalty-based boundary intersection (PBI).

The penalty-based boundary intersection(PBI) approach used in literatures [[25], [26], [27], [28], [29]] is defined as follows.

$$\text{Minimize}\text{g}\left(\mathit{x}|\mathit{w},{\mathit{z}}^{\ast }\right)=d_1+\theta d_2\text{Subject to x}\in \Omega$$

where

$$d_1=\frac{{\left(\mathit{F}\left(\mathit{x}\right)-{\mathit{z}}^{\ast }\right)}^T\mathit{w}}{\mathit{w}}$$

$$d_2=\mathit{F}\left(\mathit{x}\right)-\left({\mathit{z}}^{\ast }+d_1\frac{\mathit{w}}{\mathit{w}}\right)$$

where ${\mathit{z}}^{\ast }={\left(z_1^{\ast },z_2^{\ast }\cdots z_m^{\ast }\right)}^T$ is ideal objective vector with $z_i^{\ast }=\underset{\mathit{x}\in \text{Ω}}{\min }{f}_i\left(\mathbf{\Omega }\right).$ Fig. 1 presents an example to illustrate d₁ and d₂ of a solution x with regard to a weight vector w=(0.5, 0.5)^T. The d₁ is used to evaluate the convergence of x toward the PF, and d₂ is a measure of population diversity. $\text{g}\left(\mathit{x}|\mathit{w},{\mathit{z}}^{\ast }\right)=d_1+\theta d_2$ plays as a composite measure of both convergence and diversity. The goal of PBI approach is to push F(x) as low as possible so that it can reach the boundary of the $\text{PF}$ [29]. The major shortcomings of this approach are their applicability to MaOPs [26] and that it may result in low population diversity in some problems [27,28].

Sato [26] extends the conventional PBI approach and proposes inverted PBI (IPBI) decomposition method. In the conventional decomposition methods, such as the TCH and the PBI, solutions are evolved toward the reference point by minimizing the scalarizing function value. On the contrary, in the IPBI approach, solutions are evolved from the nadir point by maximizing the scalarizing function value [28]. The experimental study on MOKPs and WFG4 problem [30] with 2–8 objectives illustrates that the IPBI approach can better approximate widely spread PF in comparison to other scalarizing approaches. Just like the PBI, the limitation of the IPBI approach is that it involves the parameter θ which needs to be appropriately tuned. Yang et al. [31] suggest two penalty schemes, referred to as adaptive penalty scheme (APS) and subproblem-based penalty scheme (SPS), to set the value of the parameter θ. In the APS, θ is linearly increased with the number of generations from θ_min (=1) to θ_max (=10). A small value of θ is adopted initially so as to emphasize convergence and drive the search toward PF as fast as possible. The value of θ is gradually increased so as to emphasize diversity toward the later search stage. One of the main drawbacks of the PBI decomposition approach is that there is no unique setting of the penalty parameter θ that works well on different types of problems with different number of objectives [32].

As pointed out by Anupam Trivedi et al. [28], introduction of new methods to adaptively control parameter θ can be interesting research direction. This is what we do in this paper. The basic idea of the proposed PBI function based method comes from (1) Now that d₁ evaluates the convergence of x toward the PF, the time consuming Pareto domination calculation is not necessary; (2) Now that d₂ governs the diversity, the niche preservation operator and neighborhood concept for mating restriction are not necessary. Only if the value of $\theta$ is properly set, the convergence and diversity can be ensured. We provide a mathematic foundation of θ. Without the dominance calculation, the proposed method is more efficient and accurate. The major contributions of this paper are as follows.

• Precise penalty parameter calculation

The value of θ depends on the shape of PF and the number of objectives. We provide the mathematic foundation of θ calculation. This method applicable to various types of problems with different number of objectives. The calculation is simple and easy to use.

• Stable normalization

Traditionally, the nadir point ($z_1^{max},\cdots ,z_m^{max}$) is derived within the combination (of parent and child population) set. The $z_i^{max}$ denotes the extreme value of the ith objective. In the child population, there might be ‘bad’ individual which is far away from PF and has extraordinary large value of objective. It makes the search process unstable [32]. We know that the ‘bad’ individuals are eliminated after the environmental selection. Hence, the nadir and points are calculated after the environmental selection. These points are used in the next iteration of the MOEA. Besides, a damping mechanism is introduced into the normalization calculation. Initially, the normalization is less sensitive to the ideal and nadir points so as to emphasize convergence and drive the search toward PF as fast as possible. With the increase of evolutionary generations, the normalization is more and more sensitive to the ideal and nadir points so as to emphasize diversity.

This paper is organized as follows. The mathematic foundation of θ calculation is illustrated in Section 2. Our proposed algorithm is descripted in Section 3. Experimental settings and comprehensive experiments are discussed in Section 4. Finally, Section 5 concludes this paper.