Abstract
When dealing with multi-objective optimization, the proportion of non-dominated solutions increase rapidly with the increase of optimization objective. Pareto-dominance-based algorithms suffer the low selection pressure towards the true Pareto front. Decomposition-based algorithms may fail to solve the problems with highly irregular Pareto front. Based on the analysis of the two selection mechanism, a dynamic reference-vector-based many-objective evolutionary algorithm(RMaEA) is proposed. Adaptive-adjusted reference vector is used to improve the distribution of the algorithm in global area, and the improved non-dominated relationship is used to improve the convergence in a certain local area. Compared with four state-of-art algorithms on DTLZ benchmark with 5-, 10- and 15-objective, the proposed algorithm obtains 13 minimum mean IGD values and 8 minimum standard deviations among 15 test problem.
Introduction
Unlike single-objective optimization problems, multi-objective optimization problems require to optimize several conflicting objective functions simultaneously. Because of the complex relationship between these objectives, it is often impossible to achieve an optimal solution for all objectives [4, 10]. In the case of minimization optimization problem, the reduction of an objective function value tends to increase the value of other objective functions, so the solution of the multi-objective optimization problem is usually a set [9, 15].
Definition of multi-objective optimization problem
Without loss of generality, taking the objective function minimization as an example, the unconstrained, static, continuous multi-objective optimization problem is defined as follows:
For single-objective optimization problem, there is only one optimal solution (including the existence of multiple solutions with the same function value). Different form single-objective optimization problem, optimal solution for multi-objective optimization usually is a set of solutions with a certain logical relationship. The related definition can be found at our former research [12].
Generally speaking, the solution individuals present a one-to-one positional relationship in the decision space and the objective space. This positional relationship is closely related to the objective number of optimization problem. The difference between multi-objective optimization problem and single-objective optimization problem mainly focuses on the following aspects: The optimal value of the single-objective optimization problem is only one, so the optimal solution is usually a single solution; while the optimal value of the multi-objective optimization problem is usually not unique, the optimal solution is usually a solution set [25]. From the comparison of algorithms, the single-objective optimization algorithm only needs to compare the unique numerical value of the optimal solution; the multi-objective optimization algorithm needs to compare the superiority and inferior relationship between the solution sets, that is, whether the solution set converges to the optimal frontier; whether the individual is evenly distributed on the leading edge. Comparison for single-objective optimization problem is easy as the superiority and inferior relationship of the objective function are obvious; while comparison for multi-objective optimization problem is complex for the reason that the individual can not be compared only by one objective value. The difference in the superiority and inferiority of different objective function may lead to the individual’s superiority in different objective functions, that is, the superiority and inferiority cannot be isolated. For the entire solution set, the information provided by different individuals is various with its position and surroundings. Therefore, it is necessary to judge the superiority and inferior relationship according to the two evaluation criteria of convergence or diversity [14]. Selection for evolutionary parents is not same for the problem with single-objective and multi-objective. When selecting the parents from the evolutionary population for single-objective evolutionary algorithm, the fitness of objective function is the first standard or maybe the only standard. When dealing with multi-objective optimization problem, the linkage of multi-objective should be considered for the reason that complicated relationship between solutions. Because the optimal result of multi-objective problem is no longer a solo solution, the distribution of the solution set lead the complex relationship for parent selection.
As the number of objective increases, the difficulty of multi-objective optimization problems increases. Although the introduction of Pareto dominance relationship has brought great achievement in the multi-objective optimization algorithm when dealing with the optimization problem with 2-3 optimization objectives, there are still many problems in many-objective optimization(more than 3 optimization objectives).
(1) Failure of the Pareto dominance relationship
For a single individual with a limited feasible domain in space, as the number of objective increases, the area dominated by this individual decreases exponentially, while the area that does not dominated by it increases exponentially. Therefore, the ability of Pareto dominance to judge individual performance also decreases exponentially [18].
Figure 1 shows the proportion of non-dominated solutions in a population with 150 individuals which are uniformly produced in [0, 1]. As the objective number grows, the proportion of non-dominated solutions increases rapidly. When the objective dimension reaches a certain number, the ranking method based on the Pareto dominance relationship is almost completely useless. Therefore, the Pareto sorting can not be used as the selection mechanism to make population converged [12].

Ratio of non-dominated individuals in the population with different objectives.
(2) Defect of distribution retention mechanism
Most multi-objective evolutionary algorithms use "convergence first, diversity second" as the core mechanism of algorithm design: convergence performance of one individual is the priority selection criterion. For the individuals with the convergence performance (non-dominated rank), distribution is used as the basis for individual selection. When the objective dimension is high, the convergence mechanism (based on Pareto dominance relation ordering) can’t provide a sorting basis for the algorithm, and the distribution will become the dominant of individual ordering. Then, the algorithm shows a better propensity for individuals with larger crowding distances. Because the distribution mechanism overemphasizes the extensiveness of individuals in the population, it may lead to many individuals with poor convergence being preserved. At this time, it is not conducive to the convergence of the population, no matter in the stage of individual selection or individual reproduction.
(3) Increase of algorithm complexity
The increase of algorithm complexity mainly includes computational complexity (time complexity) and design complexity.
No matter Pareto dominance selection mechanism, or objective decomposition selection mechanism, the computational cost of all algorithms is directly related to the objective dimension. That is, the time complexity of the algorithm is proportional to or exponentially increases with the objective dimension. Although some algorithms propose redundant objective reduction or dimensionality reduction methods, this reduction also causes an additional computational cost.
In addition, as the objective dimension increases, some algorithms require additional parameters to provide a basis for individual selection or to control the proportion of exploration in the evolutionary process. Both IBEA [27] and HypE [2] based on the indicator require a set of reference point in advance as the calculation basis for the evaluation index. However, the selection of the reference point greatly affects the evaluation of the algorithm to the individual.
This paper is organized as follows. Section 2 gives the discussion on four types of classical methods to solve the many-objective optimization problem. Section 3 gives many-objective individual reservation strategy based on improved reference vector. Based on the above study, many-objective optimization algorithm based on reference vector is proposed. Numerical simulation and result analysis are listed in Section 4.
There are mainly four types of methods to solve the many-objective optimization problem.
Objective dimension reduction method
The method, which obtains the connection between multiple objective by statistics, omits the redundant objectives or merges the associated objectives. The objective dimension reduction method can effectively solve the problem of redundant relationship between objectives. But the method may be less effective when these objectives have relatively complex associations.
The classical method of objective dimension reduction is Principal Component Analysis (PCA). When studying mathematical problems involving high dimensional and multi-objective, the number of variables is usually used to judge the difficulty of the problem. The difficulty of the problem is exponentially related to the number of variable. However, there is a correlation between variables, such as overlap or redundancy of variables. The core of PCA is to use a small number of variables to express the relationship between more variables, and to sort the variables according to their importance.
Zheng and Shen proposed an evolutionary algorithm based on information separation(ISEA) [20, 26]. The first conversation coordinate value of the algorithm is defined as convergence information, and the remaining coordinates represent individual distribution information. The coordinate information conversion method is used to separate the individual convergence information and distribution information. Then, the algorithm adopts a neighborhood penalty mechanism based on hierarchical selection to prevent neighboring individuals from being selected into the archive set at the same time. In 2005, Deb [7] proposed PCA-NSGA-II which applies the PCA principal component analysis method to the many-objective optimization problem to achieve the goal reduction effect. The algorithm obtains the non-dominated solution set through NSGA-II, and uses the data in the non-dominated solution set as the basis of PCA analysis to obtain the relationship matrix between the objective, so as to achieve the purpose of deleting the objective number.
In essence, principal component analysis and correlation analysis based on problem characteristics are based on the mathematical model to simplify the optimization problem. Both the objective reduction and the variable reduction are performed with the exception of certain variable information. Errors caused by mathematically simplified calculations accumulate during evolution and may lead to inaccuracy of the algorithm [13].
New Pareto dominance method
For Pareto non-dominated sorting of many-objective optimization problem, the comparison cost between the two solutions will increase as the number of objective increases. In the case that the number of objectives reaches a certain number, there are few individuals which are dominated by the other. Therefore, the traditional Pareto sorting method is inefficient and can’t even compare individual.
To improve the selection pressure of Pareto dominance, relaxed Pareto dominance is used the enhance the convergence ability of the algorithm. ε-dominate amplifies the dominant area of the individual by panning objective vector. Elite individuals simultaneously dominate several solutions around them, reducing local convergence caused by individual aggregation. In ε-MOEA, the grid is used to calculate the individual’s dominant range, and the number of individuals entering the population is controlled according to the dominating parameters [6]. But the algorithm is sensitive to parameter selection and control. Excessive dominance parameters cause a decrease in the number of individuals in the population, and useful population information is eliminated. If the dominance parameter is too small, the number of individuals in the population will increase [19]. When the total number of evaluations is constant, the number of iterations of the algorithm will be affected.
Wang et al. [22] proposed a non-dominated sorting method based on angular solution. In the process of comparing different objective, the algorithm outputs its non-dominated relationship when encountering contradictions.
Method based on objective decomposition
The method based on object decomposition is also called aggregation-based algorithm. These methods convert the individuals objective vector into a single scalar as the individuals fitness. When Pareto dominates fail to solve problems with many-objective optimization problems, the aggregation-based method shows its superiority. Among them, the most representative is MOEA/D [24]. The early MOEA/D mainly includes three weight vector design methods: Weighted Sum method, Tchebycheff method and Penalty Boundary Intersection method. Compared with Pareto, this kind of performance is less affected by objective number. The selection pressure provided by weighting method is not sensitive to objective number, but the design complexity of the weight vector increases with objective number.
Decomposition techniques based on reference points or reference vectors are the most popular solution to solve the many-objective optimization problems. MOEA/D-SAS adopts decomposition-based individual ordering and angle-based individual selection, and these methods improve the convergence and diversity of evolutionary population [3]. Jiang et al. [16] revive an early developed and computationally expensive strength Pareto-based evolutionary algorithm by introducing an efficient reference direction based density estimator. Different from the traditional evolutionary idea, the algorithm introduces the selection mechanism of diversity-first-and-convergence-second to improve the pressure of non-dominated solution set selection. He and Yen gave a choice strategy based on coordinate system and environment. The algorithm considers the quality of the parent individual and the utility of the parent mix to achieve the function of improving selection and recombination efficiency [8].
Method based on combination technique
For many-objective optimization problems which objective can not be reduced, the efficiency of non-dominated sorting method is low, and the weight vector of objective decomposition method should be designed elaborately based on objective characteristic. Focusing on these issues, method based on combination technique are proposed. Convergence solution set (CA) and a diversity solution set (DA) are used to design a two-archives-based algorithm(Two-Arch2) [21]. Among them, the dual solution set uses different selection mechanisms, that is, CA is based on indicators, and DA is based on Pareto. NSGA-III uses Pareto dominance to guide population to the ture front, and uses reference vector to maintain population diversity [5]. NSGA-III emphasizes characteristics of the individuals which are Pareto non-dominated and close to reference vector. For the optimization problem whose objective number reaches a certain level, the Pareto dominance relationship does not provide sufficient selection pressure, and the algorithm emphasizes diversity. Then, the combination of different convergence and diversity technique is a promising direction to solve many-objective optimization problem [11].
Method based on evaluation indicators
Using GD and IGD as evaluation indicators to calculate the performance of the solution set, ture Pareto front needs to be given in advance, but the ture Pareto front of most engineering problem cannot be obtained before the calculation. Some algorithms use the reference front as an alternative to the true Pareto. Since the position of the reference leading edge is related to the solution set obtained by the algorithm, the replacement of reference front needs to be improved [23].
Li and Tang propose a stochastic ranking-based multi-indicator algorithm (SRA) [17]. For the problem which can not be evaluated comprehensively by a single indicator, they propose a technique to balance the influence of multiple indicators, and the elite achieve with direction is introduced to maintain diversity.
Many-objective individual reservation strategy based on improved reference vector
Selection of reference vector
Based on the analysis of existing many-objective optimization methods, it can be found that the method based on reference vector and reference point perform better than dominance-based method. With the increase of objective number, selection pressure provided by dominance relationship deceases significantly. On the other hand, selection pressure provided by reference-based on method is affected by the reference distribution at the algorithm design stage. When dealing with optimization problems with complicated front, reference vector with fixed weight did not perform well. In order to solve the defect of reference vector, an individual retention strategy based on reference vector and neighborhood is proposed.
Reference vector decomposition of many-objectives optimization and individual selection are shown in Fig. 2. The generation method of reference vector is as shown in Equation 2:

Reference vector decomposition of many-objectives optimization and individual selection.
The reference vector is normalized according to Equation 4:
For individuals in the population, the normalization method is first used as in Equation 5:
For individuals in the population, the normalization method is used as Equation 6:
In the traditional dominance-based multi-objective optimization algorithm, some solutions with better convergence may be preferred due to its specific selection mechanism, as shown in Fig. 3(a). However, from the perspective of algorithm diversity, a part of the "isolated solution" which is far from the reference vector may have individual information that leads the algorithm into the optimal front in the future evolution process. When designing the selection mechanism of multi-objective evolutionary algorithm, the convergence of the individual and the guiding effect of the individual on the population should be considered simultaneously. Therefore, the dominated-individuals with the which is beneficial to the population evolution information should be given priority.

Individual selection mechanism based on reference vectors.
Figure 3 shows two special cases of individual selection mechanisms based on reference vectors. For the case of sub-figure a) of Fig. 3, since the solutions X1, X2 are closer to the origin, the formula for calculating the fitness function based on the decomposition algorithm usually obtains a better value. Therefore, solution X3 is farther from the origin, and it will be usually removed. In this case, the selection mechanism performs good convergence. As the number of objectives increases, the diversity reservation for population exhibits an equally important position to convergence. Solution X1, X2 which perform better when calculated the fitness function on convergence, may contain the same evolutionary information which may lead the evolutionary population to local area. The selection mechanism which choose these solutions at the same time may cause the deterioration on evolutionary information, especially on diversity information in objective space.
In order to avoid the situation in Fig. 3(a), some algorithms adopt the method based on angle of neighborhood. This method sort individuals according to the angle between the solution vector and the reference vector firstly, and then sort the individual by distance sorting. This method alleviates the situation of repeated individuals to a certain extent, but still cannot solve the problem that individuals with valid information are deleted. For the case of Fig. 3(b), the distances from the solution X1, X2, X3 and the origin are relatively close. Since X1 is closer to the reference vector v1 and X2 is closer to the reference vector v2, the solution X3 is usually removed in the comparison. If only the angle between the solution vector and the reference vector is considered, the individuals below the reference vector will all be deleted.
In the above two cases, it is desirable to use the solution X1 and X2 when calculating the performance of the solution set; and use the solution X2 and X3 in individual generation for evolutionary process. Therefore, two individual selection mechanisms based on vector angle and weight decomposition are designed to elite solution sets and evolutionary populations, respectively.
Figure 4 shows the dual solution set selection mechanism based on the reference vector. In order to make better use of the global information given by the reference vector and the solution information given by individual, two different angles are defined, namely the angle θi,j between the individual i and the vector j and the angle φi,j between the individual i and the nearest individual j. When selecting elite individuals used to evaluate the ability of evolutionary algorithms, angles θi,j and individual positions are used as evaluation criteria. When selecting the parent for evolutionary computing, θi,j and φi,j are used together as a selection mechanism. The specific method is shown in Equation 7.

Dual mechanism of solution set preservation based on reference vectors.
The individual fitness value of evolutionary population parent is shown in Equation 7:
According to the Kimura neutral Hypothesis, genetic variation in a body during the evolution of the entire population usually does not affect the performance of the entire population. In evolutionary algorithms, especially in the later stages of evolutionary processes, most newly generated individuals cannot join elite populations. In the following cases, the non-directionality of the evolutionary algorithm is more obvious: (1) For some multi-modal or multiple local optimization problems, the random search algorithm may perform better in preventing the population from falling into local convergence. (2) For individuals who are far away from the objective space, the individuals randomly generated by the parent selection tend to be poor fitness values, which is mainly due to the one-to-one correspondence between the objective space and the decision space. (3) In later stage of an evolutionary algorithm, some individuals are often selected as parent individuals because they perform well in the objective space. Therefore, individuals in a population tend to focus on a particular region or a few specific values, which leads to the stagnation of evolutionary algorithms.
Therefore, using the entire population as a parent selection pool is not conducive to the direction search of the evolution process. Especially in the case of an increase in the number of objectives, the generation of new individuals requires information from individuals with good evolutionary genes. At the same time, traditional evolutionary operators are optimized for single-objective optimization. There is a certain limitation when dealing with multi-objective optimization problems, especially many-objective optimization problems.
Figure 5 shows the new individual locations that can be generated by two common crossover operators when randomly selecting a parent individual. The solid line in the figure is the true Pareto front, and the dotted line is the area where the individual is probably generated. When the position of the parent individual is not good, it may occur that no feasible solution or elite solution can be generated no matter what algorithm parameters is given.

Two crossover operators.
Based on the analysis of the above phenomenon, the individual generation strategy based on probabilistic model is adopted as the evolutionary operator. Its characteristics and principles are shown in Fig. 6.

The Pareto front and Pareto set based on Gaussian mixture model.
Figure 6 shows the Pareto front and the Pareto solution set based on the Gaussian mixture model. Usually, the solutions on the true Pareto front set show a one-to-one correspondence with the solution on the true Pareto solution set. Several points of continuous distribution on the leading edge exhibit connectivity and aggregation in the solution set. According to this characteristic, the Gaussian mixture model is adopted as one of the evolutionary operators.
The fitness value calculation formula of evolutionary population parent individual is shown in Equation 8:
The individuals in the sub-populations follow the normal distribution model as shown in Equation 9:
RMaEA algorithm adopts a region decomposition method based on reference vectors, and introduces a double decomposition selection criterion to update individuals in evolutionary populations and elite populations. The algorithm updates the individual information in the population by mixing Gaussian models and differential evolution operators. The pseudo code of the RMaEA algorithm is listed in Table 1.
Pseudo code of RMaEA
Pseudo code of RMaEA
Simulation settings
The classical test suite DTLZ is selected as test function. The DTLZ test function and true front are set based on the original literature. Since the true front of DTLZ5 and DTLZ6 are difficult to obtain when the objective is greater than 3, these two functions are not used as test functions in this paper. The objective functions are 5, 10, and 15, respectively. For comparison, the population of all algorithms is fine-tuned according to the population size of NSGA-III. The population size are 212, 276, and 136 for the objective numbers with 5, 10, and 15, respectively. All test functions were terminated with 100,000 evaluations and independently run 30 times.
In order to compare the algorithms proposed in this paper, the classical indicator-based evolutionary algorithm IBEA, the decomposition-based evolutionary algorithm MOEA/D, DBEA and NSGA-III proposed are selected. The technical characteristics and references of the algorithm are shown in Table 2. Detailed comparison algorithm parameters refer to the original literature.
Characters and references of comparison algorithms
Characters and references of comparison algorithms
For each standard test function, RMaEA, DBEA, NSGA-III, MOEA/D, and IBEA are run independently 30 times, and the Pareto front solution set solved for each run is collected. Table 3 lists the mean and standard deviation of the IGD for each run of 30 algorithms. The minimum values of each function in Table 3 are shown in bold, and the number of minimum values obtained by each algorithm is counted in the last row of the table.
Statistics result of IGD with 30 runs
In Table 3, RMaEA obtains 13 minimum mean values and 8 minimum standard deviations among the 15 test problems. The good performance of RMaEA proves the ability of the proposed algorithm on the problems with more than 5 objectives. Numerical analysis shows that the NSGA-III algorithm is close to the proposed algorithm on the IGD mean value. The indicator-based algorithm IBEA and the reference point-based algorithm DBEA do not obtain the minimum mean and standard deviation of any test functions. In order to visually represent the convergence performance of the five algorithms, Figs. 7–11 show box diagrams of IGD value for five test algorithms with different objective numbers.

Boxplot of 5 algorithms on test problem DTLZ1.

Boxplot of 5 algorithms on test problem DTLZ2.

Boxplot of 5 algorithms on test problem DTLZ3.

Boxplot of 5 algorithms on test problem DTLZ4.

Boxplot of 5 algorithms on test problem DTLZ7.
The IGD values of RMaEA is always close to the algorithm NSGA-III on most of test functions except that the RMaEA exhibits some outliers. Reference vector based strategy and reference point based strategy achieve similar result on DTLZ test suite. Usually, as the number of objective increases, the difficulty of the problem increases. So, the IGD value of each algorithm gradually increases. So, the time complexity of the algorithm is gradually improved, and some traditional multi-objective optimization techniques fail to achieve promising results.
For the test function DTLZ1, RMaEA and NSGA-III obtain the best mean and stability on IGD indicator. Compared with NSGA-III, RMaEA achieves fewer abnormal points. For the test function DTLZ2, RMaEA achieves the best mean and stability, followed by NSGA-III. With the increase of objective number, the performance of DBEA gradually improves and becomes closed to RMaEA and NSGA-III.
For the test function DTLZ3, the rank of the algorithm performance is RMaEA, NSGA-III, MOEA/D, DBEA and IBEA. For the test function DTLZ4, RMaEA and NSGA-III perform very close except for a few abnormal points. As the objective dimension increases, advantages of DBEA gradually emerge. For the test function DTLZ7, the five algorithms show different performances under different objective numbers. The reason is that the leading edge of the test function DTLZ7 is composed of several discontinuous regions.
Diagram analysis of simulation results
In order to explicitly illustrate the performance of these algorithms, the parallel coordinate system of the comparison algorithm on each test function is plotted using the experimental data of the minimum IGD in 30 runs, as shown in Figs. 12–14 and Fig. S-1 to Fig. S-12 in supplementary material. Since DBEA and IBEA did not achieve satisfactory performance on IGD values, this section no longer draws the solution set for these algorithms. From the perspective of the ordinate, RMaEA achieves more uniform distribution on most of many-objective optimization problems.

Solution sets of RMaEA, NSGA-III and MOEA/D on test problem DTLZ1 with 5 objectives.

Solution sets of RMaEA, NSGA-III and MOEA/D on test problem DTLZ1 with 10 objectives.

Solution sets of RMaEA, NSGA-III and MOEA/D on test problem DTLZ1 with 15 objectives.
Figures 12–14 show the solution sets of the test function DLTZ1 for the RMaEA, NSGA-III, and MOEA/D algorithms with the objective number of 5, 10 and 15, respectively. The performance of the three algorithms in different objective functions is more obvious. As the number of objective increases, the performance advantages of RMAEA and NSGA-III algorithms gradually emerge. The NSGA-III showed anomalies in the case with 15 objectives. In the case of MOEA/D with 15 objectives, some objective do not converge to the true front.
For the the solution sets of the test function DLTZ2, DLTZ3, DLTZ4, DLTZ7 for the RMaEA, NSGA-III, and MOEA/D algorithms with the objective number of 5, 10 and 15, please refer to the supplementary material.
Based on the analysis of the difficulties of existing many-objective optimization problems, a many-objective evolutionary algorithm based on adaptive reference vector (RMaEA) is proposed along with diversity retention strategy. Focusing on solving the failure of traditional non-dominated sort method on the selection solution when objective increases, weight decomposition method is introduced. To adapt different shape of Pareto front, reference-vector adaptive method is proposed. Under the premise that different reference vector selection mechanism are suitable for different Pareto front, a dual reference selection mechanism with adaptive dynamic adjustment is designed. The algorithm uses the cosine of the nearest angle as the diversity penalty value, and uses the cosine of the angle between the individual and the reference vector and the position of the individual’s distance origin as the basis for selection. The dual-reference selection mechanism uses different selections to select evolutionary and elite populations, and exchanges information during evolution to maintain population diversity. The simulation results show that the Pareto front solution set solved by the proposed algorithm performs a good balance between convergence and diversity.
In future research, reference-vector adaptive method will be modified for complicated many-objective optimization problem. The reference-dominated rank may be a promising way to solve the many-objective problem with the linkage between multi-variable. Furthermore, the application of many-objective optimization method will be done on rolling schedule problem for cold tandem rolling.
Footnotes
Acknowledgement
This work was supported by National Natural Science Foundation of China No.62003296 and Grant No.61703361, the Natural Science Foundation of Hebei No.F2020203031 and Grant No.E2018203162, the Post-Doctoral Research Projects of Hebei B2019003021, the Doctoral Foundation of Yanshan University No.BL18048, the Project Supported by Natural Science Foundation Steel and Iron Foundation of Hebei Province E2019105123, the Project Supported by Hebei Education Department ZD2019311.
