Differential Evolution in Robust Optimization Over Time Using a Survival Time Approach

: This study presents an empirical comparison of the standard differential evolution (DE) against three random sampling methods to solve robust optimization over time problems with a survival time approach to analyze its viability and performance capacity of solving problems in dynamic environments. A set of instances with four different dynamics, generated by two different conﬁgurations of two well-known benchmarks, are solved. This work also introduces a comparison criterion that allows the algorithm to discriminate among solutions with similar survival times to beneﬁt the selection process. The results show that the standard DE holds a good performance to ﬁnd ROOT solutions, improving the results reported by state-of-the-art approaches in the studied environments. Finally, it was found that the chaotic dynamic, disregarding the type of peak movement in the search space, is a source of difﬁculty for the proposed DE algorithm.


Introduction
Optimization is an inherent process in various areas of study and everyday life. The search to improve processes, services, and performances has originated in different solution techniques. However, there are problems in which uncertainty is present over time, given that the solution's environment can change at a specific time. These types of problems are named Dynamic Optimization Problems (DOPs) [1]. This study deals with dynamic problems where the environment of the problem changes over time. Various studies have been carried out to resolve DOPs through tracking moving optima (TMO), which is characterized by the search for and implementation of the global optimal-solution every time the environment changes [2][3][4].
Evolutionary algorithms, such as Differential Evolution (DE), have shown good performance to solve tracking problems [5][6][7]. However, the search and implementation of the optimum each time the environment changes may not be feasible due to different circumstances, such as time or cost.
The approach introduced in [8] tries to solve DOPs through a procedure known as robust optimization over time (ROOT). ROOT seeks to solve DOPs by looking for a good solution for multiple environments and preserve it for as long as possible, while its quality does not decrease from a pre-established threshold. The solution found is called Robust Solution Over Time (RSOT).
In this regard, Fu et al. [9] introduced different measures to characterize environmental changes. After that, the authors developed two definitions for robustness [10]. The first one was based on "Survival Time"-when a solution is considered acceptable (an aptitude threshold must be previously defined). The second definition was based on the "Average Fitness"-the solution's average fitness is maintained during a previously defined time window. The measurements incorporate information on the concepts of robustness and consider the values of error estimators. An algorithm performance measure was suggested to find ROOT solutions. The study was carried out using a modified version of the moving peaks benchmark (mMPB) [10].
Jin et al. [11] proposed a ROOT framework that includes three variants of the Particle Swarm Optimization algorithm (PSO). PSO with a simple restart strategy (sPSO), PSO with memory scheme (memPSO), and a variant that implements the species technique SPSO. The authors applied a radial basis function as an approximation and an autoregressive model as a predictor.
On the other hand, Wang introduced the concept of robustness in a multi-objective environment, where a framework is created to find robust Pareto fronts [12]. The author adopted the dynamic multi-objective evolutionary optimization algorithm in the experiments. At the same time, Huang et al. considered the cost derived from implementing new solutions, thus addressing the ROOT problem by a multi-objective PSO (MOPSO); The Fu metric was applied in that study [13].
Yazdani et al. introduced a new semi-ROOT algorithm that looks for a new solution when the current one is not acceptable, or if the current one is acceptable but the algorithm finds a better solution, whose implementation is preferable even with the cost of change [14].
Novoa-Hernández, Pelta, and Corona analyzed the ROOT behavior using some approximation models [15]. The authors suggested that the radial basis network model with radial basis function works better for problems with a low number of peaks. However, considering all the scenarios, the SVM model with Laplace Kernel shows notably better performance to those compared in the tests carried out.
Novoa-Hernández and Amilkar in [16] reviewed different relevant contributions to ROOT. The authors analyzed papers hosted in the SCOPUS database. Concerning new methods to solve ROOT problems, Yazdani, Nguyen, and Branke proposed a new framework using a multi-population approach where sub-populations track peaks and collect information from them [17]. Adam and Yao introduced three methods to address ROOT (Mesh, Optimal in time, and Robust). The authors mentioned that they significantly improves the results obtained for ROOT in the state-of-the-art [18]. Fox, Yang and Caraffini studied different prediction methods in ROOT, including the Linear and Quadratic Regression methods, an Autoregressive model, and Support Vector Regression [19]. Finally, Liu and Liang mapped a ROOT approach to minimize the electric bus transit center's total cost in the first stage [20].
In different studies, DE has been used to solve ROOT problems using the Average Fitness approach, achieving competitive results [21,22]. However, to the best of the authors' knowledge, there are no studies that determine DE's performance in solving ROOT problems with the Survival Time approach, and this is where this work precisely focuses. This research aims to present an empirical comparison of the standard DE against three random sampling methods to solve robust optimization over time problems with a survival time approach to analyze its viability and performance capacity of solving problems in dynamic environments.
The paper is organized as follows: Section 2 includes ROOT's definition under a survival time approach, while in Section 3 the implemented methods based on random sampling are detailed. Section 4 details the standard differential evolution and the objective function used by the algorithm in the present study. Section 5 specifies the benchmark problems to be solved. After that, Section 6 specifies the experimental settings and Section 7 shows the results obtained. Finally, Section 8 summarizes the conclusions and future work.

Survival Time Approach
Under this approach, a threshold is predefined to specify the quality that a solution must have to be considered good or suitable to survive. Once the threshold is defined, the search begins for a solution whose fitness can remain above the threshold in as many environments as possible. In this sense, the solution is maintained until its quality does not meet the predefined expectations, and then new robust solution over time must be sought.
In Equation (1), the function F s ( x, t, V) to calculate the survival time fitness of a solution x at time t is detailed. It measures the number of environments that a solution remains in above the threshold V.

Random Sampling Methods
In the study presented in [18] the authors proposed three random sampling methods to solve ROOT problems, with a better performance against the state-of-the-art algorithms. The methods are described below.
In all three methods, the best solution should be searched in the current time's solutions space, modifying the solution space when the "Robust method" is used and then using that solution next time according to the approach used (Survival Time or Average Fitness).

Mesh
This method performs random sampling in the current search space, then uses the sample with the best fitness in the current environment as a robust solution over time, using the solution found in the following times ( Figure 1).

Time-Optimal
This method performs a search similar to the "Mesh" method, with the difference that the best solution found being improved using a local search ( Figure 2). Time-optimal method. The green point is the robust solution found by the method, which is used in the next times.

Robust
This method performs a search similar to the "Mesh" method, differing in that a smoothing preprocessing of the solution space is performed before the search process. As seen in Figure 3, the solution obtained (green dot) can vary concerning the solution with better suitability in the raw environment (yellow dot).

Differential Evolution
In 1995, Storm and Price proposed an evolutionary algorithm to solve optimization problems in continuous spaces. DE is based on a population of solutions that, through simple recombination and mutation, evolves, thus improving individuals' fitness [23].
Considering the fact that this work, to the best of the authors' knowledge, is the first attempt to study DE in this type of ROOT problems (i.e., survival time), and also taking into account that the most popular DE variant (DE/rand/1/bin) has provided competitive results in ROOT problems under an average fitness approach [21,22], the algorithm used in this study is precisely the most popular variant known as DE/rand/1/bin, where ''rand" (random) refers to the base vector used in the mutation, ''bin" (binomial) refers to the crossover type used, and 1 means one vector difference computed.
The algorithm starts by randomly generating a uniformly distributed population x i , G ∈ i = 1, ..., NP, where NP is the number of individuals for each generation ''G".
After that, the algorithm enters an evolution cycle until the stop condition is reached. We applied the maximum number of evaluations allowed ''MAX Eval " as the stop condition.
Subsequently, to adapt individuals, the algorithm performs recombination, mutation, and the replacement of each one in the current generation. One of the most popular mutations is DE/rand/1 in Equation (2), where r0 = r1 = r2 = i are the indices of individuals randomly chosen from the population, 1 is the number of differences used in the mutation and F > 0 is the scale factor.
The vector obtained v i,G is known as mutant vector, which is recombined with the target (parent) vector by binomial crossover, as detailed in Equation (3).
In this study, the elements of the child vector (also called trial) u i,j,G are limited according to the pre-established maximum and minimum limits, also known as boundary constraints. Based on the study in [24], we use the boundary method (see Algorithm 1,line 14). In the selection process, the algorithm determines the vector that will prevail for the next generation between parent (target) and child (trial), as expressed in Equation (4).
Algorithm 1: "DE/rand/1/bin" Algorithm. NP, MAX Eval , CR and F are parameters defined by the user. D is the dimension of the problem.
Randomly select r0 = r1 = r2 = i: In Equation (1), the function to obtain an individual's fitness through the survival time approach is shown. However, the fitness obtained is not enough to differentiate similar individuals, i.e., individuals who have survived the same amount of environments. That is why, in the implemented algorithm, we consider an additional calculation to help identify better solutions.
We propose to obtain the average of the solution's quality in the environments that have survived. This average value helps to differentiate solutions with similar survival times. Therefore, the objective function now considers both, the number of surviving environments and the performance achieved by this solution in those environments that it has survived.
Considering the fact that the maximum height of the peaks is defined at 70 (see Table 3), the objective function for the implemented algorithm is given by the result obtained in Equation (1) multiplied by 100 plus the average fitness of the solution throughout the environments it has survived.

Benchmark Problems
The problems tackled in this study are based on Moving Peaks Benchmark (MPB) [25] and are configured in a similar way to that used in various specialized literature publications on ROOT, and specifically as used in [18].
Two modified MPBs can be highlighted, which are described in the following subsections. The dynamics used are presented in Table 1, where ∆φ is the increment from time t to time t + 1 of the φ parameter.

Moving Peaks Benchmark 1 (MPB1)
In this benchmark, environments with conical peaks of height h(t) ∈ [h min , h max ], width w(t) ∈ [w min , w max ] and center c(t) ∈ [x min , x max ] are generated, where the design variable x is bounded in [x min , x max ]. The function to generate the environment is expressed in Equation (5), where the dynamic function for height and width is given as in Table 1, while the center moves according to Equation (6). r i follows an uniform distribution of a D-dimensional sphere of radius s i , and λ ∈ [0, 1] is a fixed parameter.
In the present study, two problems generated by this benchmark are solved, with λ = 0 it implies that the movement of the peaks is random, while with λ = 1 it implies that the movement is constant in the direction v i (t).

Moving Peaks Benchmark 2 (MPB2)
The set of test functions in this benchmark is described in Equation (7), where a(t) is the environment at time step t, h i (t), w i (t), c i (t) is the height, width and center of the i-th peak function at time t, respectively; x is the decision variable and m is the total number of peaks. h i (t + 1) and w i (t + 1) vary according to Table 1. An additional technique that uses a rotation matrix is used to rotate the centers [25].

Experimental Settings
Based on the information in Section 5, different environments are generated as test problems and they are summarized in Table 2. Table 2. Summary of problems.

Benchmark Abbreviation Configuration
Dynamic The parameter settings of the problems are detailed in Table 3. The height and width of the peaks were randomly initialized in the predefined ranges. The centers were randomly initialized within the solution space.
A survival threshold V = 50 is selected, representing the most difficult cases that have been resolved in the literature under the survival approach. The higher the survival threshold, the more difficult it is to find solutions that satisfy multiple scenarios.
The DE parameters were fine-tuned using the iRace tool [26] and they are summarized in Table 4, where NP is the population size, CR is the crossover parameter, and F is the scale factor. For each problem, a solution is sought at each time t ∈ (2, ..., 100). In order to evaluate an RSOT in a specific time, approximate and predictive methods have been used in the literature so that the performance of an algorithm depends on their accuracy. However, in this study, we want to know the DE behavior when solving the ROOT environments considering they had ideal predictors to evaluate the solutions. In this regard, the process to study the algorithm's ability to find RSOT using DE at each instant of time is as follows: • A solution is sought according to the algorithm described in Section 4, and the measured solution value by Equation (1) is recorded. • Subsequently, to obtain the algorithm's performance in the following environment, the search process is performed again using the real-environments; the best solution found is newly measured by Equation (1) and is also recorded. • The described procedure is carried out at each instant of time that is being recorded. Therefore, in the present study, it is not necessary to detect environmental changes to know at what point in time a solution is no longer considered good. Each time a solution is sought, the algorithm initializes its population randomly, avoiding diversity problems.

Results and Discussion
The results for the problems generated with dynamics 1-4 are detailed in Table 7 and graphically shown in Figures 6-9, for each one of the four dynamics. In all four figures, those labeled with (a) and (b) present the results obtained in the MPB1 problems, while those labeled with (c) and (d) refer to the MPB2 problems. In all cases, the average survival values obtained by Mesh, time-optimal and robust approaches are compared against DE.
Non-parametric statistical tests [27] were applied to the corresponding numerical results presented in Table 7. The 95%-confidence Kruskal-Wallis and 95%-confidence Friedman tests were applied and their obtained p-values are reported in Table 5. Table 5. Results of the 95%-confidence Kruskal-Wallis (KW) and Friedman (F) tests. The symbol (*) after letter ''D" in the Problem column refers to the type of dynamic used according to columns Dynamic. A p-value less than 0.05 means that there are significant differences among the compared algorithms in such problems. To further determine differences among the compared algorithms, the 95%-confidence Wilcoxon test was applied to pair-wise comparisons for each problem instance. The obtained p-values are reported in Table 6, where the significant improvement with a significance level α = 0.5 is shown in boldface. We can observe that the Wilcoxon test confirmed significant differences obtained in Kruskal-Wallis and Friedman Tests, most of them comparing DE/rand/1/bin versus random sampling methods, with the exception of four corresponding to problems generated with dynamic 4 (B1D4-1 and B1D4-2 both, in the comparison of DE/rand/1/bin versus Mesh and DE/rand/1/bin versus Time-optimal).  Table 7 summarizes the mean and standard deviation statistical results obtained by the compared algorithms. It can be seen that DE/rand/1/bin obtains the highest average values in all the problems that were solved. Nevertheless, the higher standard deviation values obtained by DE/rand/1/bin in problems B1D4-1 and B1D4-2 confirm those expressed by the non-parametric tests-the differences are not significant with respect to the random sampling methods. Figures 4 and 5 have the box-plots for B1D4-1 and B1D4-2, where all the compared algorithms reach survival times between 1 and 3, with the exception of the Robust approach in B1D4-1, but such a difference was not significant.  With respect to the graphical results, when MPB1 (items (a) and (b)) is compared against MPB2 (items (c) and (d)) in Figures 6-9, it is clear that MPB1 is more difficult to solve by all four approaches. However, in all cases (MPB1 and MPB2 in the four dynamics) DE is able to provide better results against the three other algorithms. Such a performance is more evident in all MPB2 instances.

p-Value
Regarding MPB1 (items (a) and (b)), it is important to note that it is more difficult to find higher survival times when the peak movement is random (items (a), where λ = 0).
Another interesting behavior found is that all four compared methods are affected mainly by the random and chaotic dynamics in those MPB1 instances, the latter one being the most complex (chaotic dynamic). However, even in such a case DE was able to match and in some cases improve the survival values of the compared approaches. This source of difficulty now found motivates part of our future research.

Conclusions
A performance analysis of the differential evolution algorithm, with one of its original variants, called DE/rand/1/bin, when solving robust optimization over time problems with a survival time approach, was presented in this paper. Three state-of-the-art random sampling methods to solve ROOT problems were used for comparison purposes. Sixteen generated problems by two benchmarks with two configurations and four different dynamics were solved. The solutions generated by the DE were obtained using the real environments without prediction mechanisms with the aim to analyze its behavior in ideal conditions. The findings supported by the obtained results indicate that DE is a suitable algorithm to deal with this type of dynamic search space when a survival time approach is considered. Moreover, the additional criterion that was added to the DE objective function allowed the algorithm to better discriminate between similar solutions in terms of survival time. Furthermore, it was found that the combination of a chaotic dynamic with both, random and constant peak movements, is a source of difficulty that requires further analysis.
This last finding is the starting point of our future research, where more recent DE variants, such as DE/current-to-p-best, will be tested in those complex ROOT instances. Moreover, the effect of predictors in DE-based approaches will be studied. The first author acknowledges support from the Mexican National Council of Science and Technology (CONACyT) through a scholarship to pursue graduate studies at the University of Veracruz.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript:

DE
Differential Evolution MPB1 Moving Peaks Benchmark 1 MPB2 Moving Peaks Benchmark 2 ROOT Robust Optimization over Time RSOT Robust Solution over Time S.D.
Standard deviation