Particle Swarm Optimization Based on a Novel Evaluation of Diversity

: In this paper, we propose a particle swarm optimization variant based on a novel evaluation of diversity (PSO-ED). By a novel encoding of the sub-space of the search space and the hash table technique, the diversity of the swarm can be evaluated efﬁciently without any information compression. This paper proposes a notion of exploration degree based on the diversity of the swarm in the exploration, exploitation, and convergence states to characterize the degree of demand for the dispersion of the swarm. Further, a disturbance update mode is proposed to help the particles jump to the promising regions while reducing the cost of function evaluations for poor particles. The effectiveness of PSO-ED is validated on the CEC2015 test suite by comparison with seven popular PSO variants out of 12 benchmark functions; PSO-ED achieves six best results for both 10- D and 30- D .


Introduction
Inspired by the emergent motion of the foraging behavior of a flock of birds in nature, the particle swarm optimization (PSO) was first proposed by Kennedy and Eberhart [1,2] in 1995 to solve the continuous nonlinear optimization problems. Compared with other evolutionary algorithms, the PSO has attracted much attention since it was proposed because of its fewer control parameters and better convergence. Nowadays, PSO has been widely used in many fields, such as communication networks [3,4], medicine engineering [5][6][7], task scheduling [8,9], energy management [10], linguistics studies [11], supply chain management [12] and neural networks [13,14].
Despite its robustness for solving complex optimization problems, rapid convergence may cause the swarm to be easily trapped into some local optima when solving multimodal problems via PSO [15][16][17][18][19][20][21][22]. Therefore, reasonably using the swarm's exploration ability (global investigation of the search place) and exploitation ability (finer search around a local optimum) is a crucial factor for PSO's success, especially for complex multimodal problems having a large number of local minima. For this purpose, several PSO variants have been reported. In the following, we shall discuss three main types: PSO with adjusted parameters, PSO with various topology structures, and PSO with hybrid strategies.
PSO with adjusted parameters: It is evident that appropriate control parameters, such as the inertia weight ω and the acceleration coefficients c 1 and c 2 , have significant effects on the exploration and exploitation abilities of the swarm. By increasing the value of ω significantly diversity of the swarm increases while increasing the values of c 1 and c 2 accelerates the particles towards the historical optimal position or the optimal position of the whole swarm. Zhan [15] proposed an adaptive particle swarm optimization (APSO) algorithm in 2009, which updated control parameters (ω, c 1 and c 2 ) adaptively based on the distribution of positions and fitness of the swarm. Zhang et al. [23] proposed an inertia weight adjustment scheme based on Bayesian techniques to enhance the swarm's exploitation ability. Tanweer et al. [24] proposed a PSO algorithm (SRPSO) that employed a self-regulating inertia weight strategy to the best particle to enhance the exploration ability.
Diversity, a parameter for evaluating the degree of dispersion in a swarm, is widely used to improve the PSO's performance by adaptively balancing the swarm's exploration and exploitation capabilities [17,[52][53][54][55]. When diversity is small, the swarm can explore more accurate solutions; on the other hand, it can explore a large space to avoid premature convergence when the diversity is large. Therefore, it is suggested to maintain a large diversity in the early stage of evolution to explore a vast space where the optimal solution may exist, and reduce the swarm's diversity near the end of the evolution to refine the result locally.
A proper diversity value should reflect the real characteristics of the swarm. If the diversity value is incorrectly evaluated, it may mislead the swarm's movement. However, measuring the diversity of the swarm is a somewhat challenging problem.
The distance-based diversity evaluation methods can be further divided into two categories: (1) the average distance between each pair of particles in the swarm [63]; (2) the average distance between each particle and geometric center of the particle swarm [64].
The information entropy-based diversity evaluation scheme is developed from the information theory proposed by Shannon in 1948 [65]. After that, many scholars introduced this concept into particle swarm optimization algorithm. For example, Pires [60] proposed an entropy-based index to measure the diversity of the swarm. The basic idea of this method is as follows: (1) evenly dividing the search space into Q subspace; (2) counting Algorithms 2021, 14, 29 4 of 16 the number Z q of particles in each subspace; (3) calculating the probability in the subspace by Equation (3); (4) evaluating the diversity (the information entropy) by Equation (4).
where N s represents the total number of particles, q represents the q-th subspace, and E represents the information entropy and is taken as the value of diversity.

PSO-ED
This section presents the technical details of PSO-ED. Section 3.1 introduces the main idea of PSO-ED; Section 3.2 presents our evaluation strategy of computing the exploration degree; Section 3.3 presents an adaptive update of inertia weight based on the exploration degree; Section 3.4 presents the swarm reinitialization mechanism that helps the swarm escape local traps; Section 3.5 presents two update modes, i.e., normal update mode and disturbance update mode, for saving function evaluations (FEs) on poor-performed particles.

The Main Idea of PSO-ED
Given a dimension D (D = 10 or 30 for CEC2015), the flowchart of the PSO-ED is shown in Figure 1. PSO-ED has four main modules: exploration degree evaluation, an adaptive update of inertia weight, normal and disturbance update modes, and swarm reinitialization.
As mentioned above, the core problem of particle swarm optimization algorithm is how to balance the exploration and exploitation abilities of the swarm. To address this problem, in Module 1, we develop a technique of exploration degree based on the diversity of the swarm in the exploration, exploitation, and convergence states. In module 2, with the help of exploration degree, we adaptively update the inertial weight of the swarm to ensure that the swarm has a larger exploration degree in the exploration state and a smaller one in the convergence state. In Module 3, we give two update modes (normal and disturbance update modes) for each particle according to its fitness, so as to achieve the global optimum quickly. In Module 4, to help the swarm escape from local traps, a novel swarm reinitialization scheme is proposed.
For ease of description, we summarize the parameters used in the rest of the paper in Table 1. The definition of each critical parameter is also given when it first appears in the context.

Symbol Quantity t
The index of the current iteration x i The position of particle i, and x i ∈ [x lb , x ub ] where x lb and x ub are lower and upper bounds of the position v i The velocity of particle i, where v ub is the upper bound of the velocity ω, c 1 , c 2 , c 3 The inertial weight and accelerations for PSO respectively, and c 3 is defined in Section 3.4 P i , P g , P t The position of Pbest, Gbest and Tbest, where Tbest denotes the best position over all historical rounds of evolutions N s The size of the swarm D The dimension of the problem space N total The total number of iterations

Symbol Quantity
The index of the current iteration The inertial weight and accelerations for PSO respectively, and is defined in Section 3.4

, ,
The position of Pbest, Gbest and Tbest, where Tbest denotes the best position over all historical rounds of evolutions The size of the swarm

A Novel Diversity Evaluation Scheme
As mentioned above, the distance-based diversity is evaluated by computing the distances between particles. The essence of this scheme is space compression, which compresses the D-dimensional space features into one-dimensional space. For example, assume that there are three particles in the three-dimensional space, their positions are p 1 = {1,0,0}, p 2 = {0,1,0} and p 3 = {0,0,1}, then the mutual distance between any two particles is {(p 1 , p 2 ) = 2 0.5 , (p 1 , p 3 ) = 2 0.5 , (p 2 , p 3 ) = 2 0.5 }, which has only one dimension. Observe that the three points are separated by a distance of 2 0.5 , but they are regarded as being aggregated together since the information entropy of the distance information is 0 (Equations (3) and (4)). Therefore, the information entropy based on distance cannot reflect the swarm diversity very well. On the other hand, if space is divided into multiple small subspaces, the complexity of traversing all subspaces is fairly high. For example, in the D-dimensional space, if each dimension is divided into K equal parts, then the complexity of traveling all subspaces is (O K D ), which is formidable when D is larger than 10.
To calculate the information entropy of the swarm efficiently and evaluate the degree of dispersion of the swarm as accurately as possible, we propose a novel diversity Evaluation scheme. The crucial steps of our approach are as follows: For each particle in the swarm, compute the index (ID) of the subspace where the particle is located. First of all, each dimension of the search space is divided into K equal parts. Then the ID of a subspace is a string of D numbers, each of which represents the index of the part of the dimension occupied by the subspace. For example, in the 3-dimensional space in Figure 2, the ID of the blue subspace is {2, 3, 1}. Denote = { , , … , } as the ID of subspace where the i-th particle is located, where is obtained from Equation (5): Build a hash table h by traversing the IDs obtained in the Step (1). h is composed of key-value pairs {key, value}, where key is the ID of the subspace, and value is the number of particles in the subspace. We traverse the IDs of the subspaces obtained by the above steps, if a key c exists in h its value h(c) is increased by 1, else a new key-value pair {c,1} is inserted into h.
Determine the information entropy E of the swarm based on h. Note that the empty subspace without any particle contributes 0 to the calculation of information entropy since the probability p of particles falling into this subspace is 0 by Equations (3) and (4). Therefore, only the subspace corresponding to a key of h can contribute to the information entropy.
in Equation (3) represents the probability of falling into the subspace q, and = ℎ( ) ⁄ . Therefore, E can be obtained by Equations (3) and (4) based on the keys of h and can be used in evaluating the exploration degree in Equation (8).
The time complexity of the above procedure is analyzed as follows: in step 1, we need to calculate the of each particle i in each dimension d, so its time complexity is ( ), where is the size of the swarm, and D is the dimension of the search space; in step 2, inserting a key-value pair into hash table h can be done in O(1) time, the time of traversing the IDs to establish hash table h is ( ); in the last step, because only the key of the hash table h needs to be traversed, and the number of IDs is no more than , the time com- Denote as the ID of subspace where the i-th particle is located, where c d i is obtained from Equation (5): Build a hash table h by traversing the IDs obtained in the Step (1). h is composed of key-value pairs {key, value}, where key is the ID of the subspace, and value is the number of particles in the subspace. We traverse the IDs of the subspaces obtained by the above steps, if a key c exists in h its value h(c) is increased by 1, else a new key-value pair {c,1} is inserted into h.
Determine the information entropy E of the swarm based on h. Note that the empty subspace without any particle contributes 0 to the calculation of information entropy since the probability p of particles falling into this subspace is 0 by Equations (3) and (4). Therefore, only the subspace corresponding to a key of h can contribute to the information entropy. p q in Equation (3) represents the probability of falling into the subspace q, and p i = h(c)/N s . Therefore, E can be obtained by Equations (3) and (4) based on the keys of h and can be used in evaluating the exploration degree in Equation (8). The time complexity of the above procedure is analyzed as follows: in step 1, we need to calculate the c d i of each particle i in each dimension d, so its time complexity is O(N s D), where N s is the size of the swarm, and D is the dimension of the search space; in step 2, inserting a key-value pair into hash table h can be done in O(1) time, the time of traversing the IDs to establish hash table h is O(N s ); in the last step, because only the key of the hash table h needs to be traversed, and the number of IDs is no more than N s , the time complexity is also O(N s ). In summary, the running time of the above procedures is O(N s D), which is significantly smaller than O N s D .

Switch of Evolutionary States
A conventional state transition of PSO simply set the initial state, the middle state and the last state of evolution as Exploration, Exploitation and Convergence, respectively. We are more flexible in switching the evolution states, and we allow multiple state switches during the whole evolution process. The swarm starts with the Exploration state, and its switching circle is Exploration→Exploitation→Convergence→Exploration. The process terminates as the maximum number of FEs is invoked.
As for the conditions of state switches, let N ud be the number of consecutive generations when the result of the swarm evolution has not been improved, if N ud > N umax , the state should be switched at this time to enhance the chance of achieving better results (exploration→exploitation, exploitation→convergence) or jumping out of the local trap (convergence→exploration). However, if the state switches too frequently, it may not lead to any promising result while wasting time in state switch. To address this issue, let N state be the number of generations the current state has been maintained, depending on whether when N ud > N umax , we constraint N state in a range of generations (i.e., [N smin , N smax ]) to control the switching frequency as follows: Note that the state S shall be used in the Section 3.2.3 to determine the upper and lower bound of the exploration degree.
Since the swarm may not evolve in every iteration during the evolution process, we use the notion of average evolution rate to judge whether the current iteration has evolved. As shown in Equation (6), the average evolution rate of consecutive K generations is obtained as the total evolution rate of these K generations divided by K. When the evolution rate is less than a user-defined threshold α, we increase N ud by 1, otherwise, we set N ud as 0.
where the F(t) is the best result obtained in t-generation.

Exploration Degree
Let E d denote the exploration degree, which characterizes the degree of demand for the dispersion of the swarm. By trial and error, the relationship of the upper and lower bounds of E d . The three evolutionary states can lead to fairly fine results if set, as shown in Figure 3.

Exploration Degree
Let denote the exploration degree, which characterizes the degree of demand for the dispersion of the swarm. By trial and error, the relationship of the upper and lower bounds of . The three evolutionary states can lead to fairly fine results if set, as shown in Figure 3. We use the diversity E(t) to determine the exact value of . When E(t) is high, it indicates that the swarm may still explore a wider area and the corresponding value of should be larger, and vice versa. To mimic this, we develop a linear function to express the relationship between and E(t) in Equation (8). We use the diversity E(t) to determine the exact value of E d . When E(t) is high, it indicates that the swarm may still explore a wider area and the corresponding value of E d should be larger, and vice versa. To mimic this, we develop a linear function to express the relationship between E d and E(t) in Equation (8).
In the following, we show how to realize an adaptive update of the inertial weight and, therefore, the balance of the exploration and exploitation by E d (t).

Adaptive Update of Inertia Weight
The inertia weight plays a critical role in balancing the swarm's exploration and exploitation capabilities, and it also significantly affects the result's accuracy [23][24][25]. Therefore, we adaptively adjust the swarm's inertia weight to ensure that the swarm has a larger exploration degree in the exploration stage and a smaller one in the convergence stage. Finally, we use a sigmoid-like mapping that decreases from 0.9 to 0.4 [66,67] to characterize our inertia weight adjustment (Equation (9)). It is evident from this formula that when the exploration degree Ed of the population is large (close to 1.0), the inertia weight ω will also be large (close to 0.9) to ensure that the swarm can continue to explore a wider space. In contrast, when Ed becomes small, the inertia weight is also reduced to enhance the swarm's convergence ability.

Excape Local Trap by Swarm Reintialization
Note that Module 3 (Normal and Disturbance Update Modes) relies on Tbest, which is equal to Gbest when no reinitialization is conducted, therefore we shall present Module 4 ahead of Module 3 for ease of description. When the swarm falls into the trap of local optimal and cannot rely on its own to jump out, reinitializing the swarm has been proved to be an effective method of addressing this issue [68][69][70].
A crucial issue is how to judge whether the swarm has fallen into a local optimum or not. In this article, we say that the swarm has fallen into a local optimum if the iterations of undeveloped generation of the swarm N ud exceeds a maximum allowed number N umax .
Another crucial issue is how to help the swarm jump out of a local optimum. When the swarm is trapped in a local optimum, the conventional solution is to reinitialize the swarm. If we abandon the global best (Gbest) by convention, then we waste a lot of early effort; on the other hand, if we use it, the swarm has a high probability of falling into this local trap again due to its powerful attraction. To address this issue, let Gbest denote the Algorithms 2021, 14, 29 9 of 16 best position of the swarm of the current round of evolution (with a single initialization or reinitialization), and let Tbest denotes the best position of the swarm for all historical rounds of evolutions, thus Tbest is at least as good as Gbest. Finally, we add the Tbest term to the particle update in Equation (10), which ensures that the optimal position (Tbest) has only a moderate effect on the swarm with less attraction power.
where P t is the position of Tbest.
For ease of computation, we fix c 1 = c 2 = 1.49445 based on a large set of experimental results of previous PSO variants (e.g., [23,36,71]). As for c 3 , since its role is to guide the swarm to the position of Tbest in a slow manner, we can set it to a very small number c 3 = 0.01. In this way, we will update Tbest with this Gbest if the result of Gbest is better than the current Tbest; otherwise, Tbest still has the power of pulling the swarm to its position after several evolutionary generations.

Normal and Disturbance Update Modes
Note that Gbest and Pbest are more likely to be closer to the real global optimal solution in some dimensions. On the other hand, a poor particle's current position is far away from the optimal solution, and it requires more cost (e.g., FEs) to move to a promising region.
Based on this consideration, we use different update modes for particles with excellent fitness and poor fitness to speed up the search process. We take the worst N w particles in each iteration as the particles with poor fitness. As shown in Equation (11), N w is proportional to E d . The larger the value of E d is, the more dispersed the particles are. In this way, we force more particles to jump to a promising region to while reducing the waste of the FEs on poor particles.
(1) For particles with better fitness, the update mode follows the usual manner shown in Equations (2 and 10). Note that this update mode is slightly different from the traditional update due to the Tbest term in the velocity update formula. (2) For particles with poor fitness, we use a disturbance update mode as follows.
Step 1: Select the position of Pbest or Gbest as the seed (denoted as P) for generating new positions. Refer to Equation (12), the choice of Pbest position (denoted as P i ) or Gbest position (denoted as P g ) depends on the value of E d . When E d is large, we tend to choose Pbest (i.e., self-cognition) to avoid fast convergence to some local optimum; when E d is small, we choose Gbest (social-cognition) to achieve a finer result.
where r is a random number uniformly distributed in [0, 1].
Step 2: Since P is very close to the optimal solution in some dimensions, we can randomly select a dimension d of P and perturb it into a new position, and keep the other dimensions remain the same thus to maintain the property of P as much as possible.
where r n is a random number and r n ∼ N 0, E 2 d .
Step 3: Replace the position of poor particles with the new position generated by the disturbance update mode.

Experiments and Comparisons
As shown in Table 2, to validate the effectiveness of PSO-ED, we use the same comparison set provided by a recent paper on MPSO [72] and choose 12 representative test functions from CEC2015 benchmark functions [73] as the test suite. The test suite consists of 2 unimodal functions (F 1 and F 2 ), 3 multimodal functions (F 3 , F 4 and F 5 ), 3 hybrid functions (F 6 , F 7 and F 8 ) and 4 composition functions (F 9 , F 10 , F 11 and F 12 ). All the guidelines of CEC2015 have been strictly followed for the experiments. For example, every function in this test suite is conducted for 30 runs independently, the search range of each dimension is set as [−100, 100], and the maximum function evaluation FEs is equal to 10,000 D. The fitness value is F i (x) − F i (x*) after the maximum iteration is reached, where F i (x*) is just a number about the corresponding function for instruction. Hybrid Function (N = 3) 600 7 Hybrid Function (N = 4) 700 8 Hybrid Function (N = 5) 800 9 Composition

Sensitivity Analysis
In this article, two parameters (α and N s ) are very critical to performance of PSO-ED. In order to examine the influence of the parameters on precision and efficiency, a sensitivity analysis was conducted to each parameter.
As shown in Equation (7), α is a user-defined threshold used to judge whether the swarm is evolving into a better result. A proper value of α can reflect the evolutionary state of the swarm and help the swarm choose an appropriate evolution strategy. In order to determine an appropriate value of α, we have selected five values (1.0, 0.1, 0.01, 0.001 and 0.0001) and conducted for 30 runs independently for each value on 12 test functions. Note that the other parameters are set as follows: N s = 40, N total = FEs/N s , N smin = 0.01N total , N smax = 0.1N total , N umax = 0.01N total , v ub = 0.01(x ub − x lb ).
As shown in Table 3, the statistical values (Mean and Rank) are listed and it is obvious that the result is best (smallest average ranking) when α = 0.01. Another important parameter is the swarm size N s . We choose five values (20, 30, 40, 50 and 60) to conduct experiments. And the other parameters are set as follows: N total = FEs/N s , N smin = 0.01N total , N smax = 0.1N total , α = 0.01, N umax = 0.01N total , v ub = 0.01(x ub − x lb ). The results are shown in Table 4. The average ranking shows that PSO-ED performs best when N s = 40.

Comparison with Other PSO Variants
We use the same criteria provided by MPSO [72] to compare our PSO-ED with MPSO and its corresponding comparison algorithms for fair comparison. As shown in on Table 5, these popular algorithms are GPSO [49], LPSO [26], SPSO [74], CLPSO [75], FIPS [27], and DMSPSO [35]. The algorithms' performances are characterized by statistical means and standard deviations. Note that the results of the seven PSO variants are directly referenced from literature [72]. For all the algorithms, the swarm size is set to 40, therefore the maximum number of iterations is 2500 for 10-D and 7500 for 30-D. The statistical results for the 12 test functions in 10-D and 30-D are summarized in Tables 6 and 7, respectively. Following conventions, the lowest mean and standard deviation in each line are regarded as the best results and are highlighted in bold. To validate the significance of the PSO-ED, Non-parametric Wilcoxon signed-rank test is conducted between PSO-ED and other PSO variants. The symbol "+", "0" and "-" means that the statistical mean values of the proposed PSO-ED are better, equal and worse than those of the compared algorithms. Modified particle swarm optimization and its applications MPSO [72] While conducting the comparison experiments, the crucial control parameters of our PSO-ED are set as follows: N s = 40 (same as the setting of MPSO [72]), N total = FEs/N s , N smin = 0.01N total , N smax = 0.1N total , α = 0.01, N umax = 0.01N total , v ub = 0.01(x ub − x lb ). As shown in Table 6 for the 10-D case, compared with the other 7 algorithms, out of the total of 12 functions, the mean values of our algorithm achieve better solutions than the other 7 algorithms on 4 functions (F 1 , F 2 , F 6 and F 11 ), and achieve exactly the same best mean values with other 5 algorithms on 2 functions (F 9 and F 12 ). However, compared with MPSO on F 3 , F 4 , F 5 , DMSPSO on F 8 , FIPS on F 7 , CLPSO on F 10 , our PSO-ED is slightly weaker.  As shown in Table 7 for 30-D case, compared to the counterparts, out of the total of 12 functions, our algorithm obtains the best mean values on 6 functions (F 2 , F 3 , F 6 , F 8 , F 11 and F 12 ). However, PSO-ED is slight weaker than MPSO on F 4 , F 5 , F 8 , SPSO on F 1 .
As shown in Table 8, in the Non-parametric Wilcoxon signed-rank test of PSO-ED and other algorithms, PSO-ED have achieved obvious advantages on 10-D and 30-D except that MPOS and PSO-ED achieved the same optimal result on 10-D.

Conclusions
In this paper, we proposed a novel measure of diversity based on sub-space encoding for the search space; a notion of exploration degree based on the diversity in the exploration, exploitation and convergence states, which efficiently evaluates the degree of demand for the dispersion of the swarm; a technique of disturbance update mode for updating the poor-performed particles' positions to save the cost of function evaluations (FEs) on them. Since the diversity evaluation can reflect the swarm distribution very well, it can provide a better basis for adaptive parameter adjustment strategy and assist the swarm to jump out of the local traps. Therefore, this method is more suitable for the complex multimodal optimization.
The effectiveness of the developed techniques was validated through a set of benchmark functions in CEC2015. Compared with 7 popular PSO variants, out of the 12 benchmark functions, PSO-ED obtains 6 best results for both the 10-D and 30-D cases.
However, the stability of the developed PSO-ED can be further improved and is worthy of investigation in the future work. For example, for F 12 in 10-D and 30-D, although PSO-ED and SPSO both achieve the optimal mean, PSO-ED is slightly weaker in the std term, which means that PSO-ED suffers from risk of falling into some local optimum.